Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:15:28 

0001 // -*- C++ -*-
0002 //
0003 // Package:    TrackerDigiGeometryAnalyzer
0004 // Class:      TrackerDigiGeometryAnalyzer
0005 //
0006 /**\class TrackerDigiGeometryAnalyzer TrackerDigiGeometryAnalyzer.cc 
0007 
0008  Description: <one line class summary>
0009 
0010  Implementation:
0011      <Notes on implementation>
0012 */
0013 //
0014 // Original Author:  Filippo Ambroglini
0015 //         Created:  Tue Jul 26 08:47:57 CEST 2005
0016 //
0017 //
0018 
0019 // system include files
0020 #include <memory>
0021 
0022 // user include files
0023 #include "FWCore/Framework/interface/Frameworkfwd.h"
0024 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0025 
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/Framework/interface/EventSetup.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029 
0030 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0031 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0032 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0033 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
0034 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
0035 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0036 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0037 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetType.h"
0038 
0039 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0040 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0041 #include "DataFormats/GeometrySurface/interface/BoundSurface.h"
0042 #include "DataFormats/GeometrySurface/interface/MediumProperties.h"
0043 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
0044 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0045 
0046 //
0047 //
0048 // class decleration
0049 //
0050 
0051 // #define PRINT(X) edm::LogInfo(X)
0052 #define PRINT(X) std::cout << X << ':'
0053 
0054 class TrackerDigiGeometryAnalyzer : public edm::one::EDAnalyzer<> {
0055 public:
0056   explicit TrackerDigiGeometryAnalyzer(const edm::ParameterSet&);
0057 
0058   void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0059 
0060 private:
0061   void analyseTrapezoidal(const GeomDetUnit& det);
0062   void checkRotation(const GeomDetUnit& det);
0063   void checkTopology(const GeomDetUnit& det);
0064   std::ostream& cylindrical(std::ostream& os, const GlobalPoint& gp) const;
0065 
0066   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> pDDToken_;
0067 };
0068 
0069 TrackerDigiGeometryAnalyzer::TrackerDigiGeometryAnalyzer(const edm::ParameterSet& iConfig) : pDDToken_(esConsumes()) {}
0070 
0071 // ------------ method called to produce the data  ------------
0072 void TrackerDigiGeometryAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0073   PRINT("TrackerDigiGeometryAnalyzer") << "Here I am";
0074   //
0075   // get the TrackerGeom
0076   //
0077   auto const& pDD = iSetup.getData(pDDToken_);
0078   PRINT("TrackerDigiGeometryAnalyzer") << " Geometry node for TrackerGeom is  " << &pDD << '\n';
0079   PRINT("TrackerDigiGeometryAnalyzer") << " I have " << pDD.detUnits().size() << " detectors" << '\n';
0080   PRINT("TrackerDigiGeometryAnalyzer") << " I have " << pDD.detTypes().size() << " types" << '\n';
0081 
0082   for (auto const& it : pDD.detUnits()) {
0083     if (dynamic_cast<const PixelGeomDetUnit*>((it)) != nullptr) {
0084       const BoundPlane& p = (dynamic_cast<const PixelGeomDetUnit*>((it)))->specificSurface();
0085       PRINT("TrackerDigiGeometryAnalyzer") << it->geographicalId() << " RadLeng Pixel " << p.mediumProperties().radLen()
0086                                            << ' ' << " Xi Pixel " << p.mediumProperties().xi() << '\n';
0087     }
0088 
0089     if (dynamic_cast<const StripGeomDetUnit*>((it)) != nullptr) {
0090       const BoundPlane& s = (dynamic_cast<const StripGeomDetUnit*>((it)))->specificSurface();
0091       PRINT("TrackerDigiGeometryAnalyzer") << it->geographicalId() << " RadLeng Strip " << s.mediumProperties().radLen()
0092                                            << " Xi Strip " << s.mediumProperties().xi() << '\n';
0093     }
0094 
0095     //analyseTrapezoidal(*it);
0096   }
0097 
0098   for (auto const& it : pDD.detTypes()) {
0099     if (dynamic_cast<const PixelGeomDetType*>((it)) != nullptr) {
0100       const PixelTopology& p = (dynamic_cast<const PixelGeomDetType*>((it)))->specificTopology();
0101       PRINT("TrackerDigiGeometryAnalyzer") << " PIXEL Det "  // << it->geographicalId()
0102                                            << "    Rows    " << p.nrows() << "    Columns " << p.ncolumns() << '\n';
0103     } else {
0104       const StripTopology& p = (dynamic_cast<const StripGeomDetType*>((it)))->specificTopology();
0105       PRINT("TrackerDigiGeometryAnalyzer") << " STRIP Det "  // << it->geographicalId()
0106                                            << "    Strips    " << p.nstrips() << '\n';
0107     }
0108   }
0109 }
0110 void TrackerDigiGeometryAnalyzer::analyseTrapezoidal(const GeomDetUnit& det) {
0111   // checkRotation( det);
0112 
0113   const double safety = 0.9999;
0114 
0115   const Bounds& bounds = det.surface().bounds();
0116   const TrapezoidalPlaneBounds* tb = dynamic_cast<const TrapezoidalPlaneBounds*>(&bounds);
0117   if (tb == nullptr)
0118     return;  // not trapezoidal
0119 
0120   checkTopology(det);
0121 
0122   const GlobalPoint& pos = det.position();
0123   double length = tb->length();
0124   double width = tb->width();
0125 
0126   const std::array<const float, 4>& par = tb->parameters();
0127   double top = std::max(par[1], par[0]);
0128   double bot = std::min(par[1], par[0]);
0129 
0130   std::cout << std::endl;
0131   std::cout << "Det at pos " << pos << " has length " << length << " width " << width << " pars ";
0132   for (int i = 0; i < 4; i++)
0133     std::cout << par[i] << ", ";
0134   std::cout << std::endl;
0135 
0136   std::cout << "det center inside bounds? " << tb->inside(det.surface().toLocal(pos)) << std::endl;
0137 
0138   //   double outerScale = (pos.perp()+safety*length/2.) / pos.perp();
0139   //   GlobalPoint outerMiddle = GlobalPoint( outerScale*pos.x(), outerScale*pos.y(), pos.z());
0140 
0141   GlobalVector yShift = det.surface().toGlobal(LocalVector(0, safety * length / 2., 0));
0142   GlobalPoint outerMiddle = pos + yShift;
0143   GlobalPoint innerMiddle = pos + (-1. * yShift);
0144   if (outerMiddle.perp() < innerMiddle.perp())
0145     std::swap(outerMiddle, innerMiddle);
0146 
0147   GlobalVector upperShift = det.surface().toGlobal(LocalVector(safety * top, 0, 0));
0148 
0149   GlobalPoint ulc = outerMiddle + upperShift;
0150   GlobalPoint urc = outerMiddle + (-1. * upperShift);
0151   std::cout << "outerMiddle " << outerMiddle << " upperShift " << upperShift << " ulc "
0152             << "(" << ulc.perp() << "," << ulc.phi() << "," << ulc.z() << " urc "
0153             << "(" << urc.perp() << "," << urc.phi() << "," << urc.z() << std::endl;
0154 
0155   std::cout << "upper corners inside bounds? " << tb->inside(det.surface().toLocal(ulc)) << " "
0156             << tb->inside(det.surface().toLocal(urc)) << std::endl;
0157 
0158   //   double innerScale = (pos.perp()-safety*length/2.) / pos.perp();
0159   //   GlobalPoint innerMiddle = GlobalPoint( innerScale*pos.x(), innerScale*pos.y(), pos.z());
0160 
0161   GlobalVector lowerShift = det.surface().toGlobal(LocalVector(safety * bot, 0, 0));
0162 
0163   std::cout << "lower corners inside bounds? " << tb->inside(det.surface().toLocal(innerMiddle + lowerShift)) << " "
0164             << tb->inside(det.surface().toLocal(innerMiddle + (-1. * lowerShift))) << std::endl;
0165 }
0166 
0167 void TrackerDigiGeometryAnalyzer::checkRotation(const GeomDetUnit& det) {
0168   const double eps = std::numeric_limits<float>::epsilon();
0169   static int first = 0;
0170   if (first == 0) {
0171     std::cout << "numeric_limits<float>::epsilon() " << std::numeric_limits<float>::epsilon() << std::endl;
0172     first = 1;
0173   }
0174 
0175   const Surface::RotationType& rot(det.surface().rotation());
0176   GlobalVector a(rot.xx(), rot.xy(), rot.xz());
0177   GlobalVector b(rot.yx(), rot.yy(), rot.yz());
0178   GlobalVector c(rot.zx(), rot.zy(), rot.zz());
0179   GlobalVector cref = a.cross(b);
0180   GlobalVector aref = b.cross(c);
0181   GlobalVector bref = c.cross(a);
0182   if ((a - aref).mag() > eps || (b - bref).mag() > eps || (c - cref).mag() > eps) {
0183     std::cout << " Rotation not good by cross product: " << (a - aref).mag() << ", " << (b - bref).mag() << ", "
0184               << (c - cref).mag() << " for det at pos " << det.surface().position() << std::endl;
0185   }
0186   if (fabs(a.mag() - 1.) > eps || fabs(b.mag() - 1.) > eps || fabs(c.mag() - 1.) > eps) {
0187     std::cout << " Rotation not good by bector mag: " << (a).mag() << ", " << (b).mag() << ", " << (c).mag()
0188               << " for det at pos " << det.surface().position() << std::endl;
0189   }
0190 }
0191 
0192 void TrackerDigiGeometryAnalyzer::checkTopology(const GeomDetUnit& det) {
0193   const StripTopology& topol = dynamic_cast<const StripTopology&>(det.topology());
0194   const int N = 5;
0195 
0196   std::cout << std::endl << "Topology test along strip 0" << std::endl;
0197   LocalVector stripDir = topol.localPosition(MeasurementPoint(0, 0.5)) - topol.localPosition(MeasurementPoint(0, 0));
0198   std::cout << "StripDir " << stripDir << std::endl;
0199   for (int i = -N; i <= N; i++) {
0200     MeasurementPoint mp(0, float(i) / float(2 * N));
0201     LocalPoint lp = topol.localPosition(mp);
0202     double strlen = topol.localStripLength(lp);
0203     LocalError le = topol.localError(mp, MeasurementError(0.25, 0, 0.25));
0204     LocalError lep = le.rotate(stripDir.y(), stripDir.x());
0205     LocalError lem = le.rotate(stripDir.y(), -stripDir.x());
0206     GlobalPoint gp(det.surface().toGlobal(lp));
0207     std::cout << "gpos (r,phi) (" << gp.perp() << ", " << gp.phi() << ") lpos " << lp << " lerr.x (0,+,-) "
0208               << sqrt(le.xx()) << "," << sqrt(lep.xx()) << "," << sqrt(lem.xx()) << " strlen " << strlen << std::endl;
0209   }
0210 
0211   std::cout << std::endl << "Topology test along middle strip" << std::endl;
0212   float midstrip = topol.strip(LocalPoint(0, 0));
0213   for (int i = -N; i <= N; i++) {
0214     MeasurementPoint mp(midstrip, float(i) / float(2 * N));
0215     LocalPoint lp = topol.localPosition(mp);
0216     double strlen = topol.localStripLength(lp);
0217     LocalError le = topol.localError(mp, MeasurementError(0.25, 0, 0.25));
0218     GlobalPoint gp(det.surface().toGlobal(lp));
0219     std::cout << "gpos (r,phi) (" << gp.perp() << ", " << gp.phi() << ") lpos " << lp << " lerr.x " << sqrt(le.xx())
0220               << " strlen " << strlen << std::endl;
0221   }
0222 }
0223 
0224 std::ostream& TrackerDigiGeometryAnalyzer::cylindrical(std::ostream& os, const GlobalPoint& gp) const {
0225   os << "(" << gp.perp() << "," << gp.phi() << "," << gp.z();
0226   return os;
0227 }
0228 
0229 //define this as a plug-in
0230 DEFINE_FWK_MODULE(TrackerDigiGeometryAnalyzer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!