MarlinTrkProcessors  2.4.1
FPCCDSiliconTracking_MarlinTrk.h
1 #ifndef FPCCDSILICONTRACKING_MarlinTrk_H
2 #define FPCCDSILICONTRACKING_MarlinTrk_H 1
3 
4 #include "marlin/Processor.h"
5 #include <marlin/Global.h>
6 #include "lcio.h"
7 #include <string>
8 #include <vector>
9 #include <map>
10 #include <algorithm>
11 #include <cmath>
12 #include <IMPL/TrackImpl.h>
13 #include "ClusterExtended.h"
14 #include "TrackExtended.h"
15 #include "TrackerHitExtended.h"
16 #include "HelixClass_double.h"
17 #include <EVENT/SimTrackerHit.h>
18 #include <EVENT/MCParticle.h>
19 #include <IMPL/TrackImpl.h>
20 
21 #include "MarlinTrk/IMarlinTrack.h"
22 
23 #include <UTIL/BitField64.h>
24 #include <UTIL/ILDConf.h>
25 
26 #include "TTree.h"
27 #include "TNtupleD.h"
28 #include "TFile.h"
29 #include "TVector3.h"
30 #include "TStopwatch.h"
31 #include "./MCPMap.h"
32 #include "../src/GetPurity.h" //by Mori
33 
34 namespace DiagnosticsHistograms {
35  class Histograms ;
36 }
37 
38 namespace DiagnosticsHistograms2D {
39  class Histograms ;
40 }
41 
42 using namespace lcio ;
43 using namespace marlin ;
44 
45 namespace gear{
46  class GearMgr ;
47 }
48 
49 namespace MarlinTrk {
50  class HelixFit;
51  class IMarlinTrkSystem ;
52 }
53 
54 namespace UTIL{
55  class LCRelationNavigator ;
56 }
57 
58 
217 class FPCCDSiliconTracking_MarlinTrk : public Processor {
218 
219 public:
220 
221  virtual Processor* newProcessor() { return new FPCCDSiliconTracking_MarlinTrk ; }
222 
223 
225 
229  virtual void init() ;
230 
233  virtual void processRunHeader( LCRunHeader* run ) ;
234 
237  virtual void processEvent( LCEvent * evt ) ;
238 
239 
240  virtual void check( LCEvent * evt ) ;
241 
242 
245  virtual void end() ;
246 
247 
248 protected:
249 
250  int _nRun ;
251  int _nEvt ;
252  EVENT::LCEvent* _current_event;
253 
254  int _nDivisionsInPhi;
255  int _nDivisionsInTheta;
256  int _nLayers;
257 
258  MarlinTrk::HelixFit* _fastfitter;
259 
262  MarlinTrk::IMarlinTrkSystem* _trksystem ;
263  bool _runMarlinTrkDiagnostics;
264  std::string _MarlinTrkDiagnosticsName;
265 
266  bool _MSOn, _ElossOn, _SmoothOn ;
267 
268  float _initialTrackError_d0;
269  float _initialTrackError_phi0;
270  float _initialTrackError_omega;
271  float _initialTrackError_z0;
272  float _initialTrackError_tanL;
273 
274  double _maxChi2PerHit;
275  double _maxChi2PerHit2nd;
276 
277  bool _UseEventDisplay;
278  std::vector<int> _colours;
279 
280  void drawEvent();
281 
282 
283  // histogram member variables
284 
285  bool _createDiagnosticsHistograms;
287 
288 
289  int _ntriplets, _ntriplets_good, _ntriplets_2MCP, _ntriplets_3MCP, _ntriplets_1MCP_Bad, _ntriplets_bad;
290 
291 
293  LCCollection* GetCollection( LCEvent * evt, std::string colName ) ;
294 
296  LCRelationNavigator* GetRelations( LCEvent * evt, std::string RelName ) ;
297 
299  std::string _colNameMCParticles;
300 
301 
304  // n.b.: a and b should be TrackExtended const *, but the getters are not const :-(
305  bool operator()(TrackExtended *a, TrackExtended *b) const {
306  if ( a == b ) return false;
307  return (a->getChi2()/a->getNDF() < b->getChi2()/b->getNDF() );
308  }
309  };
310 
311 
312  std::string _VTXHitCollection;
313  std::string _FTDPixelHitCollection;
314  std::string _FTDSpacePointCollection;
315  std::string _SITHitCollection;
316  std::string _siTrkCollection;
317 
318  std::vector< LCCollection* > _colTrackerHits;
319  std::map< LCCollection*, std::string > _colNamesTrackerHits;
320 
321  std::vector<TrackerHitExtendedVec> _sectors;
322  std::vector<TrackerHitExtendedVec> _sectorsFTD;
323 
331  public:
333  void clear();
334 
336  inline void resize(size_t maxHits) {
337  _tracksNHits.resize(maxHits-2);
338  _maxIndex=(maxHits-3);
339  }
340 
341  // Sort all track vectors according to chi2/nDof
342  // void sort();
343 
348  inline TrackExtendedVec & getTracksWithNHitsVec( size_t nHits ) {
349  //return _tracksNHits[ std::min(nHits-3, _maxIndex) ];
350  // for debugging: with boundary check
351  return _tracksNHits.at(std::min(nHits-3, _maxIndex));
352  }
353 
354  protected:
355  std::vector< TrackExtendedVec > _tracksNHits;
356  size_t _maxIndex;
357  };
358 
359  TracksWithNHitsContainer _tracksWithNHitsContainer;
360 
361  int InitialiseVTX(LCEvent * evt);
362  int InitialiseFTD(LCEvent * evt);
363  void ProcessOneSector(int iSectorPhi, int iSectorTheta);
364  void ProcessOneSectorVer2(int iSectorPhi, int iSectorTheta);
365  void CleanUp();
366  TrackExtended * TestTriplet(TrackerHitExtended * outerHit,
367  TrackerHitExtended * middleHit,
368  TrackerHitExtended * innerHit,
369  HelixClass_double & helix,
370  int omegamode);
371 
372  int BuildTrack_KalFit(TrackerHitExtended * outerHit,
373  TrackerHitExtended * middleHit,
374  TrackerHitExtended * innerHit,
375  HelixClass_double & helix,
376  int innerlayer,
377  TrackExtended * trackAR);
378 
379  int _useBuildTrackForHighPt;
380  double _cosThetaRangeForBuildTrackForHighPt;
381  double _phiRangeForBuildTrackForHighPt;
382  void getPhiThetaRegionForHighPt(int* boundaries,TrackExtended* trackAR);
383 
384  void Sorting( TrackExtendedVec & trackVec);
385  void CreateTrack(TrackExtended * trackAR );
386  void AttachRemainingVTXHitsSlow();
387  void AttachRemainingFTDHitsSlow();
388  void AttachRemainingVTXHitsFast();
389  void AttachRemainingFTDHitsFast();
390  void AttachRemainingVTXHitsVeryFast();
391  void TrackingInFTD();
392  int BuildTrackFTD(TrackExtended* trackAR, int* nLR, int iS);
393  int AttachHitToTrack(TrackExtended * trackAR, TrackerHitExtended * hit, int iopt);
394  int AttachHitToTrack_KalFit(TrackExtended * trackAR, TrackerHitExtended * hit);
395 
396  void FinalRefit(LCCollectionVec* trk_col, LCCollectionVec* rel_col);
397 
398  float _bField;
399  float _chi2WRPhiTriplet;
400  float _chi2WRPhiQuartet;
401  float _chi2WRPhiSeptet;
402  float _chi2WZTriplet;
403  float _chi2WZQuartet;
404  float _chi2WZSeptet;
405  float _minDistCutAttachForFTD;
406  float _minDistCutAttachForVXD;
407  int _minimalLayerToAttach;
408  int _minMissAddition;
409 
410  // two pi is not a constant in cmath. Calculate it, once!
411  static const double TWOPI;
412 
413  double _dPhi;
414  double _dTheta;
415  double _dPhiFTD;
416 
417 
418 
419  std::vector<int> _Combinations;
420  std::vector<int> _CombinationsFTD;
421 
422  float _resolutionRPhiVTX;
423  float _resolutionZVTX;
424 
425  float _resolutionRPhiFTD;
426  float _resolutionZFTD;
427 
428  float _resolutionRPhiSIT;
429  float _resolutionZSIT;
430 
431  float _phiCutForMerging;
432  float _tanlambdaCutForMerging;
433  float _angleCutForMerging;
434  //float _angleCutForMerging_highPt;
435  //float _angleCutForMerging_lowPt;
436 
437  int _print;
438  int _checkForDelta;
439  float _minDistToDelta;
440 
441  float _distRPhi;
442  float _distZ;
443  float _chi2FitCut;
444  float _chi2FitCut_lowPt;
445 
446 
447  TrackExtendedVec _trackImplVec;
448 
449 
450  float _cutOnD0, _cutOnZ0, _cutOnOmegaVXD, _cutOnOmegaFTD;
451  double _cutOnPtVXD,_cutOnPtFTD;
452 
453  int _minimalHits;
454  int _nHitsChi2;
455  int _attachVXD;
456  int _attachFTD;
457 
458  int _max_hits_per_sector;
459 
460  int _nTotalVTXHits,_nTotalFTDHits,_nTotalSITHits;
461  int _useSIT;
462  int _useFTD;
463 
464 
465  // int _createMap;
466 
467  UTIL::BitField64* _encoder;
468  int getDetectorID(TrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::subdet]; }
469  int getSideID(TrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::side]; };
470  int getLayerID(TrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::layer]; };
471  int getModuleID(TrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::module]; };
472  int getSensorID(TrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::sensor]; };
473 
474  int getDetectorID(SimTrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::subdet]; }
475  int getSideID(SimTrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::side]; };
476  int getLayerID(SimTrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::layer]; };
477  int getModuleID(SimTrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::module]; };
478  int getSensorID(SimTrackerHit* hit) { _encoder->setValue(hit->getCellID0()); return (*_encoder)[lcio::ILDCellID0::sensor]; };
479 
480  void setupGearGeom( const gear::GearMgr* gearMgr ) ;
481 
482 
483  unsigned int _nLayersVTX;
484 
485  unsigned int _nLayersSIT;
486 
487 
488  std::vector<float> _zLayerFTD;
489 
490  unsigned int _nlayersFTD;
491  bool _petalBasedFTDWithOverlaps;
492  int _nPhiFTD;
493 
494  int _output_track_col_quality;
495  static const int _output_track_col_quality_GOOD;
496  static const int _output_track_col_quality_FAIR;
497  static const int _output_track_col_quality_POOR;
498 
499  int _sw_theta;//search window theta
500  float _chi2FitCut_kalman;
501  bool _useClusterRejection;
502  float _minDotOf2Clusters;
503 
504  int getIntersectionEasy(HelixClass_double& helix, TrackerHit* curInmos , int layer, double* isec, double* ref);
505  int getIntersectionEasyTest(HelixClass_double& helix, TrackerHit* basisTrkhit , int layer, std::vector<double> &vec );
506  int CheckTiltOf2Clusters(TrackerHit* A, TrackerHit* B, int level);
507  float DotOf2Clusters(TrackerHit* A, TrackerHit* B);
508  int KalFit(int& ndf, float& Chi2, TrackerHitVec trkHits,TrackerHitVec& hits_in_fit, TrackerHitVec& outliers, float* par , float* epar, HelixClass_double& helix);
509  int getPhiThetaRegion(TrackExtended* trackAR, int layer, int* Boundaries);
510 
511  struct GeoData_t {
512  int nladder;
513  double rmin; // distance of inner surface of sensitive region from IP
514  double dphi; // azimuthal angle step of each ladder
515  double phi0; // aximuthal angle offset
516  std::vector<double> cosphi; // cos[phi_ladder], cos_phi of each ladder
517  std::vector<double> sinphi; // sin[phi_ladder], sin_phi of each ladder
518  std::vector<double> phi; // phi of each ladder
519  std::vector<double> phiAtXiMin; // phiAtXiMin of each ladder
520  std::vector<double> phiAtXiMax; // phiAtXiMax of each ladder
521  std::vector<double> ladder_incline;//the tilt of the line of the ladder expressed by phi
522  double sthick; // sensitive region thickness
523  double sximin; // minimum xi of sensitive region.
524  double sximax; // maximum xi of sensitive region
525  double hlength; // ladder's half length in z
526  int num_xi_pixel; // Number of xi pixel in this ladder
527  int num_zeta_pixel; // Number of zeta pixel in this ladder
528  double rmes; //distance in R of measurement surface
529  };
530 
531  struct vxdGeoData_t {
532  int nLayer;
533  int maxLadder;
534  std::vector<GeoData_t> geodata;
535  }_vxd,_sit;
536 
537  void InitVXDGeometry();
538  void InitSITGeometry();
539  FloatVec _pixelSizeVec;
540  float _pixelheight;
541  TVector3 LocalToGlobal(TVector3 local,int layer,int ladder);
542  //purityMCP CheckOriginOf2Clusters(TrackerHit* A, TrackerHit* B);
543  void calcTrackParameterOfMCP(MCParticle* pmcp, double* par);
544 
546  public :
547  int cellid0 ;
548  int cellid1;
549  unsigned int layer ;
550  unsigned int ladder;
551  unsigned int xiwidth ;
552  unsigned int zetawidth ;
553  unsigned int nPix ;
554  unsigned int tilt ;
555  unsigned int quality ;
556 
557  ClusterStatus(TrackerHit* hit){
558  cellid0 = hit->getCellID0();
559  cellid1 = hit->getCellID1();
560  //layer = ( cellid0 >> 5 ) & 0x000001ff;//maybe bug.
561  //ladder = (cellid0 >> 15 ) & 0x000000ff;
562  xiwidth = cellid1 & 0x000001ff ;
563  zetawidth = ( cellid1 >> 9 ) & 0x000001ff;
564  nPix = ( cellid1 >> 18) & 0x000003ff;
565  tilt = ( cellid1 >> 28) & 0x00000003;
566  quality = ( cellid1 >> 30) & 0x00000003;
567  const int cellId_trk = hit->getCellID0();
568  UTIL::BitField64 encoder_trk( lcio::ILDCellID0::encoder_string ) ;
569  encoder_trk.setValue(cellId_trk) ;
570  layer = encoder_trk[lcio::ILDCellID0::layer] ;
571  ladder = encoder_trk[lcio::ILDCellID0::module] ;
572  }
573  ClusterStatus(){;}
574  };
575 
576  //Old Ver//////////////////
577  typedef std::map< std::pair< int, int >, int > RangeMap;
578  RangeMap _phiRangeForTriplet;
579 
580  double getNeededPhiSectors(double Pt, int outly , int inly); //Difference of radian is returned
582 
583  //New Ver////////////////////==under construction====
584  typedef std::map< std::vector<int> , std::vector<int> > RangeMapVer2;
585  RangeMapVer2 _phiRangeForTripletVer2;
586  void getNeededPhiSectorsVer2(double Pt, std::vector<int> layers, std::vector<double>& phiDiff);
588 
589  float _safetyPhiRange_ratio;
595  float _fudgeFactorForSITsr_rphi;
596  float _fudgeFactorForSITsr_z;
597  int _fudgePhiRange;
598  int _fudgeThetaRange;
599 
600  std::string _colNameVXDTrackerHitRelations;
601  std::string _colNameSITSpacePointRelations;
602  std::string _colNameFTDSpacePointRelations;
603  std::string _colNameFTDPixelTrackerHitRelations;
604  LCRelationNavigator* _navVXD;
605  LCRelationNavigator* _navSIT;
606  LCRelationNavigator* _navFTDsp;
607  LCRelationNavigator* _navFTDpix;
608  std::string _colNameVXDSimHit;
609  std::string _colNameSITSimHit;
610  std::string _colNameFTDspSimHit;
611  std::string _colNameFTDpixSimHit;
612  LCCollection* _simVXD;
613  LCCollection* _simSIT;
614  LCCollection* _simFTDsp;
615  LCCollection* _simFTDpix;
616 
617 
618 
619  double _nSigmaBuild_phi;
620  double _nSigmaBuild_theta;
621 
622  bool _keepCandidate;//used in AttachRemainingVTXHitsVeryFast <-- under construction for now
623 
624  moriUTIL* _moriUtil;
625  GetPurityUtil* _purityUtil;
626 
630  bool _mydebug;
631  bool _mydebugKalFit;
632  bool _mydebugIntersection;
633  bool _mydebugRoot;
634  bool _mydebugRootVXD;
635  bool _mydebugRootFTD;
636  bool _mydebugRootVXDFTD;
637  bool _mydebugVXDHits;
638  bool _mydebugSITHits;
639  bool _mydebugTriplet;
640  int _mydebugTripletMode;
641  bool _mydebugTripletVXD;
642  bool _mydebugTripletFTD;
643  bool _mydebugTripletVXDFTD;
644  bool _mydebugBuildTrack;
645  int _mydebugBuildTrackMode;
646  bool _mydebugAttachVXD;
647  bool _mydebugPrintMCP;
648  bool _stopwatch;
649  class Timer{
650  public:
651  TStopwatch inAnEvt;
652  TStopwatch InitialiseVTX;
653  TStopwatch InitialiseFTD;
654  TStopwatch ProcessOneSector;
655  TStopwatch TrackingInFTD;
656  TStopwatch Sorting;
657  TStopwatch CreateTrack;
658  TStopwatch AttachRemainingVXD;
659  TStopwatch AttachRemainingFTD;
660  TStopwatch FinalRefit;
661  void reset(){
662  inAnEvt.Reset();
663  InitialiseVTX.Reset();
664  InitialiseFTD.Reset();
665  ProcessOneSector.Reset();
666  TrackingInFTD.Reset();
667  Sorting.Reset();
668  CreateTrack.Reset();
669  AttachRemainingVXD.Reset();
670  AttachRemainingFTD.Reset();
671  FinalRefit.Reset();
672  }
673  void cout(){
674  printf("inAnEvt : RT=%.3f s, CPU=%.3f s \n",inAnEvt.RealTime(),inAnEvt.CpuTime());
675  printf("InitialiseVTX : RT=%.3f s, CPU=%.3f s \n",InitialiseVTX.RealTime(),InitialiseVTX.CpuTime());
676  printf("InitialiseFTD : RT=%.3f s, CPU=%.3f s \n",InitialiseFTD.RealTime(),InitialiseFTD.CpuTime());
677  printf("ProcessOneSector : RT=%.3f s, CPU=%.3f s \n",ProcessOneSector.RealTime(),ProcessOneSector.CpuTime());
678  printf("TrackingInFTD : RT=%.3f s, CPU=%.3f s \n",TrackingInFTD.RealTime(),TrackingInFTD.CpuTime());
679  printf("Sorting : RT=%.3f s, CPU=%.3f s \n",Sorting.RealTime(),Sorting.CpuTime());
680  printf("CreateTrack : RT=%.3f s, CPU=%.3f s \n",CreateTrack.RealTime(),CreateTrack.CpuTime());
681  printf("AttachRemainingVXD: RT=%.3f s, CPU=%.3f s \n",AttachRemainingVXD.RealTime(),AttachRemainingVXD.CpuTime());
682  printf("AttachRemainingFTD: RT=%.3f s, CPU=%.3f s \n",AttachRemainingFTD.RealTime(),AttachRemainingFTD.CpuTime());
683  printf("FinalRefit : RT=%.3f s, CPU=%.3f s \n",FinalRefit.RealTime(),FinalRefit.CpuTime());
684  }
685  Timer(){reset();}
686  }_timer;
687  TStopwatch _timer2Triplet;
688  TStopwatch _timer2Build;
689  bool _mydebugstopwatch2;
690 
691  float _currentPurity;
692  MCPMap LoadMCPMap();
693  std::map< MCParticle*, SimTrackerHitVec > _mcpVXD;
694  std::map< MCParticle*, SimTrackerHitVec > _mcpVXDFTD;
695  std::map< MCParticle*, SimTrackerHitVec > _mcpVXDSIT;
696  std::map< MCParticle*, SimTrackerHitVec > _mcpVXDFTDSIT;
697  std::map< MCParticle*, SimTrackerHitVec > _mcpFTD;
698  std::map< MCParticle*, SimTrackerHitVec > _mcpFTDSIT;
699  std::map< MCParticle*, SimTrackerHitVec > _mcpSIT;
700 
701 
702  enum MCPContributions{
703  contVXD,
704  contFTD,
705  contSIT,
706  contVXDFTD,
707  contVXDSIT,
708  contFTDSIT,
709  contVXDFTDSIT,
710  contBG,
711  contSize
712  };
713 
714  enum MCPContributions getMCPContribution(IntVec nsub);
715 
716  class Triplet : public IMPL::TrackImpl{
717  public :
718  static int numOfProcesses;
719  int id;
720  int quality_code;
721  int detID[3];//内側に近いところから詰める。
722  int layer[3];
723  float probability;
724  float purity;
725  SimTrackerHitVec simVecFromTrk;//TrackerHitVecに対応するもの
726  std::pair< MCParticle* , SimTrackerHitVec > mcTruth ;//Dominant MCPについての情報
727 
728  //new (2013_08_10)
729  TrackerHitVec truthTrkHits ;//Dominant MCPのsimthitsに対応するTrackerHitVec.
730  enum MCPContributions mcpcont;
731  IntVec nsub;//vxd:0,sit:1,ftd:2 <--for Dominant MCP
732 
733 
734  float mcp_d0;
735  float mcp_z0;
736  float mcp_omega;
737  float mcp_phi0;
738  float mcp_tanL;
739 
740  //bool availableInBuildedTrack;
741 
742  Triplet(bool incrementID = true){
743  id = numOfProcesses;
744  if(incrementID){
745  numOfProcesses++;
746  }
747  quality_code = -1;
748  detID[0] = detID[1] = detID[2] = -1;
749  layer[0] = layer[1] = layer[2] = -1;
750  probability = 1e20;
751  purity = 1e20;
752  simVecFromTrk.clear();
753  mcTruth.first = 0; mcTruth.second.clear();
754  truthTrkHits.clear();
755  mcpcont = contSize;
756  nsub.clear();
757  mcp_d0 = mcp_z0 = mcp_omega = mcp_phi0 = mcp_tanL = 1e20;
758  }
759  Triplet(const Triplet& obj){
760  *this = obj;
761  }
762  };
763  Triplet* _curtriplet;
764 
765 
766  enum BuildTrackResult{
767  NormalEnd,
768  ManyMisAssignments,
769  ManyOutliers,
770  PhiThetaError,
771  BuildTrackResultSize
772  };
773  class BuildedTrack : public IMPL::TrackImpl{
774  public :
775  static int numOfProcesses;
776  int id;
777 
778  BuildTrackResult result;
779  int nMisAssign;
780  Triplet triplet;
781  const MCParticle* truthmcp;
782  SimTrackerHitVec simvec;
783  float probability;
784  float purity;
785  SimTrackerHitVec simVecFromTrk;//TrackerHitVecに対応するもの
786  BuildedTrack():triplet(false){
787  id = numOfProcesses;
788  numOfProcesses++;
789  result = BuildTrackResultSize;
790  nMisAssign = -1;
791  truthmcp = 0;
792  simvec.clear();
793  probability = 1e20;
794  purity = 1e20;
795  simVecFromTrk.clear();
796  }
797  BuildedTrack(const BuildedTrack& obj){
798  *this = obj;
799  }
800  };
801  bool _availableInBuildedTrack;
802 
803  typedef std::vector< BuildedTrack > BuildedTrackVec;
804  BuildedTrackVec _buildedTrackContainer;
805  void SetBuildedTrack(TrackExtended* trackAR, BuildedTrackVec& btrackvec);
806  void BuildedTrackDebuger1(BuildedTrackVec::iterator begin, BuildedTrackVec::iterator end);
807  void BuildedTrackDebuger2(std::vector<const BuildedTrack*>::iterator begin, std::vector<const BuildedTrack*>::iterator end);
808 
809  //時間がアレばclass templateにしてみたい。
811  public :
812  std::pair<MCParticle*,SimTrackerHitVec> pair;
813  IntVec nsub;//nhits in sub detectors. 0 VXD, 1 FTD, 2 SIT
814  std::vector<const BuildedTrack*> BuildedTrackVec;
815  std::vector<const BuildedTrack*> BuildedTrack2ndVec;
816  std::vector<const BuildedTrack*> BuildedTrack3rdVec;
817  };
818  typedef std::map<MCParticle*,MCP_BuildedTrack> BTMap;
819  BTMap _mcpBuildedTrackContainer;
820 
821  std::vector< std::map<MCParticle*, MCP_BuildedTrack> > _mcpRemainingBTContainerA;
822  //purity 100%のトリプレットが生成されなかった生成されなかったmcp
823  std::vector< std::map<MCParticle*, MCP_BuildedTrack> > _mcpRemainingBTContainerB;
824  //purity 100%のトリプレットが生成されたものの、quality_code != 0により結局再構成されなかったmcp
825 
826 
827 
828 
829  std::vector< Triplet > _tripletContainer;
830 
831  class MCP_Triplet{
832  public :
833  std::pair<MCParticle*,SimTrackerHitVec> pair;
834  IntVec nsub;//nhits in sub detectors. 0 VXD, 1 FTD, 2 SIT
835  std::vector<const Triplet*> tripvec;
836  std::vector<const Triplet*> trip2ndvec;
837  std::vector<const Triplet*> trip3rdvec;
838  };
839  typedef std::map<MCParticle*,MCP_Triplet> TripMap;
840  TripMap _mcpTripletContainer;
841 
842 
843  std::vector< std::map<MCParticle*, MCP_Triplet> > _mcpRemainingTRContainerA;
844  //purity 100%のトリプレットが生成されなかった生成されなかったmcp
845  std::vector< std::map<MCParticle*, MCP_Triplet> > _mcpRemainingTRContainerB;
846  //purity 100%のトリプレットが生成されたものの、quality_code != 0により結局再構成されなかったmcp
847 
848  void TripletDebugerWithMCPRemain(int nth, std::vector< std::map<MCParticle*, MCP_Triplet> > A);
849 
850  MCPMap _mcpMap;
851  std::vector<LCRelationNavigator*> _naviVec;
852  IntVec getNHitsInSubDet(SimTrackerHitVec simvec);
853  void TripletDebuger1(std::vector<Triplet>::iterator begin, std::vector<Triplet>::iterator end);
854  void TripletDebuger2(std::vector<const Triplet*>::iterator begin, std::vector<const Triplet*>::iterator end);
855 
856 
857 };
858 
859 int FPCCDSiliconTracking_MarlinTrk::Triplet::numOfProcesses = 0;
860 int FPCCDSiliconTracking_MarlinTrk::BuildedTrack::numOfProcesses = 0;
861 
862 #endif
Compare tracks according to their chi2/ndf.
Definition: FPCCDSiliconTracking_MarlinTrk.h:303
void resize(size_t maxHits)
Set the size to allow a maximum of maxHit hits.
Definition: FPCCDSiliconTracking_MarlinTrk.h:336
TrackExtendedVec & getTracksWithNHitsVec(size_t nHits)
Returns the TrackExtendedVec for track with n hits.
Definition: FPCCDSiliconTracking_MarlinTrk.h:348
Definition: FPCCDSiliconTracking_MarlinTrk.h:649
Definition: FPCCDSiliconTracking_MarlinTrk.h:531
Definition: FPCCDSiliconTracking_MarlinTrk.h:773
Definition: FPCCDSiliconTracking_MarlinTrk.h:810
Definition: FPCCDSiliconTracking_MarlinTrk.h:545
MarlinTrk::IMarlinTrkSystem * _trksystem
pointer to the IMarlinTrkSystem instance
Definition: FPCCDSiliconTracking_MarlinTrk.h:262
A helper class to allow good code readability by accessing tracks with N hits.
Definition: FPCCDSiliconTracking_MarlinTrk.h:330
Definition: GetPurity.h:37
=== FPCCDSiliconTracking_MarlinTrk Processor === This processor is based on SiliconTracking_MarlinT...
Definition: FPCCDSiliconTracking_MarlinTrk.h:217
std::string _colNameMCParticles
input MCParticle collection and threshold used for Drawing
Definition: FPCCDSiliconTracking_MarlinTrk.h:299
Definition: FPCCDSiliconTracking_MarlinTrk.h:511
Definition: SiliconTracking_MarlinTrk.cc:90
Definition: FPCCDSiliconTracking_MarlinTrk.h:831
Definition: FPCCDSiliconTracking_MarlinTrk.h:716
Definition: MCPMap.h:20
int _safetyPhiRange_fix
Extra range in addition to main range used for triplet construction process and needed to find triple...
Definition: FPCCDSiliconTracking_MarlinTrk.h:594