File indexing completed on 2021-12-22 01:53:50
0001
0002
0003
0004
0005
0006
0007 #include <TROOT.h>
0008 #include <TTree.h>
0009 #include <TFile.h>
0010 #include <TF1.h>
0011 #include <TH2F.h>
0012 #include <TH1F.h>
0013
0014
0015 #include <memory>
0016 #include <string>
0017 #include <iostream>
0018
0019
0020 #include "DataFormats/Common/interface/DetSetVector.h"
0021 #include "DataFormats/Common/interface/Handle.h"
0022 #include "DataFormats/Common/interface/OwnVector.h"
0023 #include "DataFormats/Common/interface/Ref.h"
0024 #include "DataFormats/DetId/interface/DetId.h"
0025 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0026 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0027 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0028 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
0029 #include "DataFormats/TrackReco/interface/Track.h"
0030 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0031 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0032 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2DCollection.h"
0033 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2DCollection.h"
0034 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0035 #include "FWCore/Framework/interface/Event.h"
0036 #include "FWCore/Framework/interface/EventSetup.h"
0037 #include "FWCore/Framework/interface/MakerMacros.h"
0038 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0039 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0040 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0041 #include "FWCore/PluginManager/interface/ModuleDef.h"
0042 #include "FWCore/Utilities/interface/InputTag.h"
0043 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0044 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0045 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0046 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0047 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0048 #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyBuilder.h"
0049 #include "Geometry/TrackerGeometryBuilder/interface/RectangularPixelTopology.h"
0050 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0051 #include "MagneticField/Engine/interface/MagneticField.h"
0052 #include "SimDataFormats/Track/interface/SimTrack.h"
0053 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0054 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0055 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0056 #include "SimDataFormats/Vertex/interface/SimVertex.h"
0057 #include "SimDataFormats/Vertex/interface/SimVertexContainer.h"
0058 #include "SimTracker/TrackerHitAssociation/interface/TrackerHitAssociator.h"
0059
0060 class StdHitNtuplizer : public edm::one::EDAnalyzer<> {
0061 public:
0062 explicit StdHitNtuplizer(const edm::ParameterSet& conf);
0063 virtual ~StdHitNtuplizer();
0064 virtual void beginJob();
0065 virtual void endJob();
0066 virtual void analyze(const edm::Event& e, const edm::EventSetup& es);
0067
0068 protected:
0069 void fillEvt(const edm::Event&);
0070 void fillSRecHit(const int subid,
0071 SiStripRecHit2DCollection::DetSet::const_iterator pixeliter,
0072 const GeomDet* theGeom);
0073 void fillSRecHit(const int subid,
0074 SiStripMatchedRecHit2DCollection::DetSet::const_iterator pixeliter,
0075 const GeomDet* theGeom);
0076 void fillSRecHit(const int subid, const FastTrackerRecHit& hit, const GeomDet* theGeom);
0077
0078
0079 void fillPRecHit(const int subid,
0080 const int layer_num,
0081 SiPixelRecHitCollection::DetSet::const_iterator pixeliter,
0082 const int num_simhit,
0083 std::vector<PSimHit>::const_iterator closest_simhit,
0084 const GeomDet* PixGeom);
0085 void fillPRecHit(const int subid, trackingRecHit_iterator pixeliter, const GeomDet* PixGeom);
0086
0087 private:
0088 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geom_esToken;
0089 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topo_esToken;
0090 edm::ParameterSet conf_;
0091 TrackerHitAssociator::Config trackerHitAssociatorConfig_;
0092 edm::InputTag src_;
0093 edm::InputTag rphiRecHits_;
0094 edm::InputTag stereoRecHits_;
0095 edm::InputTag matchedRecHits_;
0096
0097 void init();
0098
0099
0100 struct evt {
0101 int run;
0102 int evtnum;
0103
0104 void init();
0105 } evt_, stripevt_;
0106
0107 struct RecHit {
0108 float x;
0109 float y;
0110 float xx;
0111 float xy;
0112 float yy;
0113 float row;
0114 float col;
0115 float gx;
0116 float gy;
0117 float gz;
0118 int subid;
0119 int layer;
0120 int nsimhit;
0121 float hx, hy;
0122 float tx, ty;
0123 float theta, phi;
0124
0125 void init();
0126 } recHit_, striprecHit_;
0127
0128 TFile* tfile_;
0129 TTree* pixeltree_;
0130 TTree* striptree_;
0131 TTree* pixeltree2_;
0132 };
0133
0134 using namespace std;
0135 using namespace edm;
0136 using namespace reco;
0137
0138 StdHitNtuplizer::StdHitNtuplizer(edm::ParameterSet const& conf)
0139 : geom_esToken(esConsumes()),
0140 topo_esToken(esConsumes()),
0141 conf_(conf),
0142 trackerHitAssociatorConfig_(conf, consumesCollector()),
0143 src_(conf.getParameter<edm::InputTag>("src")),
0144 rphiRecHits_(conf.getParameter<edm::InputTag>("rphiRecHits")),
0145 stereoRecHits_(conf.getParameter<edm::InputTag>("stereoRecHits")),
0146 matchedRecHits_(conf.getParameter<edm::InputTag>("matchedRecHits")),
0147 tfile_(0),
0148 pixeltree_(0),
0149 striptree_(0),
0150 pixeltree2_(0) {}
0151
0152 StdHitNtuplizer::~StdHitNtuplizer() = default;
0153
0154 void StdHitNtuplizer::endJob() {
0155 std::cout << " StdHitNtuplizer::endJob" << std::endl;
0156 tfile_->Write();
0157 tfile_->Close();
0158 }
0159
0160 void StdHitNtuplizer::beginJob() {
0161 std::cout << " StdHitNtuplizer::beginJob" << std::endl;
0162 std::string outputFile = conf_.getParameter<std::string>("OutputFile");
0163
0164 tfile_ = new TFile(outputFile.c_str(), "RECREATE");
0165 pixeltree_ = new TTree("PixelNtuple", "Pixel hit analyzer ntuple");
0166 striptree_ = new TTree("StripNtuple", "Strip hit analyzer ntuple");
0167 pixeltree2_ = new TTree("Pixel2Ntuple", "Track Pixel hit analyzer ntuple");
0168
0169 int bufsize = 64000;
0170
0171
0172 pixeltree_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0173 pixeltree_->Branch("pixel_recHit",
0174 &recHit_,
0175 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0176 bufsize);
0177 pixeltree2_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0178 pixeltree2_->Branch("pixel_recHit",
0179 &recHit_,
0180 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0181 bufsize);
0182
0183
0184 striptree_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0185 striptree_->Branch("strip_recHit",
0186 &striprecHit_,
0187 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0188 bufsize);
0189 }
0190
0191
0192 void StdHitNtuplizer::analyze(const edm::Event& e, const edm::EventSetup& es) {
0193
0194 const TrackerTopology* const tTopo = &es.getData(topo_esToken);
0195
0196
0197 const TrackerGeometry* theGeometry = &es.getData(geom_esToken);
0198
0199
0200
0201
0202
0203
0204
0205 edm::Handle<SiPixelRecHitCollection> recHitColl;
0206 e.getByLabel(src_, recHitColl);
0207
0208
0209 TrackerHitAssociator associate(e, trackerHitAssociatorConfig_);
0210
0211
0212 if ((recHitColl.product())->dataSize() > 0) {
0213 SiPixelRecHitCollection::const_iterator recHitIdIterator = (recHitColl.product())->begin();
0214 SiPixelRecHitCollection::const_iterator recHitIdIteratorEnd = (recHitColl.product())->end();
0215
0216 std::string detname;
0217 std::vector<PSimHit> matched;
0218 std::vector<PSimHit>::const_iterator closest_simhit;
0219
0220
0221 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0222 SiPixelRecHitCollection::DetSet detset = *recHitIdIterator;
0223
0224 if (detset.empty())
0225 continue;
0226 DetId detId = DetId(detset.detId());
0227
0228
0229 const GeomDet* geomDet(theGeometry->idToDet(detId));
0230
0231
0232 SiPixelRecHitCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0233 SiPixelRecHitCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0234 SiPixelRecHitCollection::DetSet::const_iterator iterRecHit;
0235 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0236
0237 matched.clear();
0238 matched = associate.associateHit(*iterRecHit);
0239 if (!matched.empty()) {
0240 float closest = 9999.9;
0241 std::vector<PSimHit>::const_iterator closestit = matched.begin();
0242 LocalPoint lp = iterRecHit->localPosition();
0243 float rechit_x = lp.x();
0244 float rechit_y = lp.y();
0245
0246 for (std::vector<PSimHit>::const_iterator m = matched.begin(); m < matched.end(); m++) {
0247 float sim_x1 = (*m).entryPoint().x();
0248 float sim_x2 = (*m).exitPoint().x();
0249 float sim_xpos = 0.5 * (sim_x1 + sim_x2);
0250 float sim_y1 = (*m).entryPoint().y();
0251 float sim_y2 = (*m).exitPoint().y();
0252 float sim_ypos = 0.5 * (sim_y1 + sim_y2);
0253
0254 float x_res = fabs(sim_xpos - rechit_x);
0255 float y_res = fabs(sim_ypos - rechit_y);
0256 float dist = sqrt(x_res * x_res + y_res * y_res);
0257 if (dist < closest) {
0258 closest = dist;
0259 closestit = m;
0260 }
0261 }
0262 closest_simhit = closestit;
0263 }
0264
0265
0266
0267
0268
0269
0270
0271
0272
0273
0274
0275
0276
0277
0278
0279
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289
0290
0291
0292
0293
0294
0295
0296
0297
0298
0299
0300
0301
0302
0303
0304
0305
0306
0307
0308
0309
0310
0311
0312
0313
0314
0315
0316 unsigned int subid = detId.subdetId();
0317 int layer_num = 0;
0318 if ((subid == 1) || (subid == 2)) {
0319
0320 if (subid == PixelSubdetector::PixelBarrel) {
0321 layer_num = tTopo->pxbLayer(detId.rawId());
0322 } else if (subid == PixelSubdetector::PixelEndcap) {
0323 layer_num = tTopo->pxfDisk(detId.rawId());
0324 }
0325 int num_simhit = matched.size();
0326 fillPRecHit(subid, layer_num, iterRecHit, num_simhit, closest_simhit, geomDet);
0327 fillEvt(e);
0328 pixeltree_->Fill();
0329 init();
0330
0331
0332
0333
0334
0335
0336
0337
0338
0339
0340
0341
0342
0343
0344
0345 }
0346 }
0347 }
0348 }
0349
0350
0351
0352 edm::Handle<View<reco::Track> > trackCollection;
0353 edm::InputTag trackProducer;
0354 trackProducer = conf_.getParameter<edm::InputTag>("trackProducer");
0355 e.getByLabel(trackProducer, trackCollection);
0356
0357
0358
0359
0360
0361
0362
0363
0364 int rT = 0;
0365 for (View<reco::Track>::size_type i = 0; i < trackCollection->size(); ++i) {
0366 ++rT;
0367 RefToBase<reco::Track> track(trackCollection, i);
0368
0369 for (trackingRecHit_iterator ih = track->recHitsBegin(); ih != track->recHitsEnd(); ++ih) {
0370 TrackingRecHit* hit = (*ih)->clone();
0371 const DetId& detId = hit->geographicalId();
0372 const GeomDet* geomDet(theGeometry->idToDet(detId));
0373
0374
0375
0376
0377
0378
0379
0380
0381
0382
0383
0384
0385
0386
0387
0388
0389
0390
0391
0392
0393
0394
0395
0396
0397
0398
0399
0400
0401
0402
0403
0404
0405
0406
0407
0408
0409
0410
0411
0412
0413
0414
0415
0416
0417
0418
0419
0420 if (hit->isValid()) {
0421 unsigned int subid = detId.subdetId();
0422 if ((subid == 1) || (subid == 2)) {
0423
0424 fillPRecHit(subid, ih, geomDet);
0425 fillEvt(e);
0426 pixeltree2_->Fill();
0427 init();
0428
0429
0430
0431
0432
0433
0434
0435
0436
0437
0438
0439
0440
0441
0442
0443
0444 }
0445 }
0446 delete hit;
0447 }
0448 }
0449
0450
0451 edm::Handle<SiStripRecHit2DCollection> rechitsrphi;
0452 edm::Handle<SiStripRecHit2DCollection> rechitsstereo;
0453 edm::Handle<SiStripMatchedRecHit2DCollection> rechitsmatched;
0454 e.getByLabel(rphiRecHits_, rechitsrphi);
0455 e.getByLabel(stereoRecHits_, rechitsstereo);
0456 e.getByLabel(matchedRecHits_, rechitsmatched);
0457
0458
0459
0460
0461 if (rechitsrphi->size() > 0) {
0462
0463
0464
0465
0466 SiStripRecHit2DCollection::const_iterator recHitIdIterator = (rechitsrphi.product())->begin();
0467 SiStripRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsrphi.product())->end();
0468
0469 std::string detname;
0470
0471
0472
0473
0474 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0475 SiStripRecHit2DCollection::DetSet detset = *recHitIdIterator;
0476
0477 if (detset.empty())
0478 continue;
0479 DetId detId = DetId(detset.detId());
0480
0481
0482 const GeomDet* geomDet(theGeometry->idToDet(detId));
0483
0484
0485 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0486 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0487 SiStripRecHit2DCollection::DetSet::const_iterator iterRecHit;
0488 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0489
0490
0491
0492
0493
0494
0495
0496
0497
0498
0499
0500
0501
0502
0503
0504
0505
0506
0507
0508
0509
0510
0511
0512
0513
0514
0515
0516
0517
0518
0519
0520
0521
0522
0523
0524
0525
0526
0527
0528
0529
0530
0531
0532
0533
0534
0535
0536
0537
0538
0539
0540 unsigned int subid = detId.subdetId();
0541 fillSRecHit(subid, iterRecHit, geomDet);
0542 fillEvt(e);
0543 striptree_->Fill();
0544 init();
0545 }
0546 }
0547 }
0548
0549
0550 if (rechitsstereo.product()->dataSize() > 0) {
0551
0552
0553
0554
0555 SiStripRecHit2DCollection::const_iterator recHitIdIterator = (rechitsstereo.product())->begin();
0556 SiStripRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsstereo.product())->end();
0557
0558 std::string detname;
0559
0560
0561 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0562 SiStripRecHit2DCollection::DetSet detset = *recHitIdIterator;
0563
0564
0565
0566 if (detset.empty())
0567 continue;
0568 DetId detId = DetId(detset.detId());
0569
0570
0571 const GeomDet* geomDet(theGeometry->idToDet(detId));
0572
0573
0574 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0575 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0576 SiStripRecHit2DCollection::DetSet::const_iterator iterRecHit;
0577 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0578
0579
0580
0581
0582
0583
0584
0585
0586
0587
0588
0589
0590
0591
0592
0593
0594
0595
0596
0597
0598
0599
0600
0601
0602
0603
0604
0605
0606
0607
0608
0609
0610
0611
0612
0613
0614
0615
0616
0617
0618
0619
0620
0621
0622
0623
0624
0625
0626
0627
0628
0629 unsigned int subid = detId.subdetId();
0630 fillSRecHit(subid, iterRecHit, geomDet);
0631 fillEvt(e);
0632 striptree_->Fill();
0633 init();
0634 }
0635 }
0636 }
0637
0638
0639 if (rechitsmatched.product()->dataSize() > 0) {
0640
0641
0642
0643
0644 SiStripMatchedRecHit2DCollection::const_iterator recHitIdIterator = (rechitsmatched.product())->begin();
0645 SiStripMatchedRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsmatched.product())->end();
0646
0647 std::string detname;
0648
0649
0650 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0651 SiStripMatchedRecHit2DCollection::DetSet detset = *recHitIdIterator;
0652
0653 if (detset.empty())
0654 continue;
0655 DetId detId = DetId(detset.detId());
0656
0657
0658
0659
0660
0661 const GeomDet* geomDet(theGeometry->idToDet(detId));
0662
0663
0664 SiStripMatchedRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0665 SiStripMatchedRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0666 SiStripMatchedRecHit2DCollection::DetSet::const_iterator iterRecHit;
0667 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0668
0669
0670
0671
0672
0673
0674
0675
0676
0677
0678
0679
0680
0681
0682
0683
0684
0685
0686
0687
0688
0689
0690
0691
0692
0693
0694
0695
0696
0697
0698
0699
0700
0701
0702
0703
0704
0705
0706
0707
0708
0709
0710
0711
0712
0713
0714
0715
0716
0717
0718
0719
0720 unsigned int subid = detId.subdetId();
0721 fillSRecHit(subid, iterRecHit, geomDet);
0722 fillEvt(e);
0723 striptree_->Fill();
0724 init();
0725 }
0726 }
0727 }
0728
0729 }
0730
0731 void StdHitNtuplizer::fillSRecHit(const int subid,
0732 SiStripRecHit2DCollection::DetSet::const_iterator pixeliter,
0733 const GeomDet* theGeom) {
0734 LocalPoint lp = pixeliter->localPosition();
0735 LocalError le = pixeliter->localPositionError();
0736
0737 striprecHit_.x = lp.x();
0738 striprecHit_.y = lp.y();
0739 striprecHit_.xx = le.xx();
0740 striprecHit_.xy = le.xy();
0741 striprecHit_.yy = le.yy();
0742 GlobalPoint GP = theGeom->surface().toGlobal(pixeliter->localPosition());
0743 striprecHit_.gx = GP.x();
0744 striprecHit_.gy = GP.y();
0745 striprecHit_.gz = GP.z();
0746 striprecHit_.subid = subid;
0747 }
0748 void StdHitNtuplizer::fillSRecHit(const int subid,
0749 SiStripMatchedRecHit2DCollection::DetSet::const_iterator pixeliter,
0750 const GeomDet* theGeom) {
0751 LocalPoint lp = pixeliter->localPosition();
0752 LocalError le = pixeliter->localPositionError();
0753
0754 striprecHit_.x = lp.x();
0755 striprecHit_.y = lp.y();
0756 striprecHit_.xx = le.xx();
0757 striprecHit_.xy = le.xy();
0758 striprecHit_.yy = le.yy();
0759 GlobalPoint GP = theGeom->surface().toGlobal(pixeliter->localPosition());
0760 striprecHit_.gx = GP.x();
0761 striprecHit_.gy = GP.y();
0762 striprecHit_.gz = GP.z();
0763 striprecHit_.subid = subid;
0764 }
0765 void StdHitNtuplizer::fillSRecHit(const int subid, const FastTrackerRecHit& hit, const GeomDet* theGeom) {
0766 LocalPoint lp = hit.localPosition();
0767 LocalError le = hit.localPositionError();
0768
0769 striprecHit_.x = lp.x();
0770 striprecHit_.y = lp.y();
0771 striprecHit_.xx = le.xx();
0772 striprecHit_.xy = le.xy();
0773 striprecHit_.yy = le.yy();
0774
0775
0776
0777 GlobalPoint GP = theGeom->surface().toGlobal(hit.localPosition());
0778 striprecHit_.gx = GP.x();
0779 striprecHit_.gy = GP.y();
0780 striprecHit_.gz = GP.z();
0781 striprecHit_.subid = subid;
0782 }
0783 void StdHitNtuplizer::fillPRecHit(const int subid,
0784 const int layer_num,
0785 SiPixelRecHitCollection::DetSet::const_iterator pixeliter,
0786 const int num_simhit,
0787 std::vector<PSimHit>::const_iterator closest_simhit,
0788 const GeomDet* PixGeom) {
0789 LocalPoint lp = pixeliter->localPosition();
0790 LocalError le = pixeliter->localPositionError();
0791
0792 recHit_.x = lp.x();
0793 recHit_.y = lp.y();
0794 recHit_.xx = le.xx();
0795 recHit_.xy = le.xy();
0796 recHit_.yy = le.yy();
0797
0798
0799
0800 GlobalPoint GP = PixGeom->surface().toGlobal(pixeliter->localPosition());
0801 recHit_.gx = GP.x();
0802 recHit_.gy = GP.y();
0803 recHit_.gz = GP.z();
0804 recHit_.subid = subid;
0805 recHit_.layer = layer_num;
0806 recHit_.nsimhit = num_simhit;
0807
0808 if (num_simhit > 0) {
0809 float sim_x1 = (*closest_simhit).entryPoint().x();
0810 float sim_x2 = (*closest_simhit).exitPoint().x();
0811 recHit_.hx = 0.5 * (sim_x1 + sim_x2);
0812 float sim_y1 = (*closest_simhit).entryPoint().y();
0813 float sim_y2 = (*closest_simhit).exitPoint().y();
0814 recHit_.hy = 0.5 * (sim_y1 + sim_y2);
0815
0816 }
0817
0818
0819
0820
0821
0822
0823
0824 }
0825 void StdHitNtuplizer::fillPRecHit(const int subid, trackingRecHit_iterator ih, const GeomDet* PixGeom) {
0826 TrackingRecHit* pixeliter = (*ih)->clone();
0827 LocalPoint lp = pixeliter->localPosition();
0828 LocalError le = pixeliter->localPositionError();
0829
0830 recHit_.x = lp.x();
0831 recHit_.y = lp.y();
0832 recHit_.xx = le.xx();
0833 recHit_.xy = le.xy();
0834 recHit_.yy = le.yy();
0835 GlobalPoint GP = PixGeom->surface().toGlobal(pixeliter->localPosition());
0836 recHit_.gx = GP.x();
0837 recHit_.gy = GP.y();
0838 recHit_.gz = GP.z();
0839 delete pixeliter;
0840 recHit_.subid = subid;
0841 }
0842
0843 void StdHitNtuplizer::fillEvt(const edm::Event& E) {
0844 evt_.run = E.id().run();
0845 evt_.evtnum = E.id().event();
0846 }
0847
0848 void StdHitNtuplizer::init() {
0849 evt_.init();
0850 recHit_.init();
0851 striprecHit_.init();
0852 }
0853
0854 void StdHitNtuplizer::evt::init() {
0855 int dummy_int = 9999;
0856 run = dummy_int;
0857 evtnum = dummy_int;
0858 }
0859
0860 void StdHitNtuplizer::RecHit::init() {
0861 float dummy_float = 9999.0;
0862
0863 x = dummy_float;
0864 y = dummy_float;
0865 xx = dummy_float;
0866 xy = dummy_float;
0867 yy = dummy_float;
0868 row = dummy_float;
0869 col = dummy_float;
0870 gx = dummy_float;
0871 gy = dummy_float;
0872 gz = dummy_float;
0873 layer = 0;
0874 nsimhit = 0;
0875 hx = dummy_float;
0876 hy = dummy_float;
0877 tx = dummy_float;
0878 ty = dummy_float;
0879 theta = dummy_float;
0880 phi = dummy_float;
0881 }
0882
0883
0884 DEFINE_FWK_MODULE(StdHitNtuplizer);