File indexing completed on 2024-09-07 04:38:11
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/Ref.h"
0023 #include "DataFormats/DetId/interface/DetId.h"
0024 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0025 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0026 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0027 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
0028 #include "DataFormats/TrackReco/interface/Track.h"
0029 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0030 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
0031 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2DCollection.h"
0032 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2DCollection.h"
0033 #include "DataFormats/TrajectorySeed/interface/TrajectorySeedCollection.h"
0034 #include "FWCore/Framework/interface/Event.h"
0035 #include "FWCore/Framework/interface/EventSetup.h"
0036 #include "FWCore/Framework/interface/MakerMacros.h"
0037 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0038 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "FWCore/PluginManager/interface/ModuleDef.h"
0041 #include "FWCore/Utilities/interface/InputTag.h"
0042 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0043 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0044 #include "Geometry/CommonDetUnit/interface/PixelGeomDetType.h"
0045 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0046 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0047 #include "Geometry/TrackerGeometryBuilder/interface/PixelTopologyBuilder.h"
0048 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0049 #include "MagneticField/Engine/interface/MagneticField.h"
0050 #include "SimDataFormats/Track/interface/SimTrack.h"
0051 #include "SimDataFormats/Track/interface/SimTrackContainer.h"
0052 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0053 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0054 #include "SimDataFormats/Vertex/interface/SimVertex.h"
0055 #include "SimDataFormats/Vertex/interface/SimVertexContainer.h"
0056 #include "SimTracker/TrackerHitAssociation/interface/TrackerHitAssociator.h"
0057
0058 class StdHitNtuplizer : public edm::one::EDAnalyzer<> {
0059 public:
0060 explicit StdHitNtuplizer(const edm::ParameterSet& conf);
0061 virtual ~StdHitNtuplizer();
0062 virtual void beginJob();
0063 virtual void endJob();
0064 virtual void analyze(const edm::Event& e, const edm::EventSetup& es);
0065
0066 protected:
0067 void fillEvt(const edm::Event&);
0068 void fillSRecHit(const int subid,
0069 SiStripRecHit2DCollection::DetSet::const_iterator pixeliter,
0070 const GeomDet* theGeom);
0071 void fillSRecHit(const int subid,
0072 SiStripMatchedRecHit2DCollection::DetSet::const_iterator pixeliter,
0073 const GeomDet* theGeom);
0074 void fillSRecHit(const int subid, const FastTrackerRecHit& hit, const GeomDet* theGeom);
0075
0076
0077 void fillPRecHit(const int subid,
0078 const int layer_num,
0079 SiPixelRecHitCollection::DetSet::const_iterator pixeliter,
0080 const int num_simhit,
0081 std::vector<PSimHit>::const_iterator closest_simhit,
0082 const GeomDet* PixGeom);
0083 void fillPRecHit(const int subid, trackingRecHit_iterator pixeliter, const GeomDet* PixGeom);
0084
0085 private:
0086 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> geom_esToken;
0087 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> topo_esToken;
0088 edm::ParameterSet conf_;
0089 TrackerHitAssociator::Config trackerHitAssociatorConfig_;
0090 edm::InputTag src_;
0091 edm::InputTag rphiRecHits_;
0092 edm::InputTag stereoRecHits_;
0093 edm::InputTag matchedRecHits_;
0094
0095 void init();
0096
0097
0098 struct evt {
0099 int run;
0100 int evtnum;
0101
0102 void init();
0103 } evt_, stripevt_;
0104
0105 struct RecHit {
0106 float x;
0107 float y;
0108 float xx;
0109 float xy;
0110 float yy;
0111 float row;
0112 float col;
0113 float gx;
0114 float gy;
0115 float gz;
0116 int subid;
0117 int layer;
0118 int nsimhit;
0119 float hx, hy;
0120 float tx, ty;
0121 float theta, phi;
0122
0123 void init();
0124 } recHit_, striprecHit_;
0125
0126 TFile* tfile_;
0127 TTree* pixeltree_;
0128 TTree* striptree_;
0129 TTree* pixeltree2_;
0130 };
0131
0132 using namespace std;
0133 using namespace edm;
0134 using namespace reco;
0135
0136 StdHitNtuplizer::StdHitNtuplizer(edm::ParameterSet const& conf)
0137 : geom_esToken(esConsumes()),
0138 topo_esToken(esConsumes()),
0139 conf_(conf),
0140 trackerHitAssociatorConfig_(conf, consumesCollector()),
0141 src_(conf.getParameter<edm::InputTag>("src")),
0142 rphiRecHits_(conf.getParameter<edm::InputTag>("rphiRecHits")),
0143 stereoRecHits_(conf.getParameter<edm::InputTag>("stereoRecHits")),
0144 matchedRecHits_(conf.getParameter<edm::InputTag>("matchedRecHits")),
0145 tfile_(0),
0146 pixeltree_(0),
0147 striptree_(0),
0148 pixeltree2_(0) {}
0149
0150 StdHitNtuplizer::~StdHitNtuplizer() = default;
0151
0152 void StdHitNtuplizer::endJob() {
0153 std::cout << " StdHitNtuplizer::endJob" << std::endl;
0154 tfile_->Write();
0155 tfile_->Close();
0156 }
0157
0158 void StdHitNtuplizer::beginJob() {
0159 std::cout << " StdHitNtuplizer::beginJob" << std::endl;
0160 std::string outputFile = conf_.getParameter<std::string>("OutputFile");
0161
0162 tfile_ = new TFile(outputFile.c_str(), "RECREATE");
0163 pixeltree_ = new TTree("PixelNtuple", "Pixel hit analyzer ntuple");
0164 striptree_ = new TTree("StripNtuple", "Strip hit analyzer ntuple");
0165 pixeltree2_ = new TTree("Pixel2Ntuple", "Track Pixel hit analyzer ntuple");
0166
0167 int bufsize = 64000;
0168
0169
0170 pixeltree_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0171 pixeltree_->Branch("pixel_recHit",
0172 &recHit_,
0173 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0174 bufsize);
0175 pixeltree2_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0176 pixeltree2_->Branch("pixel_recHit",
0177 &recHit_,
0178 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0179 bufsize);
0180
0181
0182 striptree_->Branch("evt", &evt_, "run/I:evtnum/I", bufsize);
0183 striptree_->Branch("strip_recHit",
0184 &striprecHit_,
0185 "x/F:y:xx:xy:yy:row:col:gx:gy:gz:subid/I:layer:nsimhit:hx/F:hy:tx:ty:theta:phi",
0186 bufsize);
0187 }
0188
0189
0190 void StdHitNtuplizer::analyze(const edm::Event& e, const edm::EventSetup& es) {
0191
0192 const TrackerTopology* const tTopo = &es.getData(topo_esToken);
0193
0194
0195 const TrackerGeometry* theGeometry = &es.getData(geom_esToken);
0196
0197
0198
0199
0200
0201
0202
0203 edm::Handle<SiPixelRecHitCollection> recHitColl;
0204 e.getByLabel(src_, recHitColl);
0205
0206
0207 TrackerHitAssociator associate(e, trackerHitAssociatorConfig_);
0208
0209
0210 if ((recHitColl.product())->dataSize() > 0) {
0211 SiPixelRecHitCollection::const_iterator recHitIdIterator = (recHitColl.product())->begin();
0212 SiPixelRecHitCollection::const_iterator recHitIdIteratorEnd = (recHitColl.product())->end();
0213
0214 std::string detname;
0215 std::vector<PSimHit> matched;
0216 std::vector<PSimHit>::const_iterator closest_simhit;
0217
0218
0219 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0220 SiPixelRecHitCollection::DetSet detset = *recHitIdIterator;
0221
0222 if (detset.empty())
0223 continue;
0224 DetId detId = DetId(detset.detId());
0225
0226
0227 const GeomDet* geomDet(theGeometry->idToDet(detId));
0228
0229
0230 SiPixelRecHitCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0231 SiPixelRecHitCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0232 SiPixelRecHitCollection::DetSet::const_iterator iterRecHit;
0233 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0234
0235 matched.clear();
0236 matched = associate.associateHit(*iterRecHit);
0237 if (!matched.empty()) {
0238 float closest = 9999.9;
0239 std::vector<PSimHit>::const_iterator closestit = matched.begin();
0240 LocalPoint lp = iterRecHit->localPosition();
0241 float rechit_x = lp.x();
0242 float rechit_y = lp.y();
0243
0244 for (std::vector<PSimHit>::const_iterator m = matched.begin(); m < matched.end(); m++) {
0245 float sim_x1 = (*m).entryPoint().x();
0246 float sim_x2 = (*m).exitPoint().x();
0247 float sim_xpos = 0.5 * (sim_x1 + sim_x2);
0248 float sim_y1 = (*m).entryPoint().y();
0249 float sim_y2 = (*m).exitPoint().y();
0250 float sim_ypos = 0.5 * (sim_y1 + sim_y2);
0251
0252 float x_res = fabs(sim_xpos - rechit_x);
0253 float y_res = fabs(sim_ypos - rechit_y);
0254 float dist = sqrt(x_res * x_res + y_res * y_res);
0255 if (dist < closest) {
0256 closest = dist;
0257 closestit = m;
0258 }
0259 }
0260 closest_simhit = closestit;
0261 }
0262
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 unsigned int subid = detId.subdetId();
0315 int layer_num = 0;
0316 if ((subid == 1) || (subid == 2)) {
0317
0318 if (subid == PixelSubdetector::PixelBarrel) {
0319 layer_num = tTopo->pxbLayer(detId.rawId());
0320 } else if (subid == PixelSubdetector::PixelEndcap) {
0321 layer_num = tTopo->pxfDisk(detId.rawId());
0322 }
0323 int num_simhit = matched.size();
0324 fillPRecHit(subid, layer_num, iterRecHit, num_simhit, closest_simhit, geomDet);
0325 fillEvt(e);
0326 pixeltree_->Fill();
0327 init();
0328
0329
0330
0331
0332
0333
0334
0335
0336
0337
0338
0339
0340
0341
0342
0343 }
0344 }
0345 }
0346 }
0347
0348
0349
0350 edm::Handle<View<reco::Track> > trackCollection;
0351 edm::InputTag trackProducer;
0352 trackProducer = conf_.getParameter<edm::InputTag>("trackProducer");
0353 e.getByLabel(trackProducer, trackCollection);
0354
0355
0356
0357
0358
0359
0360
0361
0362 for (View<reco::Track>::size_type i = 0; i < trackCollection->size(); ++i) {
0363 RefToBase<reco::Track> track(trackCollection, i);
0364 for (trackingRecHit_iterator ih = track->recHitsBegin(); ih != track->recHitsEnd(); ++ih) {
0365 TrackingRecHit* hit = (*ih)->clone();
0366 const DetId& detId = hit->geographicalId();
0367 const GeomDet* geomDet(theGeometry->idToDet(detId));
0368
0369
0370
0371
0372
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 if (hit->isValid()) {
0416 unsigned int subid = detId.subdetId();
0417 if ((subid == 1) || (subid == 2)) {
0418
0419 fillPRecHit(subid, ih, geomDet);
0420 fillEvt(e);
0421 pixeltree2_->Fill();
0422 init();
0423
0424
0425
0426
0427
0428
0429
0430
0431
0432
0433
0434
0435
0436
0437
0438
0439 }
0440 }
0441 delete hit;
0442 }
0443 }
0444
0445
0446 edm::Handle<SiStripRecHit2DCollection> rechitsrphi;
0447 edm::Handle<SiStripRecHit2DCollection> rechitsstereo;
0448 edm::Handle<SiStripMatchedRecHit2DCollection> rechitsmatched;
0449 e.getByLabel(rphiRecHits_, rechitsrphi);
0450 e.getByLabel(stereoRecHits_, rechitsstereo);
0451 e.getByLabel(matchedRecHits_, rechitsmatched);
0452
0453
0454
0455
0456 if (rechitsrphi->size() > 0) {
0457
0458
0459
0460
0461 SiStripRecHit2DCollection::const_iterator recHitIdIterator = (rechitsrphi.product())->begin();
0462 SiStripRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsrphi.product())->end();
0463
0464 std::string detname;
0465
0466
0467
0468
0469 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0470 SiStripRecHit2DCollection::DetSet detset = *recHitIdIterator;
0471
0472 if (detset.empty())
0473 continue;
0474 DetId detId = DetId(detset.detId());
0475
0476
0477 const GeomDet* geomDet(theGeometry->idToDet(detId));
0478
0479
0480 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0481 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0482 SiStripRecHit2DCollection::DetSet::const_iterator iterRecHit;
0483 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0484
0485
0486
0487
0488
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 unsigned int subid = detId.subdetId();
0536 fillSRecHit(subid, iterRecHit, geomDet);
0537 fillEvt(e);
0538 striptree_->Fill();
0539 init();
0540 }
0541 }
0542 }
0543
0544
0545 if (rechitsstereo.product()->dataSize() > 0) {
0546
0547
0548
0549
0550 SiStripRecHit2DCollection::const_iterator recHitIdIterator = (rechitsstereo.product())->begin();
0551 SiStripRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsstereo.product())->end();
0552
0553 std::string detname;
0554
0555
0556 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0557 SiStripRecHit2DCollection::DetSet detset = *recHitIdIterator;
0558
0559
0560
0561 if (detset.empty())
0562 continue;
0563 DetId detId = DetId(detset.detId());
0564
0565
0566 const GeomDet* geomDet(theGeometry->idToDet(detId));
0567
0568
0569 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0570 SiStripRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0571 SiStripRecHit2DCollection::DetSet::const_iterator iterRecHit;
0572 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0573
0574
0575
0576
0577
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 unsigned int subid = detId.subdetId();
0625 fillSRecHit(subid, iterRecHit, geomDet);
0626 fillEvt(e);
0627 striptree_->Fill();
0628 init();
0629 }
0630 }
0631 }
0632
0633
0634 if (rechitsmatched.product()->dataSize() > 0) {
0635
0636 SiStripMatchedRecHit2DCollection::const_iterator recHitIdIterator = (rechitsmatched.product())->begin();
0637 SiStripMatchedRecHit2DCollection::const_iterator recHitIdIteratorEnd = (rechitsmatched.product())->end();
0638
0639 std::string detname;
0640
0641
0642 for (; recHitIdIterator != recHitIdIteratorEnd; recHitIdIterator++) {
0643 SiStripMatchedRecHit2DCollection::DetSet detset = *recHitIdIterator;
0644
0645 if (detset.empty())
0646 continue;
0647 DetId detId = DetId(detset.detId());
0648
0649
0650
0651
0652
0653 const GeomDet* geomDet(theGeometry->idToDet(detId));
0654
0655
0656 SiStripMatchedRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorBegin = detset.begin();
0657 SiStripMatchedRecHit2DCollection::DetSet::const_iterator rechitRangeIteratorEnd = detset.end();
0658 SiStripMatchedRecHit2DCollection::DetSet::const_iterator iterRecHit;
0659 for (iterRecHit = rechitRangeIteratorBegin; iterRecHit != rechitRangeIteratorEnd; ++iterRecHit) {
0660
0661
0662
0663
0664
0665
0666
0667
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 unsigned int subid = detId.subdetId();
0713 fillSRecHit(subid, iterRecHit, geomDet);
0714 fillEvt(e);
0715 striptree_->Fill();
0716 init();
0717 }
0718 }
0719 }
0720
0721 }
0722
0723 void StdHitNtuplizer::fillSRecHit(const int subid,
0724 SiStripRecHit2DCollection::DetSet::const_iterator pixeliter,
0725 const GeomDet* theGeom) {
0726 LocalPoint lp = pixeliter->localPosition();
0727 LocalError le = pixeliter->localPositionError();
0728
0729 striprecHit_.x = lp.x();
0730 striprecHit_.y = lp.y();
0731 striprecHit_.xx = le.xx();
0732 striprecHit_.xy = le.xy();
0733 striprecHit_.yy = le.yy();
0734 GlobalPoint GP = theGeom->surface().toGlobal(pixeliter->localPosition());
0735 striprecHit_.gx = GP.x();
0736 striprecHit_.gy = GP.y();
0737 striprecHit_.gz = GP.z();
0738 striprecHit_.subid = subid;
0739 }
0740 void StdHitNtuplizer::fillSRecHit(const int subid,
0741 SiStripMatchedRecHit2DCollection::DetSet::const_iterator pixeliter,
0742 const GeomDet* theGeom) {
0743 LocalPoint lp = pixeliter->localPosition();
0744 LocalError le = pixeliter->localPositionError();
0745
0746 striprecHit_.x = lp.x();
0747 striprecHit_.y = lp.y();
0748 striprecHit_.xx = le.xx();
0749 striprecHit_.xy = le.xy();
0750 striprecHit_.yy = le.yy();
0751 GlobalPoint GP = theGeom->surface().toGlobal(pixeliter->localPosition());
0752 striprecHit_.gx = GP.x();
0753 striprecHit_.gy = GP.y();
0754 striprecHit_.gz = GP.z();
0755 striprecHit_.subid = subid;
0756 }
0757 void StdHitNtuplizer::fillSRecHit(const int subid, const FastTrackerRecHit& hit, const GeomDet* theGeom) {
0758 LocalPoint lp = hit.localPosition();
0759 LocalError le = hit.localPositionError();
0760
0761 striprecHit_.x = lp.x();
0762 striprecHit_.y = lp.y();
0763 striprecHit_.xx = le.xx();
0764 striprecHit_.xy = le.xy();
0765 striprecHit_.yy = le.yy();
0766
0767
0768
0769 GlobalPoint GP = theGeom->surface().toGlobal(hit.localPosition());
0770 striprecHit_.gx = GP.x();
0771 striprecHit_.gy = GP.y();
0772 striprecHit_.gz = GP.z();
0773 striprecHit_.subid = subid;
0774 }
0775 void StdHitNtuplizer::fillPRecHit(const int subid,
0776 const int layer_num,
0777 SiPixelRecHitCollection::DetSet::const_iterator pixeliter,
0778 const int num_simhit,
0779 std::vector<PSimHit>::const_iterator closest_simhit,
0780 const GeomDet* PixGeom) {
0781 LocalPoint lp = pixeliter->localPosition();
0782 LocalError le = pixeliter->localPositionError();
0783
0784 recHit_.x = lp.x();
0785 recHit_.y = lp.y();
0786 recHit_.xx = le.xx();
0787 recHit_.xy = le.xy();
0788 recHit_.yy = le.yy();
0789
0790
0791
0792 GlobalPoint GP = PixGeom->surface().toGlobal(pixeliter->localPosition());
0793 recHit_.gx = GP.x();
0794 recHit_.gy = GP.y();
0795 recHit_.gz = GP.z();
0796 recHit_.subid = subid;
0797 recHit_.layer = layer_num;
0798 recHit_.nsimhit = num_simhit;
0799
0800 if (num_simhit > 0) {
0801 float sim_x1 = (*closest_simhit).entryPoint().x();
0802 float sim_x2 = (*closest_simhit).exitPoint().x();
0803 recHit_.hx = 0.5 * (sim_x1 + sim_x2);
0804 float sim_y1 = (*closest_simhit).entryPoint().y();
0805 float sim_y2 = (*closest_simhit).exitPoint().y();
0806 recHit_.hy = 0.5 * (sim_y1 + sim_y2);
0807
0808 }
0809
0810
0811
0812
0813
0814
0815
0816 }
0817 void StdHitNtuplizer::fillPRecHit(const int subid, trackingRecHit_iterator ih, const GeomDet* PixGeom) {
0818 TrackingRecHit* pixeliter = (*ih)->clone();
0819 LocalPoint lp = pixeliter->localPosition();
0820 LocalError le = pixeliter->localPositionError();
0821
0822 recHit_.x = lp.x();
0823 recHit_.y = lp.y();
0824 recHit_.xx = le.xx();
0825 recHit_.xy = le.xy();
0826 recHit_.yy = le.yy();
0827 GlobalPoint GP = PixGeom->surface().toGlobal(pixeliter->localPosition());
0828 recHit_.gx = GP.x();
0829 recHit_.gy = GP.y();
0830 recHit_.gz = GP.z();
0831 delete pixeliter;
0832 recHit_.subid = subid;
0833 }
0834
0835 void StdHitNtuplizer::fillEvt(const edm::Event& E) {
0836 evt_.run = E.id().run();
0837 evt_.evtnum = E.id().event();
0838 }
0839
0840 void StdHitNtuplizer::init() {
0841 evt_.init();
0842 recHit_.init();
0843 striprecHit_.init();
0844 }
0845
0846 void StdHitNtuplizer::evt::init() {
0847 int dummy_int = 9999;
0848 run = dummy_int;
0849 evtnum = dummy_int;
0850 }
0851
0852 void StdHitNtuplizer::RecHit::init() {
0853 float dummy_float = 9999.0;
0854
0855 x = dummy_float;
0856 y = dummy_float;
0857 xx = dummy_float;
0858 xy = dummy_float;
0859 yy = dummy_float;
0860 row = dummy_float;
0861 col = dummy_float;
0862 gx = dummy_float;
0863 gy = dummy_float;
0864 gz = dummy_float;
0865 layer = 0;
0866 nsimhit = 0;
0867 hx = dummy_float;
0868 hy = dummy_float;
0869 tx = dummy_float;
0870 ty = dummy_float;
0871 theta = dummy_float;
0872 phi = dummy_float;
0873 }
0874
0875
0876 DEFINE_FWK_MODULE(StdHitNtuplizer);