InttSeedTrackPerformance-Jingyu/src/InttSeedTrackPerformance.cxx

0001 #ifndef InttSeedTrackPerformance_cxx
0002 #define InttSeedTrackPerformance_cxx
0003
0004 #include "InttSeedTracking.h"
0005 #include "InttSeedTrackPerformance.h"
0006
0007
0008 #include <TStyle.h>
0009 #include <TCanvas.h>
0010
0011 #include "SPHTracKuma.h"
0012
0013 bool bCaloClu = false;
0014 Int_t strockEvents = 1;
0015
0016 void InttSeedTrackPerformance::Loop(Int_t runNum)
0017 {
0018    if (fChain == 0) return;
0019    HistInit();
0020
0021    Long64_t nentries = fChain->GetEntriesFast();
0022    Long64_t nbytes = 0, nb = 0;
0023
0024    // Using the EventJudge function, you can get the event IDs finally.
0025    // ex) If you want to searching for the event -2 < dpt < -1,
0026    // you need to set EventJudge(dpt, -2, -1, true);
0027    // The target event IDs are straged in the vector m_vTargetEvents, and finally they are shown.
0028    // You can run only target events using bTargetEV.
0029    // (1) Turn on the bTargetEV to truth.
0030    // (2) Set the events in vTargetEvents.
0031    // (3) Change the loop the second one ("jentry< vTargetEvents.size()")
0032    bool bTargetEV = false;
0033    std::vector<Int_t > vTargetEvents = {};
0034    // 9, 13, 17, 19, 29, 53, 73, 84, 92, 367
0035    // 35, 102, 108, 192, 213, 272, 283, 315, 324, 327, 341, 350
0036
0037    // == s0 == Event Loop 000000000000000000000000000000000000000000000000000000
0038    // for (Long64_t jentry=0; jentry<100;jentry++) {
0039    // for (Long64_t jentry=0; jentry< vTargetEvents.size();jentry++) {
0040    for (Long64_t jentry=(1000*runNum); jentry<(1000*runNum + 9900); jentry++) {
0041    // for (Long64_t jentry=0; jentry<nentries;jentry++) {
0042
0043       Long64_t tempJEntry = jentry;
0044       if(bTargetEV) jentry =  vTargetEvents.at(jentry);
0045       pubEvNum = jentry;
0046
0047       Long64_t ientry = LoadTree(jentry);
0048       if (ientry < 0) break;
0049       nb = fChain->GetEntry(jentry);   nbytes += nb;
0050
0051       if(jentry%100 == 0) std::cout << " == event" << jentry << " =====" << std::endl;
0052
0053       CheckPrimP(m_TruthParticle);
0054       // std::cout << "truth size = " << m_TruthParticle.size() << std::endl;
0055
0056       // MVTX f: first, s: second, t: third
0057       ReadInttHitting(m_fMvtxHits, m_sMvtxHits, m_tMvtxHits, m_iInttHits, m_oInttHits, m_TpcHits);
0058       if(bCaloClu) ReadCalCluHitting(m_emcalHits, m_iHCalHits, m_oHCalHits);
0059       else ReadCalHitting(m_emcalHits, m_iHCalHits, m_oHCalHits);
0060
0061       // std::cout << "== event No." << jentry << "  =========================== "<< std::endl;
0062       // ShowEventInfo();
0063
0064       m_vTruthParticle.push_back(m_TruthParticle);
0065       m_viInttHits.push_back(m_iInttHits);
0066       m_voInttHits.push_back(m_oInttHits);
0067       m_vTpcHits.push_back(m_TpcHits);
0068       m_vemcalHits.push_back(m_emcalHits);
0069       m_viHCalHits.push_back(m_iHCalHits);
0070       m_voHCalHits.push_back(m_oHCalHits);
0071       if(jentry%strockEvents == 0){
0072          for(Int_t iEvent = 0; iEvent < strockEvents; iEvent++){
0073             m_iInttHits = m_viInttHits.at(iEvent);
0074             m_oInttHits = m_voInttHits.at(iEvent);
0075             m_TpcHits = m_vTpcHits.at(iEvent);
0076             m_emcalHits = m_vemcalHits.at(iEvent);
0077             InttSeedTracking *TracKumaContainer = new InttSeedTracking(m_tracks, m_fMvtxHits, m_sMvtxHits, m_tMvtxHits, m_viInttHits.at(iEvent), m_voInttHits.at(iEvent),\
0078                m_vemcalHits.at(iEvent), m_viHCalHits.at(iEvent), m_voHCalHits.at(iEvent));
0079             CaloEnergyQA(0, m_vemcalHits.at(iEvent));
0080             CaloEnergyQA(1, m_viHCalHits.at(iEvent));
0081             CaloEnergyQA(2, m_voHCalHits.at(iEvent));
0082             TrackQA(m_vTruthParticle.at(iEvent), m_vemcalHits.at(iEvent));
0083             m_bDecayEvent = false;
0084          }
0085          AllResetValuesForEvent();
0086       }
0087
0088       PartResetValuesForEvent();
0089
0090       if(bTargetEV) jentry = tempJEntry;
0091    }
0092    // == e0 == Event Loop 000000000000000000000000000000000000000000000000000000
0093
0094    WrightHists();
0095
0096    std::cout << "std::vector<Int_t > vTargetEvents = {";
0097    for(Int_t i = 0; i < m_vTargetEvents.size(); i++) std::cout << m_vTargetEvents.at(i) << ", ";
0098    std::cout << "};" << std::endl;
0099 }
0100
0101 void InttSeedTrackPerformance::HistInit(){
0102    m_HINTTHitMap = new TH2D( "m_HINTTHitMap", "Global Hit Map;x [cm];y [cm]",
0103                 2000, -100., 100., 2000, -100., 100. );
0104
0105    m_HGlobalHitMap = new TH2D( "m_HGlobalHitMap", "Global Hit Map;x [cm];y [cm]",
0106                 5000, -250., 250., 2500, 0., 250. );
0107
0108
0109    m_HMatchCalEVsHighestCalE = new TH1D( "m_HMatchCalEVsHighestCalE", \
0110       "EMCal Energy disparsion map;phi [rad];E [MeV]", \
0111       320, -1.6, 1.6);
0112
0113    m_HdPhiInttdPhiECal = new TH2D( "m_HdPhiInttdPhiECal", \
0114       "EMCal Energy disparsion map;d#it{#phi}_{INTT} [rad];d#it{#phi}_{ECal} [rad]", \
0115          320, -1.6, 1.6, 320, -1.6, 1.6);
0116    m_HTruthPtVsdPhiInttdPhiECal = new TH2D( "m_HTruthPtVsdPhiInttdPhiECal", \
0117       "truth pT vs dPhi (INTT to ECal);#it{p}_{T, truth} [cm];dPhi [rad]", \
0118          150, 0., 15.,  6, -3., 3. );
0119
0120    m_HTruthPtVsdPhiIIntt = new TH2D( "m_HTruthPtVsdPhiIIntt", \
0121       "truth pT vs dPhi (vertex to iINTT);#it{p}_{T, truth} [cm];dPhi [rad]",
0122         150, 0., 15., 6, -3., 3. );
0123    m_HTruthPtVsdPhiOIntt = new TH2D( "m_HTruthPtVsdPhiOIntt", \
0124       "truth pT vs dPhi (vertex to oINTT);#it{p}_{T, truth} [cm];dPhi [rad]",
0125         150, 0., 15.,  6, -3., 3. );
0126
0127
0128    // === pT resolution histograms  ==========
0129    m_HTruthPtVsSagittaPt_rough = new TH3D( "m_HTruthPtVsSagittaPt_rough", \
0130       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0131         150, 0., 15., 250, -2., 2., 20, -1., 1);
0132    m_HTruthPtVsSagittaPt_VtxIntt = new TH3D( "m_HTruthPtVsSagittaPt_VtxIntt", \
0133       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0134         150, 0., 15., 250, -2., 2., 20, -1., 1);
0135    m_HTruthPtVsSagittaPt_InttEmcal = new TH3D( "m_HTruthPtVsSagittaPt_InttEmcal", \
0136       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0137         150, 0., 15., 250, -2., 2., 20, -1., 1);
0138    m_HTruthPtVsSagittaPt_VtxInttEmcal = new TH3D( "m_HTruthPtVsSagittaPt_VtxInttEmcal", \
0139       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0140         150, 0., 15., 250, -2., 2., 20, -1., 1);
0141    m_HTruthPtVsSagittaPt_MvtxInttEmcal = new TH3D( "m_HTruthPtVsSagittaPt_MvtxInttEmcal", \
0142       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0143         150, 0., 15., 250, -2., 2., 20, -1., 1);
0144    m_HTruthPtVsSagittaPt_VtxMvtxInttEmcal = new TH3D( "m_HTruthPtVsSagittaPt_VtxMvtxInttEmcal", \
0145       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0146         150, 0., 15., 250, -2., 2., 20, -1., 1);
0147
0148    m_HTruthPtVsFitFuncPt_IInttOInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_IInttOInttEmcal", \
0149       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0150         150, 0., 15., 250, -2., 2., 20, -1., 1);
0151
0152    m_HTruthPtVsFitFuncPt_VtxIInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_VtxIInttEmcal", \
0153       "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0154         150, 0., 15., 250, -2., 2., 20, -1., 1);
0155    m_HTruthPtVsFitFuncPt_VtxOInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_VtxOInttEmcal", \
0156    "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0157     150, 0., 15., 250, -2., 2., 20, -1., 1);
0158    m_HTruthPtVsFitFuncPt_VtxInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_VtxInttEmcal", \
0159    "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0160     150, 0., 15., 250, -2., 2., 20, -1., 1);
0161
0162    m_HTruthPtVsFitFuncPt_MvtxIInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_MVtxIInttEmcal", \
0163    "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0164     150, 0., 15., 250, -2., 2., 20, -1., 1);
0165    m_HTruthPtVsFitFuncPt_MvtxOInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_MVtxOInttEmcal", \
0166    "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0167     150, 0., 15., 250, -2., 2., 20, -1., 1);
0168    m_HTruthPtVsFitFuncPt_MvtxInttEmcal = new TH3D( "m_HTruthPtVsFitFuncPt_MVtxInttEmcal", \
0169    "truth pT vs delta pt;#it{p}_{T} [GeV/#it{c}];#it{#delta}#it{p}_{T}/#it{p}_{T};#eta",
0170     150, 0., 15., 250, -2., 2., 20, -1., 1);
0171
0172
0173    // === theta resolution histograms  ==========
0174    m_HTruthPVsTheta_InttEmcal = new TH2D( "m_HTruthPVsTheta_InttEmcal", \
0175       "truth p vs delta pt;#it{p} [GeV/#it{c}];#it{#delta}#it{#theta}",
0176         150, 0., 15., 100, -0.5, 0.5 );
0177    m_HTruthPVsTheta_MvtxInttEmcal = new TH2D( "m_HTruthPVsTheta_MvtxInttEmcal", \
0178       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{#delta}#it{#theta}",
0179         150, 0., 15., 100, -0.5, 0.5 );
0180
0181    // === phi resolution histograms  ==========
0182    m_HTruthPVsPhi_InttEmcal = new TH2D( "m_HTruthPVsPhi_InttEmcal", \
0183       "truth p vs delta pt;#it{p} [GeV/#it{c}];#it{#delta}#it{#phi}",
0184         150, 0., 15., 100, -0.5, 0.5 );
0185    m_HTruthPVsPhi_MvtxInttEmcal = new TH2D( "m_HTruthPVsPhi_MvtxInttEmcal", \
0186       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{#delta}#it{#phi}",
0187         150, 0., 15., 100, -0.5, 0.5 );
0188
0189    // === dp histograms  ==========
0190    m_HTruthPVsRecoP_InttEmcal = new TH2D( "m_HTruthPVsRecoP_InttEmcal", \
0191       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{#delta}#it{p}/#it{p}",
0192         150, 0., 15., 50, -2.5, 2.5 );
0193    m_HTruthPVsRecoP_MvtxInttEmcal = new TH2D( "m_HTruthPVsRecoP_MvtxInttEmcal", \
0194       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{#delta}#it{p}/#it{p}",
0195         150, 0., 15., 50, -2.5, 2.5 );
0196    m_HTruthPVsRecoP_FitFunc = new TH2D( "m_HTruthPVsRecoP_FitFunc", \
0197       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{#delta}#it{p}/#it{p}",
0198         150, 0., 15., 50, -2.5, 2.5 );
0199
0200    // === dE histograms  ==========
0201    m_HDE = new TH3D( "m_HDE", "truth E vs Reco E;#it{E} [MeV];#it{#delta}#it{E}/#it{E};#eta",
0202     200, 0., 20., 50, -2.5, 2.5, 20, -1., 1);
0203
0204    // === E/p histograms  ==========
0205    m_HTruthPVsEOverP_InttEmcal = new TH2D( "m_HTruthPVsEOverP_InttEmcal", \
0206       "truth p vs delta pt;#it{p} [GeV/#it{c}];#it{E}/#it{p}",
0207         150, 0., 15., 60, 0., 3. );
0208    m_HTruthPVsEOverP_MvtxInttEmcal = new TH2D( "m_HTruthPVsEOverP_MvtxInttEmcal", \
0209       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{E}/#it{p}",
0210         150, 0., 15., 60, 0., 3. );
0211    m_HTruthPVsEOverP_FitFunc = new TH2D( "m_HTruthPVsEOverP_FitFunc", \
0212       "truth p vs delta p;#it{p} [GeV/#it{c}];#it{E}/#it{p}",
0213         150, 0., 15., 60, 0., 3. );
0214
0215    // === dVtx histograms  ==========
0216    m_dVtxXY_InttEmcal = new TH2D( "m_dVtxXY_InttEmcal", "vertexReso;d#it{X}_{(reco - truth)} [cm];d#it{Y}_{(reco - truth)} [cm];count", 1000, -50., 50.,  1000, -50., 50.);
0217    m_dVtxXY_MvtxInttEmcal = new TH2D( "m_dVtxXY_MvtxInttEmcal", "vertexReso;d#it{X}_{(reco - truth)} [cm];d#it{Y}_{(reco - truth)} [cm];count", 1000, -50., 50., 1000, -50., 50.);
0218
0219    m_dVtxR_InttEmcal = new TH1D( "m_dVtxR_InttEmcal", "vertexReso;d#it{R}_{(reco - truth)} [cm];count", 500, 0., 50.);
0220    m_dVtxR_MvtxInttEmcal = new TH1D( "m_dVtxR_MvtxInttEmcal", "vertexReso;d#it{R}_{(reco - truth)} [cm];count", 500, 0., 50.);
0221
0222    m_dVtxZ_InttEmcal = new TH1D( "m_dVtxZ_InttEmcal", "vertexReso;d#it{Z}_{(reco - truth)} [cm];count", 1000, -50., 50.);
0223    m_dVtxZ_MvtxInttEmcal = new TH1D( "m_dVtxZ_MvtxInttEmcal", "vertexReso;d#it{Z}_{(reco - truth)} [cm];count", 1000, -50., 50.);
0224
0225    // == s == Calo QA Histograms ==================================
0226    m_HECalPhiVsE  = new TProfile2D( "m_HECalPhiVsE",  "EMCal Energy disparsion map;#phi [rad];#eta;E [MeV]", 320, -1.6, 1.6, 300, -1.5, 1.5);
0227    m_HIHCalPhiVsE = new TProfile2D( "m_HIHCalPhiVsE", "iHCal Energy disparsion map;#phi [rad];#eta;E [MeV]", 320, -1.6, 1.6, 300, -1.5, 1.5);
0228    m_HOHCalPhiVsE = new TProfile2D( "m_HOHCalPhiVsE", "oHCal Energy disparsion map;#phi [rad];#eta;E [MeV]", 320, -1.6, 1.6, 300, -1.5, 1.5);
0229    // == e == Calo QA  Histograms ==================================
0230
0231    // == s == Magnetic shift Histograms ==================================
0232    m_HDPhiVsDStraightVsTruPt = new TH3D("m_HDPhiVsDStraightVsTruPt", "dPhi/dl between vertex and inner INTT;d#it{l} [cm];d#it{phi} [rad];#it{p}_{T} [GeV/#it{c}];", 150, 0., 150., 100, -0.05, 0.05, 150, 0., 15.);
0233
0234    m_HDPhiDStraight_VtxIIntt = new TH2D("m_HDPhiDStraight_VtxIIntt", "dPhi/dl between vertex and inner INTT;#it{p}_{T} [GeV/#it{c}];d#it{phi}/d#it{l} [rad/cm]", 150, 0., 15., 100, -0.05, 0.05);
0235    m_HDPhiDStraight1D_VtxIIntt = new TH1D("m_HDPhiDStraight1D_VtxIIntt", "dPhi/dl between vertex and inner INTT;d#it{phi}/d#it{l} [rad/cm];count", 100, -0.05, 0.05);
0236
0237    m_HDPhiDStraight_IInttOIntt = new TH2D("m_HDPhiDStraight_IInttOIntt", "dPhi/dl between vertex and inner INTT;#it{p}_{T} [GeV/#it{c}];d#it{phi}/d#it{l} [rad/cm]", 150, 0., 15., 100, -0.05, 0.05);
0238    m_HDPhiDStraight1D_IInttOIntt = new TH1D("m_HDPhiDStraight1D_IInttOIntt", "dPhi/dl between vertex and inner INTT;d#it{phi}/d#it{l} [rad/cm];count", 100, -0.05, 0.05);
0239
0240    m_HTruthPtVsDdPhiddStraight_IInttOInttEmcal = new TH2D("m_HTruthPtVsDdPhiddStraight_IInttOInttEmcal", "dPhi/dl between vertex and inner INTT;#it{p}_{T} [GeV/#it{c}];d#it{phi}/d#it{l} [rad/cm]", 150, 0., 15., 100, -0.05, 0.05);
0241
0242    m_HDPhiDStraight_OInttEmcal = new TH2D("m_HDPhiDStraight_OInttEmcal", "dPhi/dl between vertex and inner INTT;#it{p}_{T} [GeV/#it{c}];d#it{phi}/d#it{l} [rad/cm]", 150, 0., 15., 100, -0.05, 0.05);
0243    m_HDPhiDStraight1D_OInttEmcal = new TH1D("m_HDPhiDStraight1D_OInttEmcal", "dPhi/dl between vertex and inner INTT;d#it{phi}/d#it{l} [rad/cm];count", 100, -0.05, 0.05);
0244
0245
0246    m_HDLVsTruthPtVsEta_VtxIIntt = new TH3D( "m_HDLVsTruthPtVsEta_VtxIIntt", \
0247       "truth pT vs shift distanse (L);L [cm];#it{p}_{T} [GeV/#it{c}];#eta",
0248         1200, -10., 10., 1000, -50., 50., 20, -1., 1);
0249    m_HDL1D_VtxIIntt = new TH1D( "m_HDL1D_VtxIIntt", \
0250       "truth pT vs shift distanse (L);L [cm]; count", 200, -10., 10.);
0251    m_HDLVsTruthPtVsEta_IInttOIntt = new TH3D( "m_HDLVsTruthPtVsEta_IInttOIntt", \
0252       "truth pT vs shift distanse (L);L [cm];#it{p}_{T} [GeV/#it{c}];#eta",
0253         200, -10., 10., 1000, -50., 50., 20, -1., 1);
0254    m_HDL1D_IInttOIntt = new TH1D( "m_HDL1D_IInttOIntt", \
0255       "truth pT vs shift distanse (L);L [cm]; count", 200, -10., 10.);
0256    m_HDLVsTruthPtVsEta_OInttEmcal = new TH3D( "m_HDLVsTruthPtVsEta_OInttEmcal", \
0257       "truth pT vs shift distanse (L);L [cm];#it{p}_{T} [GeV/#it{c}];#eta",
0258         1000, -50., 50., 1000, -50., 50., 20, -1., 1);
0259    m_HDL1D_OInttEmcal = new TH1D( "m_HDL1D_OInttEmcal", \
0260       "truth pT vs shift distanse (L);L [cm]; count", 1000, -50., 50.);
0261
0262
0263    m_HDPhiVsTruthPtVsEta_VtxIIntt = new TH3D( "m_HDPhiVsTruthPtVsEta_VtxIIntt", \
0264       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0265         200, -1., 1., 150, 0., 15., 10, -1., 1);
0266    m_HDPhiVsTruthPtVsEta_IInttOIntt = new TH3D( "m_HDPhiVsTruthPtVsEta_IInttOIntt", \
0267       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0268         200, -1., 1., 150, 0., 15., 10, -1., 1);
0269    m_HDPhiVsTruthPtVsEta_OInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_OInttEmcal", \
0270       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0271         200, -1., 1., 150, 0., 15., 10, -1., 1);
0272
0273
0274    m_HDPhiVsTruthPtVsEta_IInttOInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_IInttOInttEmcal", \
0275       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0276         200, -1., 1., 150, 0., 15., 10, -1., 1);
0277
0278    m_HDPhiVsTruthPtVsEta_VtxIInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_VtxIInttEmcal", \
0279       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0280         200, -1., 1., 150, 0., 15., 10, -1., 1);
0281    m_HDPhiVsTruthPtVsEta_VtxOInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_VtxOInttEmcal", \
0282       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0283         200, -1., 1., 150, 0., 15., 10, -1., 1);
0284    m_HDPhiVsTruthPtVsEta_VtxInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_VtxInttEmcal", \
0285       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0286         200, -1., 1., 150, 0., 15., 10, -1., 1);
0287
0288    m_HDPhiVsTruthPtVsEta_MvtxIInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_MvtxIInttEmcal", \
0289       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0290         200, -1., 1., 150, 0., 15., 10, -1., 1);
0291    m_HDPhiVsTruthPtVsEta_MvtxOInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_MvtxOInttEmcal", \
0292       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0293         200, -1., 1., 150, 0., 15., 10, -1., 1);
0294    m_HDPhiVsTruthPtVsEta_MvtxInttEmcal = new TH3D( "m_HDPhiVsTruthPtVsEta_MvtxInttEmcal", \
0295       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];#eta",
0296         200, -1., 1., 150, 0., 15., 10, -1., 1);
0297
0298    m_HDPhiVsTruthPtProf_IInttOInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_IInttOInttEmcal", \
0299       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];", 200, -1., 1.);
0300    m_HDPhiVsTruthPtProf_VtxIIntt = new TProfile( "m_HDPhiVsTruthPtProf_VtxIIntt", \
0301       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}];", 200, -1., 1.);
0302    m_HDPhiVsTruthPtProf_IInttOIntt = new TProfile( "m_HDPhiVsTruthPtProf_IInttOIntt", \
0303       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0304    m_HDPhiVsTruthPtProf_OInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_OInttEmcal", \
0305       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0306
0307    m_HDPhiVsTruthPtProf_VtxIInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_VtxIInttEmcal", \
0308       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0309    m_HDPhiVsTruthPtProf_VtxOInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_VtxOInttEmcal", \
0310       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]",200, -1., 1.);
0311    m_HDPhiVsTruthPtProf_VtxInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_VtxInttEmcal", \
0312       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0313
0314    m_HDPhiVsTruthPtProf_MvtxIInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_MvtxIInttEmcal", \
0315       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]",200, -1., 1.);
0316    m_HDPhiVsTruthPtProf_MvtxOInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_MvtxOInttEmcal", \
0317       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0318    m_HDPhiVsTruthPtProf_MvtxInttEmcal = new TProfile( "m_HDPhiVsTruthPtProf_MvtxInttEmcal", \
0319       "truth pT vs shift distanse (L);#phi [rad];#it{p}_{T} [GeV/#it{c}]", 200, -1., 1.);
0320
0321
0322   m_HPtEfficiency = new TH1D( "m_HPtEfficiency", \
0323       "reconstrucntion efficinecy for pT;#it{p}_{T} [GeV/#it{c}];Efficiency [%]", 150, 0., 15);
0324   m_HTruTrackNum = new TH1D( "m_HTruTrackNum", \
0325       "reconstrucntion efficinecy for pT;#it{p}_{T} [GeV/#it{c}];Efficiency [%]", 150, 0., 15);
0326 }
0327
0328
0329 void InttSeedTrackPerformance::CheckPrimP(std::vector<hitStruct >& vTruthParticle){
0330    Int_t numOfP = PrimaryG4P_PID->size();
0331    if(numOfP == 0) return;
0332
0333    for(Int_t iP = 0; iP < numOfP; iP++){
0334       hitStruct truthParticle;
0335       Double_t pPhi = PrimaryG4P_Phi->at(iP);
0336       Double_t pEta = PrimaryG4P_Eta->at(iP);
0337       Double_t pE = PrimaryG4P_E->at(iP);
0338       Double_t pPt = PrimaryG4P_Pt->at(iP);
0339
0340       truthParticle.phi = pPhi;
0341       truthParticle.eta = pEta;
0342       truthParticle.pt = pPt;
0343       truthParticle.energy = pE;
0344
0345       vTruthParticle.push_back(truthParticle);
0346    }
0347 }
0348
0349 void InttSeedTrackPerformance::ReadInttHitting(std::vector<hitStruct >& vFMvtxHits,\
0350     std::vector<hitStruct >& vSMvtxHits, std::vector<hitStruct >& vTMvtxHits,\
0351     std::vector<hitStruct >& vIInttHits, std::vector<hitStruct >& vOInttHits,\
0352     std::vector<hitStruct >& vTpcHits)
0353 {
0354     Int_t numOfClu = trk_system->size();
0355     if(numOfClu == 0) return;
0356
0357     for(Int_t iClu = 0; iClu < numOfClu; iClu++)
0358     {
0359         Int_t CluSysId = trk_system->at(iClu);
0360         Int_t CluLayId = trk_layer->at(iClu);
0361
0362         // if(CluSysId != 1) continue;
0363
0364         Double_t CluX = trk_X->at(iClu);
0365         Double_t CluY = trk_Y->at(iClu);
0366         Double_t CluZ = trk_Z->at(iClu);
0367
0368         Double_t CluR = sqrt(CluX*CluX + CluY*CluY);
0369         Double_t CluPhi = std::atan(CluY/CluX); // -pi/2 to pi/2
0370         if((CluPhi < 0)&&(CluX < 0)) CluPhi += TMath::Pi();
0371         else if((CluPhi > 0)&&(CluX < 0)) CluPhi -= TMath::Pi();
0372
0373         Double_t CluTheta = std::atan(CluR/CluZ);// -pi - pi
0374         Double_t CluEta = (CluZ/std::abs(CluZ)) * (- log(std::abs(std::tan(CluTheta/2))));
0375
0376         CluX = CluR*std::cos(CluPhi);
0377         CluY = CluR*std::sin(CluPhi);
0378         m_HINTTHitMap->Fill(CluX, CluY);
0379
0380         hitStruct cluHit;
0381         cluHit.r = CluR;
0382         cluHit.z = CluZ;
0383         cluHit.phi = CluPhi;
0384         cluHit.eta = CluEta;
0385
0386         if(CluSysId == 0) // mvtx
0387         {
0388             if(CluLayId == 0) vFMvtxHits.push_back(cluHit);
0389             else if(CluLayId == 1) vSMvtxHits.push_back(cluHit);
0390             else if(CluLayId == 2) vTMvtxHits.push_back(cluHit);
0391         }
0392         else if(CluSysId == 1) // intt
0393         {
0394             if((CluLayId == 3) || (CluLayId == 4)) vIInttHits.push_back(cluHit);
0395             else if((CluLayId == 5) || (CluLayId == 6)) vOInttHits.push_back(cluHit);
0396         }
0397         else vTpcHits.push_back(cluHit); // tpc
0398     }
0399 }
0400
0401 // calo tower
0402 void InttSeedTrackPerformance::ReadCalHitting(std::vector<hitStruct >& vEmcalHits,\
0403     std::vector<hitStruct >& vIHCalHits, std::vector<hitStruct >& vOHcalHits)
0404 {
0405     Int_t numOfCalTower = tower_Phi->size();
0406     if(numOfCalTower == 0) return;
0407
0408     for(Int_t iCalTower = 0; iCalTower < numOfCalTower; iCalTower++)
0409     {
0410         Int_t calId = tower_system->at(iCalTower);
0411
0412         Double_t caloX = tower_X->at(iCalTower);
0413         Double_t caloY = tower_Y->at(iCalTower);
0414         Double_t caloZ = tower_Z->at(iCalTower);
0415         Double_t caloPhi = tower_Phi->at(iCalTower);
0416
0417         Double_t caloEta = tower_Eta->at(iCalTower);
0418         Double_t caloE = tower_edep->at(iCalTower);
0419
0420         hitStruct caloHit;
0421         caloHit.r = std::sqrt(caloX*caloX + caloY*caloY);
0422         caloHit.z = caloZ;
0423         caloHit.phi = caloPhi;
0424         caloHit.eta = caloEta;
0425         caloHit.energy = caloE;
0426
0427         if(calId == 0) vEmcalHits.push_back(caloHit);
0428         else if(calId == 1) vIHCalHits.push_back(caloHit);
0429         else if(calId == 2) vOHcalHits.push_back(caloHit);
0430
0431         // if(caloE > 0.1){
0432         //    std::cout << "calE = " << caloE << ", bin(phi, eta) = " << tower_Phi_bin->at(iCalTower) << ", " << tower_Eta_bin->at(iCalTower) << std::endl;
0433         // }
0434     }
0435 }
0436
0437 // calo cluster
0438 void InttSeedTrackPerformance::ReadCalCluHitting(std::vector<hitStruct>& vEmcalHits,\
0439     std::vector<hitStruct>& vIHCalHits, std::vector<hitStruct>& vOHcalHits)
0440 {
0441     Int_t numOfCaloCluster = caloClus_Phi->size();
0442     if(numOfCaloCluster == 0) return;
0443
0444     for(Int_t iCaloCluster = 0; iCaloCluster < numOfCaloCluster; iCaloCluster++)
0445     {
0446         Int_t calId = caloClus_system->at(iCaloCluster);
0447
0448         Double_t caloX = caloClus_X->at(iCaloCluster);
0449         Double_t caloY = caloClus_Y->at(iCaloCluster);
0450         Double_t caloZ = caloClus_Z->at(iCaloCluster);
0451
0452         Double_t caloPhi = caloClus_Phi->at(iCaloCluster);
0453         Double_t caloE = caloClus_edep->at(iCaloCluster);
0454
0455         hitStruct caloHit;
0456         caloHit.r = std::sqrt(caloX*caloX + caloY*caloY);
0457         caloHit.z = caloZ;
0458         caloHit.phi = caloPhi;
0459
0460         Double_t caloTheta = std::atan(caloHit.r/caloHit.z);
0461         Double_t caloEta = (caloZ/std::abs(caloZ)) * (- log(std::abs(std::tan(caloTheta/2))));
0462         caloHit.eta = caloEta;
0463         caloHit.energy = caloE;
0464
0465         if(calId == 0) vEmcalHits.push_back(caloHit);
0466         else if(calId == 1) vIHCalHits.push_back(caloHit);
0467         else if(calId == 2) vOHcalHits.push_back(caloHit);
0468     }
0469 }
0470
0471
0472 // == s == Track QA Functions  ###############################################
0473 void InttSeedTrackPerformance::TrackQA(std::vector<hitStruct > vTruthPs, std::vector<hitStruct > vEmcalHits){
0474    std::vector<Int_t > matchiedRecoTrkId = {};
0475    Int_t numOfP = PrimaryG4P_PID->size();
0476    for(Int_t iP = 0; iP < numOfP; iP++){
0477       m_HTruTrackNum->Fill(vTruthPs.at(iP).pt);
0478       if(m_tracks.size() == 0) continue;
0479       Int_t trkId = TruRecoMatching(vTruthPs.at(iP), m_tracks, matchiedRecoTrkId);
0480       if(trkId == 9999) continue;
0481       // Int_t trkId = 0;
0482
0483       // m_bDecayEvent = ParticleDecayCheck(m_tracks.at(trkId), vEmcalHits);
0484
0485       DeltaPtPerform(vTruthPs.at(iP), m_tracks.at(trkId));
0486       EstimateMagneticShift(vTruthPs.at(iP), m_tracks.at(trkId));
0487       Double_t truE = vTruthPs.at(iP).energy;
0488       Double_t recoE = m_tracks.at(trkId).getTrackE();
0489       m_HDE->Fill(truE, (recoE - truE)/truE, vTruthPs.at(iP).eta);
0490       m_HPtEfficiency->Fill(vTruthPs.at(iP).pt);
0491    }
0492
0493    // std::cout << "size tru, reco = " << vTruthPs.size() << ", " << m_tracks.size() << std::endl;
0494 }
0495
0496 Int_t InttSeedTrackPerformance::TruRecoMatching(hitStruct truthP, std::vector<tracKuma > vRecoTrk,\
0497    std::vector<Int_t > vMatchiedRecoTrkId){
0498    Int_t matchiedRecoTrkId = 9999;
0499    Double_t closestDPt = 9999.;
0500    for(Int_t iRecoTrk = 0; iRecoTrk < vRecoTrk.size(); iRecoTrk++){
0501       Double_t tempTruTheta = 2*atan(std::exp(-truthP.eta));
0502       Double_t tempRecoTheta = vRecoTrk.at(iRecoTrk).getTrackTheta();
0503       Double_t tempRecoPhi = 0.;
0504       if(vRecoTrk.at(iRecoTrk).getHitIs(4)) tempRecoPhi = vRecoTrk.at(iRecoTrk).getHitPhi(4);
0505       else if(vRecoTrk.at(iRecoTrk).getHitIs(5)) tempRecoPhi = vRecoTrk.at(iRecoTrk).getHitPhi(5);
0506       // std::cout << "recoTheta, truTheta, dTheta = " << tempRecoTheta << ", " << tempTruTheta << ", " << std::abs(tempRecoTheta - tempTruTheta) << std::endl;
0507       if((std::abs(tempRecoTheta - tempTruTheta) < 0.1)){
0508          // std::cout << "dPhi = " << std::abs(tempRecoPhi - truthP.phi) << std::endl;
0509          if((std::abs(tempRecoPhi - truthP.phi) < 0.3)){
0510             // std::cout << "dPt = " << std::abs(vRecoTrk.at(iRecoTrk).getTrackPt() - truthP.pt) << std::endl;
0511             if(closestDPt > std::abs(vRecoTrk.at(iRecoTrk).getTrackPt() - truthP.pt)){
0512                bool matchedIs = std::find( vMatchiedRecoTrkId.begin(),  vMatchiedRecoTrkId.end(), \
0513                   iRecoTrk) !=  vMatchiedRecoTrkId.end();
0514                if(matchedIs) continue;
0515                closestDPt = std::abs(vRecoTrk.at(iRecoTrk).getTrackPt() - truthP.pt);
0516                matchiedRecoTrkId = matchedIs;
0517             }
0518          }
0519       }
0520    }
0521    vMatchiedRecoTrkId.push_back(matchiedRecoTrkId);
0522
0523    return matchiedRecoTrkId;
0524 }
0525
0526 void InttSeedTrackPerformance::DeltaPtPerform(hitStruct truthP, tracKuma trk){
0527    InttSeedTracking truckF;
0528
0529    Double_t truthPt = truthP.pt;
0530    Double_t truthEta = truthP.eta;
0531
0532    Double_t sagittaR = 9999.;
0533    Double_t centerX = 9999.;
0534    Double_t centerY = 9999.;
0535
0536    Double_t HitsXY[3][2];
0537    truckF.Set3PointsXY(HitsXY, trk, 0);
0538    truckF.RoughEstiSagittaCenter3Point(sagittaR, centerX, centerY, HitsXY);
0539    Double_t sagittaPt = truckF.CalcSagittaPt(sagittaR);
0540    Double_t tempdPt = (sagittaPt - truthPt)/truthPt;
0541    m_HTruthPtVsSagittaPt_rough->Fill(truthPt, tempdPt, truthEta);
0542
0543    std::vector<Int_t > subDetIds_VtxIntt = {0, 4, 5};
0544    std::vector<Double_t > hitsR_VtxIntt = {};
0545    std::vector<Double_t > hitsPhi_VtxIntt = {};
0546    if(truckF.ReturnHitsRPhiVect(hitsR_VtxIntt, hitsPhi_VtxIntt, subDetIds_VtxIntt, trk)){
0547       Double_t tempHitOInttPhi = trk.getHitPhi(4);
0548       Double_t tempHitEmcalPhi = trk.getHitPhi(6);
0549       truckF.SagittaRByCircleFit(centerX, centerY, sagittaR, hitsR_VtxIntt, hitsPhi_VtxIntt, tempHitOInttPhi, tempHitEmcalPhi);
0550       Double_t recoPt = truckF.CalcSagittaPt(sagittaR);
0551       Double_t dPt = (sagittaPt - truthPt)/truthPt;
0552       m_HTruthPtVsSagittaPt_VtxIntt->Fill(truthPt, dPt, truthEta);
0553    }
0554
0555    // Double_t HitsXY2[3][2];
0556    // truckF.Set3PointsXY(HitsXY2, trk, 0);
0557    // truckF.RoughEstiSagittaCenter3Point(sagittaR, centerX, centerY, HitsXY);
0558
0559    std::vector<Int_t > subDetIds_InttEmcal = {4, 5, 6};
0560    std::vector<Double_t > hitsR_InttEmcal = {};
0561    std::vector<Double_t > hitsPhi_InttEmcal = {};
0562    if(truckF.ReturnHitsRPhiVect(hitsR_InttEmcal, hitsPhi_InttEmcal, subDetIds_InttEmcal, trk)){
0563       Double_t tempHitOInttPhi = trk.getHitPhi(4);
0564       Double_t tempHitEmcalPhi = trk.getHitPhi(6);
0565       truckF.SagittaRByCircleFit(centerX, centerY, sagittaR, hitsR_InttEmcal, hitsPhi_InttEmcal, tempHitOInttPhi, tempHitEmcalPhi);
0566       Double_t recoPt = truckF.CalcSagittaPt(sagittaR);
0567       Double_t dPt = (recoPt - truthPt)/truthPt;
0568       m_HTruthPtVsSagittaPt_InttEmcal->Fill(truthPt, dPt, truthEta);
0569
0570       TrackOtherPropertiesWTruth(truthP, trk, sagittaR, centerX, centerY, recoPt,\
0571          m_HTruthPVsTheta_InttEmcal, m_HTruthPVsPhi_InttEmcal, m_HTruthPVsRecoP_InttEmcal, m_HTruthPVsEOverP_InttEmcal,\
0572             m_dVtxXY_InttEmcal, m_dVtxR_InttEmcal, m_dVtxZ_InttEmcal);
0573    }
0574
0575    std::vector<Int_t > subDetIds_VtxInttEmcal = {0, 4, 5, 6};
0576    std::vector<Double_t > hitsR_VtxInttEmcal = {};
0577    std::vector<Double_t > hitsPhi_VtxInttEmcal = {};
0578    if(truckF.ReturnHitsRPhiVect(hitsR_VtxInttEmcal, hitsPhi_VtxInttEmcal, subDetIds_VtxInttEmcal, trk)){
0579       Double_t tempHitOInttPhi = trk.getHitPhi(4);
0580       Double_t tempHitEmcalPhi = trk.getHitPhi(6);
0581       truckF.SagittaRByCircleFit(centerX, centerY, sagittaR, hitsR_VtxInttEmcal, hitsPhi_VtxInttEmcal, tempHitOInttPhi, tempHitEmcalPhi);
0582       Double_t recoPt = truckF.CalcSagittaPt(sagittaR);
0583       Double_t dPt = (recoPt - truthPt)/truthPt;
0584       m_HTruthPtVsSagittaPt_VtxInttEmcal->Fill(truthPt, dPt, truthEta);
0585    }
0586
0587    std::vector<Int_t > subDetIds_MvtxInttEmcal = {1, 2, 3, 4, 5, 6};
0588    std::vector<Double_t > hitsR_MvtxInttEmcal = {};
0589    std::vector<Double_t > hitsPhi_MvtxInttEmcal = {};
0590    if(truckF.ReturnHitsRPhiVect(hitsR_MvtxInttEmcal, hitsPhi_MvtxInttEmcal, subDetIds_MvtxInttEmcal, trk)){
0591       Double_t tempHitOInttPhi = trk.getHitPhi(4);
0592       Double_t tempHitEmcalPhi = trk.getHitPhi(6);
0593       truckF.SagittaRByCircleFit(centerX, centerY, sagittaR, hitsR_MvtxInttEmcal, hitsPhi_MvtxInttEmcal, tempHitOInttPhi, tempHitEmcalPhi);
0594       Double_t recoPt = truckF.CalcSagittaPt(sagittaR);
0595       Double_t dPt = (recoPt - truthPt)/truthPt;
0596       m_HTruthPtVsSagittaPt_MvtxInttEmcal->Fill(truthPt, dPt, truthEta);
0597       // std::cout << "pubEvNum = " << pubEvNum << std::endl;
0598       // EventJudge(pubEvNum, dPt, -2., 2., false); //CheckumaDaYo!!!
0599
0600       TrackOtherPropertiesWTruth(truthP, trk, sagittaR, centerX, centerY, recoPt,\
0601          m_HTruthPVsTheta_MvtxInttEmcal, m_HTruthPVsPhi_MvtxInttEmcal,\
0602          m_HTruthPVsRecoP_MvtxInttEmcal, m_HTruthPVsEOverP_MvtxInttEmcal,\
0603          m_dVtxXY_MvtxInttEmcal, m_dVtxR_MvtxInttEmcal, m_dVtxZ_MvtxInttEmcal);
0604
0605       // std::cout << "dPt = " << dPt << std::endl;
0606       // CheckumaDaYo!!!
0607       // ShowTrackInfo(trk, dPt, centerX, centerY, sagittaR);
0608    }
0609
0610    std::vector<Int_t > subDetIds_VtxMvtxInttEmcal = {0, 1, 2, 3, 4, 5, 6};
0611    std::vector<Double_t > hitsR_VtxMvtxInttEmcal = {};
0612    std::vector<Double_t > hitsPhi_VtxMvtxInttEmcal = {};
0613    if(truckF.ReturnHitsRPhiVect(hitsR_VtxMvtxInttEmcal, hitsPhi_VtxMvtxInttEmcal, subDetIds_VtxMvtxInttEmcal, trk)){
0614       Double_t tempHitOInttPhi = trk.getHitPhi(4);
0615       Double_t tempHitEmcalPhi = trk.getHitPhi(6);
0616       truckF.SagittaRByCircleFit(centerX, centerY, sagittaR, hitsR_VtxMvtxInttEmcal, hitsPhi_VtxMvtxInttEmcal, tempHitOInttPhi, tempHitEmcalPhi);
0617       Double_t recoPt = truckF.CalcSagittaPt(sagittaR);
0618       Double_t dPt = (recoPt - truthPt)/truthPt;
0619       m_HTruthPtVsSagittaPt_VtxMvtxInttEmcal->Fill(truthPt, dPt, truthEta);
0620    }
0621
0622    Double_t dPhiOInttEmcal = truckF.dPhiOInttEmcalEsti(trk);
0623    Double_t OriFunTrackPt = truckF.FitFunctionPt(dPhiOInttEmcal);
0624    Double_t dPtOriFun = (OriFunTrackPt - truthPt)/truthPt;
0625    m_HTruthPtVsFitFuncPt_IInttOInttEmcal->Fill(truthPt, dPtOriFun, truthEta);
0626    EventJudge(pubEvNum, dPtOriFun ,-2, -1, true);
0627    // std::cout << "111 dPt = " << dPtOriFun << std::endl;
0628    // CheckumaDaYo!! oooooooo
0629
0630    Double_t dPhiVtxIInttEmcal = truckF.dPhiVtxIInttEmcalEsti(trk);
0631    Double_t OriFunTrackPt_VtxIInttEmcal = truckF.FitFunctionPt_VtxIInttEmcal(dPhiOInttEmcal);
0632    Double_t dPtOriFun_VtxIInttEmcal = (OriFunTrackPt_VtxIInttEmcal - truthPt)/truthPt;
0633    m_HTruthPtVsFitFuncPt_VtxIInttEmcal->Fill(truthPt, dPtOriFun_VtxIInttEmcal, truthEta);
0634
0635    Double_t dPhiVtxOInttEmcal = truckF.dPhiVtxOInttEmcalEsti(trk);
0636    Double_t OriFunTrackPt_VtxOInttEmcal = truckF.FitFunctionPt_VtxOInttEmcal(dPhiVtxOInttEmcal);
0637    Double_t dPtOriFun_VtxOInttEmcal = (OriFunTrackPt_VtxOInttEmcal - truthPt)/truthPt;
0638    m_HTruthPtVsFitFuncPt_VtxOInttEmcal->Fill(truthPt, dPtOriFun_VtxOInttEmcal, truthEta);
0639
0640    Double_t dPhiVtxInttEmcal = truckF.dPhiVtxInttEmcalEsti(trk);
0641    Double_t OriFunTrackPt_VtxInttEmcal = truckF.FitFunctionPt_VtxInttEmcal(dPhiVtxInttEmcal);
0642    Double_t dPtOriFun_VtxInttEmcal = (OriFunTrackPt_VtxInttEmcal - truthPt)/truthPt;
0643    m_HTruthPtVsFitFuncPt_VtxInttEmcal->Fill(truthPt, dPtOriFun_VtxInttEmcal, truthEta);
0644
0645    Double_t dPhiMvtxIInttEmcal = truckF.dPhiMvtxIInttEmcalEsti(trk);
0646    Double_t OriFunTrackPt_MvtxIInttEmcal = truckF.FitFunctionPt_MVtxIInttEmcal(dPhiMvtxIInttEmcal);
0647    Double_t dPtOriFun_MvtxIInttEmcal = (OriFunTrackPt_MvtxIInttEmcal - truthPt)/truthPt;
0648    m_HTruthPtVsFitFuncPt_MvtxIInttEmcal->Fill(truthPt, dPtOriFun_MvtxIInttEmcal, truthEta);
0649
0650    Double_t dPhiMvtxOInttEmcal = truckF.dPhiMvtxOInttEmcalEsti(trk);
0651    Double_t OriFunTrackPt_MvtxOInttEmcal = truckF.FitFunctionPt_MVtxInttEmcal(dPhiMvtxOInttEmcal);
0652    Double_t dPtOriFun_MvtxOInttEmcal = (OriFunTrackPt_MvtxOInttEmcal - truthPt)/truthPt;
0653    m_HTruthPtVsFitFuncPt_MvtxOInttEmcal->Fill(truthPt, dPtOriFun_MvtxOInttEmcal, truthEta);
0654
0655    Double_t dPhiMvtxInttEmcal = truckF.dPhiMvtxInttEmcalEsti(trk);
0656    Double_t OriFunTrackPt_MvtxInttEmcal = truckF.FitFunctionPt_MVtxInttEmcal(dPhiMvtxInttEmcal);
0657    Double_t dPtOriFun_MvtxInttEmcal = (OriFunTrackPt_MvtxInttEmcal - truthPt)/truthPt;
0658    m_HTruthPtVsFitFuncPt_MvtxInttEmcal->Fill(truthPt, dPtOriFun_MvtxInttEmcal, truthEta);
0659 }
0660
0661
0662
0663 void InttSeedTrackPerformance::TrackOtherPropertiesWTruth(hitStruct truthP, tracKuma trk,\
0664    Double_t sagittaR, Double_t centerX, Double_t centerY, Double_t recoPt,\
0665    TH2D* hTruthPVsTheta, TH2D* hTruthPVsPhi, TH2D* hTruthPVsRecoP,TH2D* hTruthPVsEOverP,\
0666    TH2D* hDVtxXY, TH1D* hDVtxR, TH1D* hDVtxZ){
0667    InttSeedTracking truckF;
0668
0669    Double_t truthPt = truthP.pt;
0670    Double_t truthPhi = truthP.phi;
0671    Double_t truthTheta = 2*atan(std::exp(-truthP.eta));
0672    Double_t truVtxX = 0.; //TruthPV_trig_x
0673    Double_t truVtxY = 0.; //TruthPV_trig_y
0674    Double_t truVtxZ = 0.; //TruthPV_trig_y
0675
0676    Double_t recoTheta = truckF.EstimateRecoTheta(trk, 0);
0677
0678    Double_t trackE = trk.getTrackE();
0679    Double_t recoP = recoPt/sin(recoTheta);
0680
0681    Double_t vtxX = 0.;
0682    Double_t vtxY = 0.;
0683    truckF.CrossLineCircle(vtxX, vtxY, centerX, centerY, sagittaR);
0684    Double_t vtxR = std::sqrt(vtxX*vtxX + vtxY*vtxY);
0685    Double_t vtxPhi = atan(vtxY/vtxX);
0686    if((vtxPhi < 0)&&(vtxX < 0)) vtxPhi += TMath::Pi();
0687    else if((vtxPhi > 0)&&(vtxY < 0)) vtxPhi -= TMath::Pi();
0688
0689    Double_t trkPhi = -1/std::tan(vtxPhi);
0690    Double_t vtxZ = 0.;
0691    if(trk.getHitIs(6)){
0692       vtxZ = trk.getHitZ(6) - trk.getHitZ(6)/std::tan(recoTheta);
0693    }
0694
0695
0696    Double_t dVtxX = 9999.;
0697    Double_t dVtxY = 9999.;
0698    Double_t dVtxR = 9999.;
0699    Double_t dVtxZ = 9999.;
0700    dVtxX = vtxX - truVtxX;
0701    dVtxY = vtxY - truVtxY;
0702    dVtxR = std::sqrt(dVtxX*dVtxX + dVtxY*dVtxY);
0703    dVtxZ = vtxZ - truVtxZ;
0704    // std::cout << "dvtx r = " << dVtxR << std::endl;
0705
0706    hTruthPVsTheta->Fill(truthPt, recoTheta - truthTheta);
0707
0708    // hTruthPVsRecoP->Fill(recoP, (recoP - )/recoP);
0709    hTruthPVsEOverP->Fill(truthPt, trackE/recoP);
0710    hTruthPVsPhi->Fill(truthPt, trkPhi - truthPhi);
0711
0712    hDVtxXY->Fill(dVtxX, dVtxY);
0713    hDVtxR->Fill(dVtxR);
0714    hDVtxZ->Fill(dVtxZ);
0715
0716 }
0717
0718 void InttSeedTrackPerformance::EstimateMagneticShift(hitStruct truthP, tracKuma trk){
0719    InttSeedTracking truckF;
0720
0721    Double_t truthPt  = truthP.pt;
0722    Double_t truthEta = truthP.eta;
0723    Double_t truthPhi = truthP.phi;
0724
0725    Double_t refPhi = truthPhi;
0726    refPhi = TMath::Pi()/2 - refPhi;
0727    Double_t tempTruPhi = truthPhi + refPhi;
0728
0729
0730    Double_t fMvtxHitR = trk.getHitR(1);
0731    Double_t fMvtxHitPhi = trk.getHitPhi(1) + refPhi;
0732    Double_t sMvtxHitR = trk.getHitR(2);
0733    Double_t sMvtxHitPhi = trk.getHitPhi(2) + refPhi;
0734    Double_t tMvtxHitR = trk.getHitR(3);
0735    Double_t tMvtxHitPhi = trk.getHitPhi(3) + refPhi;
0736    Double_t MvtxHitR = (fMvtxHitR + sMvtxHitR + tMvtxHitR)/3;
0737    Double_t MvtxHitPhi = (fMvtxHitPhi + sMvtxHitPhi + tMvtxHitPhi)/3;
0738
0739    Double_t iInttHitR = trk.getHitR(4);
0740    Double_t iInttHitPhi = trk.getHitPhi(4) + refPhi;
0741    Double_t oInttHitR = trk.getHitR(5);
0742    Double_t oInttHitPhi = trk.getHitPhi(5) + refPhi;
0743    Double_t InttHitR = (iInttHitR + oInttHitR)/2;
0744    Double_t InttHitPhi = (iInttHitPhi + oInttHitPhi)/2;
0745
0746    Double_t emcalHitR = trk.getHitR(6);
0747    Double_t emcalHitPhi = trk.getHitPhi(6) + refPhi;
0748
0749
0750    Double_t dStraight_VtxIIntt = iInttHitR*sin(iInttHitPhi);
0751    Double_t dShiftL_VtxIIntt   = iInttHitR*cos(iInttHitPhi);
0752
0753    Double_t slopeIInttOIntt = 0.;
0754    Double_t sectionIInttOIntt = 0.;
0755    truckF.CalcLineEq(slopeIInttOIntt, sectionIInttOIntt, \
0756       iInttHitR*cos(iInttHitPhi), iInttHitR*sin(iInttHitPhi),\
0757          oInttHitR*cos(oInttHitPhi), oInttHitR*sin(oInttHitPhi));
0758
0759    Double_t crossX1 = 0.;
0760    Double_t crossY1 = 0.;
0761    Double_t slopeIINTT = dStraight_VtxIIntt/dShiftL_VtxIIntt;
0762    Double_t sectionIINTT = - slopeIINTT*iInttHitR*cos(iInttHitPhi) + iInttHitR*sin(iInttHitPhi);
0763    truckF.CalcCrossPoint(crossX1, crossY1, slopeIINTT, sectionIINTT, \
0764       oInttHitR*cos(oInttHitPhi), oInttHitR*sin(oInttHitPhi));
0765    Double_t dShiftL_VtxOIntt   = oInttHitR * cos(oInttHitPhi);
0766    Double_t dShiftL_IInttOIntt = dShiftL_VtxOIntt - dShiftL_VtxIIntt;
0767    Double_t dStraight_IInttOIntt = \
0768    std::sqrt((crossX1 - iInttHitR*cos(iInttHitPhi))*(crossX1 - iInttHitR*cos(iInttHitPhi))\
0769             + (crossY1 - iInttHitR*sin(iInttHitPhi))*(crossY1 - iInttHitR*sin(iInttHitPhi)));
0770
0771
0772    Double_t slopeOInttEmcal = 0.;
0773    Double_t sectionOInttEmcal = 0.;
0774    truckF.CalcLineEq(slopeOInttEmcal, sectionOInttEmcal, oInttHitR*cos(oInttHitPhi), oInttHitR*sin(oInttHitPhi),\
0775       emcalHitR*cos(emcalHitPhi), emcalHitR*sin(emcalHitPhi));
0776    Double_t crossX2 = 0.;
0777    Double_t crossY2 = 0.;
0778    truckF.CalcCrossPoint(crossX2, crossY2, slopeIInttOIntt, sectionIInttOIntt, \
0779       emcalHitR*cos(emcalHitPhi), emcalHitR*sin(emcalHitPhi));
0780    Double_t dStraight_OInttEmcal \
0781       = std::sqrt((crossX2 - oInttHitR*cos(oInttHitPhi))*(crossX2 - oInttHitR * cos(oInttHitPhi))\
0782             + (crossY2 - oInttHitR*sin(oInttHitPhi))*(crossY2 - oInttHitR * sin(oInttHitPhi)));
0783
0784
0785    Double_t dShiftL_VtxEmcal   = emcalHitR*cos(emcalHitPhi);
0786    Double_t dStraight_VtxEmcal = emcalHitR*sin(emcalHitPhi);
0787    Double_t dShiftL_OInttEmcal = dShiftL_VtxEmcal - dShiftL_VtxOIntt;
0788
0789    Double_t dPhi_VtxIIntt = TMath::Pi()/2 - iInttHitPhi;
0790    Double_t temp_dPhi_IInttOIntt = atan(dStraight_IInttOIntt/dShiftL_IInttOIntt);
0791    if(temp_dPhi_IInttOIntt < 0) temp_dPhi_IInttOIntt += TMath::Pi();
0792    Double_t dPhi_IInttOIntt = TMath::Pi()/2 - temp_dPhi_IInttOIntt - dPhi_VtxIIntt;
0793    Double_t temp_dPhi_OInttEmcal = atan(dStraight_OInttEmcal/dShiftL_OInttEmcal);
0794    if(temp_dPhi_OInttEmcal < 0) temp_dPhi_OInttEmcal += TMath::Pi();
0795    Double_t dPhi_OInttEmcal = TMath::Pi()/2 - temp_dPhi_OInttEmcal - dPhi_VtxIIntt - dPhi_IInttOIntt;
0796
0797    Double_t dPhiDStraight_VtxIIntt = dPhi_VtxIIntt/dStraight_VtxIIntt;
0798    m_HDPhiDStraight_VtxIIntt->Fill(truthPt, dPhiDStraight_VtxIIntt);
0799    m_HDPhiDStraight1D_VtxIIntt->Fill(dPhiDStraight_VtxIIntt);
0800
0801    Double_t dPhiDStraight_IInttOIntt = dPhi_IInttOIntt/dStraight_IInttOIntt;
0802    m_HDPhiDStraight_IInttOIntt->Fill(truthPt,  dPhiDStraight_IInttOIntt);
0803    m_HDPhiDStraight1D_IInttOIntt->Fill(dPhiDStraight_IInttOIntt);
0804
0805    Double_t dPhiDStraight_OInttEmcal = dPhi_OInttEmcal/dStraight_OInttEmcal;
0806    m_HDPhiDStraight_OInttEmcal->Fill(truthPt, dPhiDStraight_OInttEmcal);
0807    m_HDPhiDStraight1D_OInttEmcal->Fill(dPhiDStraight_OInttEmcal);
0808
0809    Double_t ddPhiddStraight_IInttOInttEmcal = dPhiDStraight_OInttEmcal - dPhiDStraight_IInttOIntt;
0810    m_HTruthPtVsDdPhiddStraight_IInttOInttEmcal->Fill(truthPt, ddPhiddStraight_IInttOInttEmcal);
0811
0812    m_HDPhiVsDStraightVsTruPt->Fill(dStraight_VtxIIntt, dPhi_VtxIIntt, truthPt);
0813    m_HDPhiVsDStraightVsTruPt->Fill(dStraight_IInttOIntt, dPhi_IInttOIntt, truthPt);
0814    m_HDPhiVsDStraightVsTruPt->Fill(dPhiDStraight_OInttEmcal, dPhi_OInttEmcal, truthPt);
0815
0816    m_HDLVsTruthPtVsEta_VtxIIntt->Fill(dShiftL_VtxIIntt, truthPt, truthEta);
0817    m_HDL1D_VtxIIntt->Fill(dShiftL_VtxIIntt);
0818    m_HDLVsTruthPtVsEta_IInttOIntt->Fill(dShiftL_IInttOIntt, truthPt, truthEta);
0819    m_HDL1D_IInttOIntt->Fill(dShiftL_IInttOIntt);
0820    m_HDLVsTruthPtVsEta_OInttEmcal->Fill(dShiftL_OInttEmcal, truthPt, truthEta);
0821    m_HDL1D_OInttEmcal->Fill(dShiftL_OInttEmcal);
0822
0823    Double_t tempMvtxX = MvtxHitR*std::cos(MvtxHitPhi);
0824    Double_t tempMvtxY = MvtxHitR*std::sin(MvtxHitPhi);
0825    Double_t tempIInttX = iInttHitR*std::cos(iInttHitPhi);
0826    Double_t tempIInttY = iInttHitR*std::sin(iInttHitPhi);
0827    Double_t tempOInttX = oInttHitR*std::cos(oInttHitPhi);
0828    Double_t tempOInttY = oInttHitR*std::sin(oInttHitPhi);
0829    Double_t tempInttX = InttHitR*std::cos(InttHitPhi);
0830    Double_t tempInttY = InttHitR*std::sin(InttHitPhi);
0831    Double_t tempEmcalX = emcalHitR*std::cos(emcalHitPhi);
0832    Double_t tempEmcalY = emcalHitR*std::sin(emcalHitPhi);
0833
0834    Double_t dPhiMvtxIIntt = std::atan((tempIInttY - tempMvtxY)/(tempIInttX - tempMvtxX));
0835    if(dPhiMvtxIIntt < 0) dPhiMvtxIIntt += TMath::Pi();
0836    Double_t dPhiVtxMvtxIIntt = dPhiMvtxIIntt - MvtxHitPhi;
0837
0838    Double_t dPhiMvtxOIntt = std::atan((tempOInttY - tempMvtxY)/(tempOInttX - tempMvtxX));
0839    if(dPhiMvtxOIntt < 0) dPhiMvtxOIntt += TMath::Pi();
0840    Double_t dPhiVtxMvtxOIntt = dPhiMvtxOIntt - MvtxHitPhi;
0841
0842    Double_t dPhiMvtxIntt = std::atan((tempInttY - tempMvtxY)/(tempInttX - tempMvtxX));
0843    if(dPhiMvtxIntt < 0) dPhiMvtxIntt += TMath::Pi();
0844    Double_t dPhiVtxMvtxIntt = dPhiMvtxIntt - MvtxHitPhi;
0845
0846
0847    // Double_t tempHitsXY[3][2];
0848    // truckF.Set3PointsXY(tempHitsXY, trk, 3);
0849    // Double_t tempR = 0.;
0850    // Double_t tempXc = 0.;
0851    // Double_t tempYc = 0.;
0852    // truckF.RoughEstiSagittaCenter3Point(tempR, tempXc, tempYc, tempHitsXY);
0853    // Double_t crossXcal = 0.;
0854    // Double_t crossYcal = 0.;
0855    // Double_t calPhi = truckF.CrossCircleCircle(crossXcal, crossYcal, tempXc, tempYc, tempR, trk.getHitPhi(5));
0856
0857    // if(std::abs(trk.getHitPhi(6) - calPhi) > 0.3) bDecayEvent = true;
0858    // std::cout << "x, y = " << crossXcal << ",  " << crossYcal << std::endl;
0859    // std::cout << "dPhi = " << std::abs(trk.getHitPhi(6) - calPhi) << std::endl;
0860    if(!m_bDecayEvent){
0861       Double_t dPhiIInttEmcal = std::atan((tempEmcalY - tempIInttY)/(tempEmcalX - tempIInttX));
0862       if(dPhiIInttEmcal < 0) dPhiIInttEmcal += TMath::Pi();
0863       Double_t dPhiVtxIInttEmcal = dPhiIInttEmcal - iInttHitPhi;
0864       Double_t dPhiMvtxIInttEmcal = dPhiIInttEmcal - dPhiMvtxIIntt;
0865
0866       Double_t dPhiOInttEmcal = std::atan((tempEmcalY - tempOInttY)/(tempEmcalX - tempOInttX));
0867       if(dPhiOInttEmcal < 0) dPhiOInttEmcal += TMath::Pi();
0868       Double_t dPhiVtxOInttEmcal = dPhiOInttEmcal - oInttHitPhi;
0869       Double_t dPhiMvtxOInttEmcal = dPhiOInttEmcal - dPhiMvtxOIntt;
0870
0871       Double_t dPhiInttEmcal = std::atan((tempEmcalY - tempInttY)/(tempEmcalX - tempInttX));
0872       if(dPhiInttEmcal < 0) dPhiInttEmcal += TMath::Pi();
0873       Double_t dPhiVtxInttEmcal = dPhiInttEmcal - InttHitPhi;
0874       Double_t dPhiMvtxInttEmcal = dPhiInttEmcal - dPhiMvtxIntt;
0875
0876       Double_t dPhiIInttOIntt = std::atan((tempOInttY - tempIInttY)/(tempOInttX - tempIInttX));
0877       if(dPhiIInttOIntt < 0) dPhiIInttOIntt += TMath::Pi();
0878       Double_t dPhi_IInttOInttEmcal = dPhiOInttEmcal - dPhiIInttOIntt;
0879
0880       m_HDPhiVsTruthPtVsEta_VtxIIntt->Fill(dPhi_VtxIIntt, truthPt, truthEta);
0881       m_HDPhiVsTruthPtVsEta_IInttOIntt->Fill(dPhi_IInttOIntt, truthPt, truthEta);
0882       m_HDPhiVsTruthPtVsEta_OInttEmcal->Fill(dPhi_OInttEmcal, truthPt, truthEta);
0883
0884       m_HDPhiVsTruthPtProf_VtxIIntt->Fill(dPhi_VtxIIntt, truthPt);
0885       m_HDPhiVsTruthPtProf_IInttOIntt->Fill(dPhi_IInttOIntt, truthPt);
0886       m_HDPhiVsTruthPtProf_OInttEmcal->Fill(dPhi_OInttEmcal, truthPt);
0887
0888       m_HDPhiVsTruthPtVsEta_IInttOInttEmcal->Fill(dPhi_IInttOInttEmcal, truthPt, truthEta);
0889       m_HDPhiVsTruthPtProf_IInttOInttEmcal->Fill(dPhi_IInttOInttEmcal, truthPt);
0890
0891       m_HDPhiVsTruthPtVsEta_VtxIInttEmcal->Fill(dPhiVtxIInttEmcal, truthPt, truthEta);
0892       m_HDPhiVsTruthPtVsEta_VtxOInttEmcal->Fill(dPhiVtxOInttEmcal, truthPt, truthEta);
0893       m_HDPhiVsTruthPtVsEta_VtxInttEmcal->Fill(dPhiVtxInttEmcal, truthPt, truthEta);
0894
0895       m_HDPhiVsTruthPtVsEta_MvtxIInttEmcal->Fill(dPhiMvtxIInttEmcal, truthPt, truthEta);
0896       m_HDPhiVsTruthPtVsEta_MvtxOInttEmcal->Fill(dPhiMvtxOInttEmcal, truthPt, truthEta);
0897       m_HDPhiVsTruthPtVsEta_MvtxInttEmcal->Fill(dPhiMvtxInttEmcal, truthPt, truthEta);
0898
0899       m_HDPhiVsTruthPtProf_VtxIInttEmcal->Fill(dPhiVtxIInttEmcal, truthPt);
0900       m_HDPhiVsTruthPtProf_VtxOInttEmcal->Fill(dPhiVtxOInttEmcal, truthPt);
0901       m_HDPhiVsTruthPtProf_VtxInttEmcal->Fill(dPhiVtxInttEmcal, truthPt);
0902
0903       m_HDPhiVsTruthPtProf_MvtxIInttEmcal->Fill(dPhiMvtxIInttEmcal, truthPt);
0904       m_HDPhiVsTruthPtProf_MvtxOInttEmcal->Fill(dPhiMvtxOInttEmcal, truthPt);
0905       m_HDPhiVsTruthPtProf_MvtxInttEmcal->Fill(dPhiMvtxInttEmcal, truthPt);
0906    }
0907
0908
0909 }
0910
0911 void InttSeedTrackPerformance::CaloEnergyQA(Int_t CalId, std::vector<hitStruct > calHits){
0912    Int_t highestEID = 9999;
0913    Double_t highestE = 0.;
0914    Int_t numOfCalT = calHits.size();
0915    for(Int_t iCalTower = 0; iCalTower < numOfCalT; iCalTower++){
0916       if(highestE < calHits.at(iCalTower).energy){
0917          highestE = calHits.at(iCalTower).energy;
0918          highestEID = iCalTower;
0919       }
0920    }
0921    if(highestEID == 9999) return;
0922    Double_t refCalPhi = calHits.at(highestEID).phi;
0923    Double_t refCalEta = calHits.at(highestEID).eta;
0924    for(Int_t iCalT = 0; iCalT < numOfCalT; iCalT++){
0925       Double_t phi = calHits.at(iCalT).phi - refCalPhi;
0926       Double_t eta = calHits.at(iCalT).eta - refCalEta;
0927       Double_t E = calHits.at(iCalT).energy;
0928       if(CalId == 0) m_HECalPhiVsE->Fill(phi, eta, E);
0929       else if(CalId == 1) m_HIHCalPhiVsE->Fill(phi, eta, E);
0930       else if(CalId == 2) m_HOHCalPhiVsE->Fill(phi, eta, E);
0931    }
0932 }
0933
0934
0935 bool InttSeedTrackPerformance::ParticleDecayCheck(tracKuma trk, std::vector<hitStruct > vEmcalHits){
0936    bool bIsDecay = false;
0937
0938    Int_t numOfCal = vEmcalHits.size();
0939    for(Int_t iCal = 0; iCal < numOfCal; iCal++){
0940       if((vEmcalHits.at(iCal).energy > 0.1)\
0941          && (std::abs(vEmcalHits.at(iCal).phi - trk.getHitPhi(6)) > 0.04)){
0942          bIsDecay = true;
0943       }
0944    }
0945
0946    return bIsDecay;
0947 }
0948 // == e == Track QA Functions  ###############################################
0949
0950
0951
0952 void InttSeedTrackPerformance::ResetValuesForEvent(){
0953    m_tracks.clear();
0954    m_tracks.shrink_to_fit();
0955
0956    m_TruthParticle.clear();
0957    m_TruthParticle.shrink_to_fit();
0958
0959    m_fMvtxHits.clear();
0960    m_sMvtxHits.clear();
0961    m_tMvtxHits.clear();
0962
0963    m_iInttHits.clear();
0964    m_oInttHits.clear();
0965
0966    m_TpcHits.clear();
0967
0968    m_emcalHits.clear();
0969    m_iHCalHits.clear();
0970    m_oHCalHits.clear();
0971
0972    // ==
0973    m_fMvtxHits.shrink_to_fit();
0974    m_sMvtxHits.shrink_to_fit();
0975    m_tMvtxHits.shrink_to_fit();
0976
0977    m_iInttHits.shrink_to_fit();
0978    m_oInttHits.shrink_to_fit();
0979
0980    m_TpcHits.shrink_to_fit();
0981
0982    m_emcalHits.shrink_to_fit();
0983    m_iHCalHits.shrink_to_fit();
0984    m_oHCalHits.shrink_to_fit();
0985
0986 }
0987
0988 void InttSeedTrackPerformance::AllResetValuesForEvent(){
0989    m_tracks.clear();
0990    m_tracks.shrink_to_fit();
0991
0992    m_vTruthParticle.clear();
0993    m_viInttHits.clear();
0994    m_voInttHits.clear();
0995    m_vTpcHits.clear();
0996    m_vemcalHits.clear();
0997    m_viHCalHits.clear();
0998    m_voHCalHits.clear();
0999
1000    m_vTruthParticle.shrink_to_fit();
1001    m_viInttHits.shrink_to_fit();
1002    m_voInttHits.shrink_to_fit();
1003    m_vTpcHits.shrink_to_fit();
1004    m_vemcalHits.shrink_to_fit();
1005    m_viHCalHits.shrink_to_fit();
1006    m_voHCalHits.shrink_to_fit();
1007
1008
1009    m_TruthParticle.clear();
1010    m_TruthParticle.shrink_to_fit();
1011
1012    m_fMvtxHits.clear();
1013    m_sMvtxHits.clear();
1014    m_tMvtxHits.clear();
1015
1016    m_iInttHits.clear();
1017    m_oInttHits.clear();
1018
1019    m_TpcHits.clear();
1020
1021    m_emcalHits.clear();
1022    m_iHCalHits.clear();
1023    m_oHCalHits.clear();
1024
1025    // ==
1026    m_fMvtxHits.shrink_to_fit();
1027    m_sMvtxHits.shrink_to_fit();
1028    m_tMvtxHits.shrink_to_fit();
1029
1030    m_iInttHits.shrink_to_fit();
1031    m_oInttHits.shrink_to_fit();
1032
1033    m_TpcHits.shrink_to_fit();
1034
1035    m_emcalHits.shrink_to_fit();
1036    m_iHCalHits.shrink_to_fit();
1037    m_oHCalHits.shrink_to_fit();
1038
1039 }
1040
1041 void InttSeedTrackPerformance::PartResetValuesForEvent(){
1042    m_TruthParticle.clear();
1043    m_TruthParticle.shrink_to_fit();
1044
1045    m_iInttHits.clear();
1046    m_oInttHits.clear();
1047
1048    m_TpcHits.clear();
1049
1050    m_emcalHits.clear();
1051    m_iHCalHits.clear();
1052    m_oHCalHits.clear();
1053
1054    // ==
1055    m_iInttHits.shrink_to_fit();
1056    m_oInttHits.shrink_to_fit();
1057
1058    m_TpcHits.shrink_to_fit();
1059
1060    m_emcalHits.shrink_to_fit();
1061    m_iHCalHits.shrink_to_fit();
1062    m_oHCalHits.shrink_to_fit();
1063 }
1064
1065
1066
1067 void InttSeedTrackPerformance::WrightHists(){
1068    oFile = new TFile(fOutput.c_str(), "recreate");
1069    oFile->cd();
1070    m_HINTTHitMap->Write();
1071
1072    m_HGlobalHitMap->Write();
1073
1074    m_HECalPhiVsE->Write();
1075    m_HIHCalPhiVsE->Write();
1076    m_HOHCalPhiVsE->Write();
1077
1078    m_HMatchCalEVsHighestCalE->Write();
1079
1080    m_HdPhiInttdPhiECal->Write();
1081    m_HTruthPtVsdPhiInttdPhiECal->Write();
1082
1083    m_HTruthPtVsSagittaPt_rough->Write();
1084    m_HTruthPtVsSagittaPt_VtxIntt->Write();
1085    m_HTruthPtVsSagittaPt_InttEmcal->Write();
1086    m_HTruthPtVsSagittaPt_VtxInttEmcal->Write();
1087    m_HTruthPtVsSagittaPt_MvtxInttEmcal->Write();
1088    m_HTruthPtVsSagittaPt_VtxMvtxInttEmcal->Write();
1089
1090    m_HTruthPtVsFitFuncPt_IInttOInttEmcal->Write();
1091
1092    m_HTruthPtVsFitFuncPt_VtxIInttEmcal->Write();
1093    m_HTruthPtVsFitFuncPt_VtxOInttEmcal->Write();
1094    m_HTruthPtVsFitFuncPt_VtxInttEmcal->Write();
1095
1096    m_HTruthPtVsFitFuncPt_MvtxIInttEmcal->Write();
1097    m_HTruthPtVsFitFuncPt_MvtxOInttEmcal->Write();
1098    m_HTruthPtVsFitFuncPt_MvtxInttEmcal->Write();
1099
1100    m_HTruthPVsRecoP_InttEmcal->Write();
1101    m_HTruthPVsRecoP_MvtxInttEmcal->Write();
1102    m_HTruthPVsRecoP_FitFunc->Write();
1103
1104    m_HTruthPVsEOverP_InttEmcal->Write();
1105    m_HTruthPVsEOverP_MvtxInttEmcal->Write();
1106    m_HTruthPVsEOverP_FitFunc->Write();
1107
1108    m_HTruthPVsTheta_InttEmcal->Write();
1109    m_HTruthPVsTheta_MvtxInttEmcal->Write();
1110
1111    m_HTruthPVsPhi_InttEmcal->Write();
1112    m_HTruthPVsPhi_MvtxInttEmcal->Write();
1113
1114
1115    m_dVtxXY_InttEmcal->Write();
1116    m_dVtxXY_MvtxInttEmcal->Write();
1117
1118    m_dVtxR_InttEmcal->Write();
1119    m_dVtxR_MvtxInttEmcal->Write();
1120
1121    m_dVtxZ_InttEmcal->Write();
1122    m_dVtxZ_MvtxInttEmcal->Write();
1123
1124    m_HDE->Write();
1125
1126    m_HDPhiVsDStraightVsTruPt->Write();
1127
1128    m_HDPhiDStraight_VtxIIntt->Write();
1129    m_HDPhiDStraight1D_VtxIIntt->Write();
1130    m_HDPhiDStraight_IInttOIntt->Write();
1131    m_HDPhiDStraight1D_IInttOIntt->Write();
1132    m_HDPhiDStraight_OInttEmcal->Write();
1133    m_HDPhiDStraight1D_OInttEmcal->Write();
1134    m_HTruthPtVsDdPhiddStraight_IInttOInttEmcal->Write();
1135
1136
1137    m_HDLVsTruthPtVsEta_VtxIIntt->Write();
1138    m_HDL1D_VtxIIntt->Write();
1139    m_HDLVsTruthPtVsEta_IInttOIntt->Write();
1140    m_HDL1D_IInttOIntt->Write();
1141    m_HDLVsTruthPtVsEta_OInttEmcal->Write();
1142    m_HDL1D_OInttEmcal->Write();
1143
1144    m_HDPhiVsTruthPtVsEta_VtxIIntt->Write();
1145    m_HDPhiVsTruthPtVsEta_IInttOIntt->Write();
1146    m_HDPhiVsTruthPtVsEta_OInttEmcal->Write();
1147
1148    m_HDPhiVsTruthPtVsEta_IInttOInttEmcal->Write();
1149
1150    m_HDPhiVsTruthPtVsEta_VtxIInttEmcal->Write();
1151    m_HDPhiVsTruthPtVsEta_VtxOInttEmcal->Write();
1152    m_HDPhiVsTruthPtVsEta_VtxInttEmcal->Write();
1153    m_HDPhiVsTruthPtVsEta_MvtxIInttEmcal->Write();
1154    m_HDPhiVsTruthPtVsEta_MvtxOInttEmcal->Write();
1155    m_HDPhiVsTruthPtVsEta_MvtxInttEmcal->Write();
1156    m_HDPhiVsTruthPtProf_VtxIIntt->Write();
1157
1158    m_HDPhiVsTruthPtProf_IInttOInttEmcal->Write();
1159    m_HDPhiVsTruthPtProf_IInttOIntt->Write();
1160    m_HDPhiVsTruthPtProf_OInttEmcal->Write();
1161    m_HDPhiVsTruthPtProf_VtxIInttEmcal->Write();
1162    m_HDPhiVsTruthPtProf_VtxOInttEmcal->Write();
1163    m_HDPhiVsTruthPtProf_VtxInttEmcal->Write();
1164    m_HDPhiVsTruthPtProf_MvtxIInttEmcal->Write();
1165    m_HDPhiVsTruthPtProf_MvtxOInttEmcal->Write();
1166    m_HDPhiVsTruthPtProf_MvtxInttEmcal->Write();
1167
1168    m_HPtEfficiency->Divide(m_HPtEfficiency, m_HTruTrackNum, 1.0, 1.0, "b");
1169    m_HPtEfficiency->Write();
1170    m_HTruTrackNum->Write();
1171
1172    oFile->Close();
1173
1174    std::cout << "root " << fOutput << std::endl;
1175 }
1176
1177 // MakeClass Original Functions  ###############################################
1178 InttSeedTrackPerformance::InttSeedTrackPerformance(TTree *tree, std::string fInputName, std::string fOutputName, Int_t runNum) : fChain(0)
1179 {
1180 // if parameter tree is not specified (or zero), connect the file
1181 // used to generate this class and read the Tree.
1182    fInput = fInputName;
1183    fOutput = fOutputName + std::to_string(runNum) + ".root";
1184
1185    if (tree == 0) {
1186       TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject(fInput.c_str());
1187       if (!f || !f->IsOpen()) {
1188          f = new TFile(fInput.c_str());
1189       }
1190       f->GetObject("tree",tree);
1191
1192    }
1193    Init(tree);
1194 }
1195
1196 InttSeedTrackPerformance::~InttSeedTrackPerformance()
1197 {
1198    if (!fChain) return;
1199    delete fChain->GetCurrentFile();
1200 }
1201
1202 Int_t InttSeedTrackPerformance::GetEntry(Long64_t entry)
1203 {
1204 // Read contents of entry.
1205    if (!fChain) return 0;
1206    return fChain->GetEntry(entry);
1207 }
1208
1209 Long64_t InttSeedTrackPerformance::LoadTree(Long64_t entry)
1210 {
1211 // Set the environment to read one entry
1212    if (!fChain) return -5;
1213    Long64_t centry = fChain->LoadTree(entry);
1214    if (centry < 0) return centry;
1215    if (fChain->GetTreeNumber() != fCurrent) {
1216       fCurrent = fChain->GetTreeNumber();
1217    }
1218    return centry;
1219 }
1220
1221 void InttSeedTrackPerformance::Init(TTree *tree)
1222 {
1223    // The Init() function is called when the selector needs to initialize
1224    // a new tree or chain. Typically here the branch addresses and branch
1225    // pointers of the tree will be set.
1226    // It is normally not necessary to make changes to the generated
1227    // code, but the routine can be extended by the user if needed.
1228    // Init() will be called many times when running on PROOF
1229    // (once per file to be processed).
1230
1231    // Set object pointer
1232    trk_system = 0;
1233    trk_layer = 0;
1234    trk_adc = 0;
1235    trk_X = 0;
1236    trk_Y = 0;
1237    trk_Z = 0;
1238    trk_size = 0;
1239    trk_phi_size = 0;
1240    trk_Z_size = 0;
1241
1242    tower_system = 0;
1243    tower_X = 0;
1244    tower_Y = 0;
1245    tower_Z = 0;
1246    tower_Eta = 0;
1247    tower_Phi = 0;
1248    tower_Eta_test = 0;
1249    tower_Phi_test = 0;
1250    tower_Eta_bin = 0;
1251    tower_Phi_bin = 0;
1252    tower_edep = 0;
1253
1254    caloClus_system = 0;
1255    caloClus_X = 0;
1256    caloClus_Y = 0;
1257    caloClus_Z = 0;
1258    caloClus_R = 0;
1259    caloClus_Phi = 0;
1260    caloClus_edep = 0;
1261
1262    PrimaryG4P_Pt = 0;
1263    PrimaryG4P_Eta = 0;
1264    PrimaryG4P_Phi = 0;
1265    PrimaryG4P_E = 0;
1266    PrimaryG4P_PID = 0;
1267    PrimaryG4P_isChargeHadron = 0;
1268    // Set branch addresses and branch pointers
1269    if (!tree) return;
1270    fChain = tree;
1271    fCurrent = -1;
1272    fChain->SetMakeClass(1);
1273
1274    fChain->SetBranchAddress("trk_NClus", &trk_NClus, &b_trk_NClus);
1275    fChain->SetBranchAddress("trk_system", &trk_system, &b_trk_system);
1276    fChain->SetBranchAddress("trk_layer", &trk_layer, &b_trk_layer);
1277    fChain->SetBranchAddress("trk_adc", &trk_adc, &b_trk_adc);
1278    fChain->SetBranchAddress("trk_X", &trk_X, &b_trk_X);
1279    fChain->SetBranchAddress("trk_Y", &trk_Y, &b_trk_Y);
1280    fChain->SetBranchAddress("trk_Z", &trk_Z, &b_trk_Z);
1281    fChain->SetBranchAddress("trk_size", &trk_size, &b_trk_size);
1282    fChain->SetBranchAddress("trk_phi_size", &trk_phi_size, &b_trk_phi_size);
1283    fChain->SetBranchAddress("trk_Z_size", &trk_Z_size, &b_trk_Z_size);
1284
1285    fChain->SetBranchAddress("nTowers", &nTowers, &b_nTowers);
1286    fChain->SetBranchAddress("tower_system", &tower_system, &b_tower_system);
1287    fChain->SetBranchAddress("tower_X", &tower_X, &b_tower_X);
1288    fChain->SetBranchAddress("tower_Y", &tower_Y, &b_tower_Y);
1289    fChain->SetBranchAddress("tower_Z", &tower_Z, &b_tower_Z);
1290    fChain->SetBranchAddress("tower_Eta", &tower_Eta, &b_tower_Eta);
1291    fChain->SetBranchAddress("tower_Phi", &tower_Phi, &b_tower_Phi);
1292    fChain->SetBranchAddress("tower_Eta_test", &tower_Eta_test, &b_tower_Eta_test);
1293    fChain->SetBranchAddress("tower_Phi_test", &tower_Phi_test, &b_tower_Phi_test);
1294    fChain->SetBranchAddress("tower_Eta_bin", &tower_Eta_bin, &b_tower_Eta_bin);
1295    fChain->SetBranchAddress("tower_Phi_bin", &tower_Phi_bin, &b_tower_Phi_bin);
1296    fChain->SetBranchAddress("tower_edep", &tower_edep, &b_tower_edep);
1297
1298    //checkumaDaYooo!!
1299    if(bCaloClu){
1300       fChain->SetBranchAddress("nCaloClus", &nCaloClus, &b_nCaloClus);
1301       fChain->SetBranchAddress("caloClus_system", &caloClus_system, &b_caloClus_system);
1302       fChain->SetBranchAddress("caloClus_X", &caloClus_X, &b_caloClus_X);
1303       fChain->SetBranchAddress("caloClus_Y", &caloClus_Y, &b_caloClus_Y);
1304       fChain->SetBranchAddress("caloClus_Z", &caloClus_Z, &b_caloClus_Z);
1305       fChain->SetBranchAddress("caloClus_R", &caloClus_R, &b_caloClus_R);
1306       fChain->SetBranchAddress("caloClus_Phi", &caloClus_Phi, &b_caloClus_Phi);
1307       fChain->SetBranchAddress("caloClus_edep", &caloClus_edep, &b_caloClus_edep);
1308    }
1309
1310
1311    fChain->SetBranchAddress("NTruthVtx", &NTruthVtx, &b_NTruthVtx);
1312    fChain->SetBranchAddress("TruthPV_trig_x", &TruthPV_trig_x, &b_TruthPV_trig_x);
1313    fChain->SetBranchAddress("TruthPV_trig_y", &TruthPV_trig_y, &b_TruthPV_trig_y);
1314    fChain->SetBranchAddress("TruthPV_trig_z", &TruthPV_trig_z, &b_TruthPV_trig_z);
1315    fChain->SetBranchAddress("NPrimaryG4P", &NPrimaryG4P, &b_NPrimaryG4P);
1316    fChain->SetBranchAddress("PrimaryG4P_Pt", &PrimaryG4P_Pt, &b_PrimaryG4P_Pt);
1317    fChain->SetBranchAddress("PrimaryG4P_Eta", &PrimaryG4P_Eta, &b_PrimaryG4P_Eta);
1318    fChain->SetBranchAddress("PrimaryG4P_Phi", &PrimaryG4P_Phi, &b_PrimaryG4P_Phi);
1319    fChain->SetBranchAddress("PrimaryG4P_E", &PrimaryG4P_E, &b_PrimaryG4P_E);
1320    fChain->SetBranchAddress("PrimaryG4P_PID", &PrimaryG4P_PID, &b_PrimaryG4P_PID);
1321    fChain->SetBranchAddress("PrimaryG4P_isChargeHadron", &PrimaryG4P_isChargeHadron, &b_PrimaryG4P_isChargeHadron);
1322
1323 }
1324 // MakeClass Original Functions  ###############################################
1325
1326
1327 void InttSeedTrackPerformance::DrawEventDisplay(Int_t eventId, std::vector<tracKuma > trks){
1328    if(qaEvent != 9999){
1329       TCanvas* c = new TCanvas("c","c: pT",800,850);
1330       c->cd();
1331       TPad* pad1 = new TPad("pad1", "pad1", 0.1, 0.1, 1, 1.0);
1332       pad1->SetBottomMargin(0.15);
1333       pad1->SetLeftMargin(0.15);
1334       pad1->SetRightMargin(0.05);
1335       pad1->SetTopMargin(0.05);
1336       pad1->Draw();
1337       pad1->cd();
1338
1339       TH2D* hHitMapED = new TH2D( "hHitMapED", "hHitMapED", 200, -100, 100, 200, -100, 100);
1340       std::vector<TArc* > trkLines;
1341       for(Int_t iTrk = 0; iTrk < trks.size(); iTrk++){
1342          tracKuma trk = trks.at(iTrk);
1343          for(Int_t iHit = 1; iHit < 7; iHit++){
1344             Double_t tempX = trk.getHitR(iHit)*std::cos(trk.getHitPhi(iHit));
1345             Double_t tempY = trk.getHitR(iHit)*std::sin(trk.getHitPhi(iHit));
1346             hHitMapED->Fill(tempX, tempY);
1347          }
1348
1349          Double_t tempR = trk.getTrackSagR();
1350          Double_t tempXc = trk.getTrackSagX();
1351          Double_t tempYc = trk.getTrackSagY();
1352
1353          Double_t edgeX1 = 0.;
1354          Double_t edgeY1 = 0.;
1355          Double_t edgePhi1 = std::atan((edgeY1 - tempYc)/(edgeX1 - tempXc));
1356          if((edgePhi1 < 0)&&((edgeX1 - tempXc) < 0)) edgePhi1 += TMath::Pi();
1357          else if((edgePhi1 > 0)&&((edgeX1 - tempXc) < 0))  edgePhi1 += TMath::Pi();
1358
1359          Double_t edgeX2 = trk.getHitR(6)*std::cos(trk.getHitPhi(6));
1360          Double_t edgeY2 = trk.getHitR(6)*std::sin(trk.getHitPhi(6));
1361          Double_t edgePhi2 = std::atan((edgeY2 - tempYc)/(edgeX2 - tempXc));
1362          if((edgePhi2 < 0)&&((edgeX2 - tempXc) < 0)) edgePhi2 += TMath::Pi();
1363          else if((edgePhi2 > 0)&&((edgeX2 - tempXc) < 0))  edgePhi2 += TMath::Pi();
1364
1365          Double_t diffPhi1 = std::abs(edgePhi2 - edgePhi1);
1366          Double_t diffPhi2 = std::abs((edgePhi2 - 2*TMath::Pi()) - edgePhi1);
1367          if(diffPhi2 < diffPhi1) edgePhi2 -= 2*TMath::Pi();
1368
1369          Double_t tempMinPhi = 0.;
1370          Double_t tempMaxPhi = 360.;
1371          if(edgePhi1 < edgePhi2){
1372             tempMinPhi = edgePhi1 - 10.;
1373             tempMaxPhi = edgePhi2 + 10.;
1374          }else{
1375             tempMinPhi = edgePhi2 - 10.;
1376             tempMaxPhi = edgePhi1 + 10.;
1377          }
1378
1379          auto arc = new TArc(tempXc, tempYc, tempR, tempMinPhi, tempMaxPhi);
1380          arc->SetFillStyle(0);
1381
1382          trkLines.push_back(arc);
1383       }
1384
1385       hHitMapED->SetFillColor(3);
1386       hHitMapED->Draw("BOX");
1387       for(Int_t iTrk = 0; iTrk < trkLines.size(); iTrk++){
1388          trkLines.at(iTrk)->Draw("only");
1389       }
1390
1391       c->SaveAs("hQaSagittaFit.root");
1392       c->Close();
1393       // oFile->cd();
1394       // hHitMap->Write();
1395       // fCircleQA->Write();
1396       qaEvent = 9999;
1397    }
1398 }
1399
1400 void InttSeedTrackPerformance::ShowEventInfo(){
1401    std::cout << "TH2D* hHitMapMvtx = new TH2D( \"hHitMapMvtx\", \"hHitMap\", 2000, -100, 100, 2000, -100, 100);" << std::endl;
1402    std::cout << "TH2D* hHitMapIntt = new TH2D( \"hHitMapIntt\", \"hHitMap\", 2000, -100, 100, 2000, -100, 100);" << std::endl;
1403    std::cout << "TH2D* hHitMapTpc = new TH2D( \"hHitMapTpc\", \"hHitMap\", 2000, -100, 100, 2000, -100, 100);" << std::endl;
1404    std::cout << "TH2D* hHitMapEmcal = new TH2D( \"hHitMapEmcal\", \"hHitMap\", 2000, -100, 100, 2000, -100, 100);" << std::endl;
1405    std::cout << "hHitMapMvtx->SetLineWidth(8);" << std::endl;
1406    std::cout << "hHitMapIntt->SetLineWidth(8);" << std::endl;
1407    std::cout << "hHitMapTpc->SetLineWidth(8);" << std::endl;
1408    std::cout << "hHitMapEmcal->SetLineWidth(8);" << std::endl;
1409    std::cout << "hHitMapMvtx->SetLineColor(860);" << std::endl;
1410    std::cout << "hHitMapIntt->SetLineColor(875);" << std::endl;
1411    std::cout << "hHitMapTpc->SetLineColor(870);" << std::endl;
1412    std::cout << "hHitMapEmcal->SetLineColor(807);" << std::endl;
1413
1414
1415    for(Int_t i=0; i<m_fMvtxHits.size(); i++){
1416       std::cout << "hHitMapMvtx->Fill(" << m_fMvtxHits.at(i).r*std::cos(m_fMvtxHits.at(i).phi) \
1417          << ", " << m_fMvtxHits.at(i).r*std::sin(m_fMvtxHits.at(i).phi) << ");" << std::endl;
1418    }
1419    for(Int_t i=0; i<m_sMvtxHits.size(); i++){
1420       std::cout << "hHitMapMvtx->Fill(" << m_sMvtxHits.at(i).r*std::cos(m_sMvtxHits.at(i).phi) \
1421          << ", " << m_sMvtxHits.at(i).r*std::sin(m_sMvtxHits.at(i).phi) << ");" << std::endl;
1422    }
1423    for(Int_t i=0; i<m_tMvtxHits.size(); i++){
1424       std::cout << "hHitMapMvtx->Fill(" << m_tMvtxHits.at(i).r*std::cos(m_tMvtxHits.at(i).phi) \
1425          << ", " << m_tMvtxHits.at(i).r*std::sin(m_tMvtxHits.at(i).phi) << ");" << std::endl;
1426    }
1427    for(Int_t i=0; i< m_iInttHits.size(); i++){
1428       std::cout << "hHitMapIntt->Fill(" << m_iInttHits.at(i).r*std::cos(m_iInttHits.at(i).phi) \
1429          << ", " << m_iInttHits.at(i).r*std::sin(m_iInttHits.at(i).phi) << ");" << std::endl;
1430    }
1431    for(Int_t i=0; i< m_oInttHits.size(); i++){
1432       std::cout << "hHitMapIntt->Fill(" << m_oInttHits.at(i).r*std::cos(m_oInttHits.at(i).phi) \
1433          << ", " << m_oInttHits.at(i).r*std::sin(m_oInttHits.at(i).phi) << ");" << std::endl;
1434    }
1435    for(Int_t i=0; i< m_TpcHits.size(); i++){
1436       std::cout << "hHitMapTpc->Fill(" << m_TpcHits.at(i).r*std::cos(m_TpcHits.at(i).phi) \
1437          << ", " << m_TpcHits.at(i).r*std::sin(m_TpcHits.at(i).phi) << ");" << std::endl;
1438    }
1439    for(Int_t i=0; i< m_emcalHits.size(); i++){
1440       if(m_emcalHits.at(i).energy > 0.1) std::cout << "hHitMapEmcal->Fill(" << m_emcalHits.at(i).r*std::cos(m_emcalHits.at(i).phi) \
1441          << ", " << m_emcalHits.at(i).r*std::sin(m_emcalHits.at(i).phi) << ");" << std::endl;
1442    }
1443
1444    std::cout << "hHitMapMvtx->Draw(\"BOX\");" << std::endl;
1445    std::cout << "hHitMapIntt->Draw(\"BOX same\");" << std::endl;
1446    std::cout << "hHitMapTpc->Draw(\"BOX same\");" << std::endl;
1447    std::cout << "hHitMapEmcal->Draw(\"BOX same\");" << std::endl;
1448 }
1449
1450
1451 void InttSeedTrackPerformance::ShowTrackInfo(tracKuma trk, Double_t dPt, Double_t centerX, Double_t centerY, Double_t sagittaR){
1452    std::cout << "222: dPt, R, Xc, Yc = " << dPt << ", " << sagittaR << ", " << centerX << ", " << centerY << std::endl;
1453    std::cout << "auto arc = new  TArc(" << centerX << ", " << centerY << ", " << sagittaR << ");"  << std::endl;
1454    std::cout << "arc->SetFillStyle(0);" << std::endl;
1455    std::cout << "arc->Draw(\"only\");" << std::endl;
1456
1457
1458    std::cout << "TH2D* hHitMapTrk = new TH2D( \"hHitMapTrk\", \"hHitMapTrk\", 2000, -100, 100, 2000, -100, 100);" << std::endl;
1459    std::cout << "hHitMapTrk->SetLineWidth(3);" << std::endl;
1460    std::cout << "hHitMapTrk->SetLineColor(632);" << std::endl;
1461
1462    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(1)*std::cos(trk.getHitPhi(1)) \
1463             << ", " << trk.getHitR(1)*std::sin(trk.getHitPhi(1)) << ");" << std::endl;
1464    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(2)*std::cos(trk.getHitPhi(2)) \
1465             << ", " << trk.getHitR(2)*std::sin(trk.getHitPhi(2)) << ");" << std::endl;
1466    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(3)*std::cos(trk.getHitPhi(3)) \
1467             << ", " << trk.getHitR(3)*std::sin(trk.getHitPhi(3)) << ");" << std::endl;
1468    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(4)*std::cos(trk.getHitPhi(4)) \
1469             << ", " << trk.getHitR(4)*std::sin(trk.getHitPhi(4)) << ");" << std::endl;
1470    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(5)*std::cos(trk.getHitPhi(5)) \
1471             << ", " << trk.getHitR(5)*std::sin(trk.getHitPhi(5)) << ");" << std::endl;
1472    std::cout << "hHitMapTrk->Fill(" << trk.getHitR(6)*std::cos(trk.getHitPhi(6)) \
1473             << ", " << trk.getHitR(6)*std::sin(trk.getHitPhi(6)) << ");" << std::endl;
1474
1475    std::cout << "hHitMapTrk->Draw(\"BOX same\");" << std::endl;
1476 }
1477
1478 void InttSeedTrackPerformance::EventJudge(Int_t eventNum, Double_t targetVal, Double_t minLim, Double_t maxLim, bool bIn){
1479    if(bIn){
1480       if((targetVal > minLim)&&(targetVal < maxLim)){
1481          m_vTargetEvents.push_back(eventNum);
1482       }
1483    }else{
1484       if((targetVal < minLim)||(targetVal > maxLim)){
1485          m_vTargetEvents.push_back(eventNum);
1486       }
1487    }
1488
1489 }
1490
1491
1492 void InttSeedTrackPerformance::ChecKuma(std::string checkNo){
1493     std::cout <<"CheeeeeeeeeeeeeeeeeeeeeeeeeeeeeeecKuma"<< checkNo << std::endl;
1494 }
1495
1496 #endif // #ifdef InttSeedTrackPerformance_cxx