TrackletHistogram/old_code/TrackletHistogramFill.cpp

0001 #include "TrackletHistogramFill.h"
0002
0003 TrackletHistogramFill::TrackletHistogramFill(
0004     int process_id_in,
0005     int runnumber_in,
0006
0007     std::string output_directory_in,
0008     std::string output_file_name_suffix_in
0009 ):
0010     process_id(process_id_in),
0011     runnumber(runnumber_in),
0012     output_directory(output_directory_in),
0013     output_file_name_suffix(output_file_name_suffix_in),
0014     isWithRotate(true)
0015 {
0016     nCentrality_bin = centrality_edges.size() - 1;
0017     vtxZReweightFactor = 1.;
0018 }
0019
0020
0021 void TrackletHistogramFill::PrepareOutPutFileName()
0022 {
0023     std::string job_index = std::to_string( process_id );
0024     int job_index_len = 5;
0025     job_index.insert(0, job_index_len - job_index.size(), '0');
0026
0027     std::string runnumber_str = std::to_string( runnumber );
0028     if (runnumber != -1){
0029         int runnumber_str_len = 8;
0030         runnumber_str.insert(0, runnumber_str_len - runnumber_str.size(), '0');
0031     }
0032
0033     if (output_file_name_suffix.size() > 0 && output_file_name_suffix[0] != '_') {
0034         output_file_name_suffix = "_" + output_file_name_suffix;
0035     }
0036
0037     output_filename = "TrackletHistogram";
0038     output_filename = (runnumber != -1) ? "Data_" + output_filename : "MC_" + output_filename;
0039     output_filename += (vtxZReweight) ? "_vtxZReweight" : "";
0040     output_filename += (BcoFullDiffCut && runnumber != -1) ? "_BcoFullDiffCut" : "";
0041     output_filename += (INTT_vtxZ_QA) ? "_INTT_vtxZ_QA" : "";
0042     output_filename += (isClusQA.first) ? Form("_ClusQAAdc%.0fPhiSize%.0f",isClusQA.second.first,isClusQA.second.second) : "";
0043     output_filename += output_file_name_suffix;
0044     output_filename += (runnumber != -1) ? Form("_%s_%s.root",runnumber_str.c_str(),job_index.c_str()) : Form("_%s.root",job_index.c_str());
0045 }
0046
0047 void TrackletHistogramFill::PrepareOutPutRootFile()
0048 {
0049     file_out = new TFile(Form("%s/%s", output_directory.c_str(), output_filename.c_str()), "RECREATE");
0050     tree_out_par = new TTree("tree_par", "used parameters");
0051
0052     tree_out_par -> Branch("process_id", &process_id);
0053     tree_out_par -> Branch("runnumber", &runnumber);
0054     tree_out_par -> Branch("vtxZReweight", &vtxZReweight);
0055     tree_out_par -> Branch("BcoFullDiffCut", &BcoFullDiffCut);
0056     tree_out_par -> Branch("INTT_vtxZ_QA", &INTT_vtxZ_QA);
0057     tree_out_par -> Branch("isWithRotate", &isWithRotate);
0058
0059     tree_out_par -> Branch("centrality_edges", &centrality_edges);
0060     tree_out_par -> Branch("nCentrality_bin", &nCentrality_bin);
0061
0062     tree_out_par -> Branch("CentralityFineEdge_min", &CentralityFineEdge_min);
0063     tree_out_par -> Branch("CentralityFineEdge_max", &CentralityFineEdge_max);
0064     tree_out_par -> Branch("nCentralityFineBin", &nCentralityFineBin);
0065
0066     tree_out_par -> Branch("EtaEdge_min", &EtaEdge_min);
0067     tree_out_par -> Branch("EtaEdge_max", &EtaEdge_max);
0068     tree_out_par -> Branch("nEtaBin", &nEtaBin);
0069
0070     tree_out_par -> Branch("VtxZEdge_min", &VtxZEdge_min);
0071     tree_out_par -> Branch("VtxZEdge_max", &VtxZEdge_max);
0072     tree_out_par -> Branch("nVtxZBin", &nVtxZBin);
0073
0074     tree_out_par -> Branch("DeltaPhiEdge_min", &DeltaPhiEdge_min);
0075     tree_out_par -> Branch("DeltaPhiEdge_max", &DeltaPhiEdge_max);
0076     tree_out_par -> Branch("nDeltaPhiBin", &nDeltaPhiBin);
0077
0078     tree_out_par -> Branch("typeA_sensorZID", &typeA_sensorZID);
0079     tree_out_par -> Branch("cut_GlobalMBDvtxZ", &cut_GlobalMBDvtxZ);
0080     tree_out_par -> Branch("cut_vtxZDiff", &cut_vtxZDiff);
0081     tree_out_par -> Branch("cut_TrapezoidalFitWidth", &cut_TrapezoidalFitWidth);
0082     tree_out_par -> Branch("cut_TrapezoidalFWHM", &cut_TrapezoidalFWHM);
0083     tree_out_par -> Branch("cut_InttBcoFullDIff_next", &cut_InttBcoFullDIff_next);
0084     tree_out_par -> Branch("cut_bestPair_DeltaPhi", &cut_bestPair_DeltaPhi);
0085     tree_out_par -> Branch("cut_GoodProtoTracklet_DeltaPhi", &cut_GoodProtoTracklet_DeltaPhi);
0086 }
0087
0088 void TrackletHistogramFill::PrepareHistograms()
0089 {
0090     h1D_centrality_bin = new TH1D("h1D_centrality_bin","h1D_centrality_bin;Centrality [%];Entries",nCentrality_bin,&centrality_edges[0]); // note : the 0-5%
0091
0092     h1D_eta_template = new TH1D("h1D_eta_template", "h1D_eta_template", nEtaBin, EtaEdge_min, EtaEdge_max); // note : coarse
0093     h1D_vtxz_template = new TH1D("h1D_vtxz_template", "h1D_vtxz_template", nVtxZBin, VtxZEdge_min, VtxZEdge_max); // note : coarse
0094
0095     h1D_PairDeltaEta_inclusive = new TH1D("h1D_PairDeltaEta_inclusive", "h1D_PairDeltaEta_inclusive;Pair #eta;Entries", nDeltaEtaBin, DeltaEtaEdge_min, DeltaEtaEdge_max);
0096     h1D_PairDeltaPhi_inclusive = new TH1D("h1D_PairDeltaPhi_inclusive", "h1D_PairDeltaPhi_inclusive;Pair #eta;Entries", nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max);
0097
0098     std::vector<int> insert_check; insert_check.clear();
0099     bool isInserted = false;
0100
0101     // note : for 1D
0102     h1D_map.clear();
0103
0104     isInserted = h1D_map.insert( std::make_pair(
0105             "h1D_centrality",
0106             new TH1D("h1D_centrality", "h1D_centrality;Centrality [%];Entries", nCentralityFineBin, CentralityFineEdge_min, CentralityFineEdge_max)
0107         )
0108     ).second; insert_check.push_back(isInserted);
0109
0110     for (int vtxz_bin = 0; vtxz_bin < nVtxZBin; vtxz_bin++)
0111     {
0112         isInserted = h1D_map.insert( std::make_pair(
0113                 Form("h1D_centrality_InttVtxZ%d", vtxz_bin),
0114                 new TH1D(Form("h1D_centrality_InttVtxZ%d", vtxz_bin), Form("h1D_centrality_InttVtxZ%d;Centrality [%];Entries",vtxz_bin), nCentralityFineBin, CentralityFineEdge_min, CentralityFineEdge_max)
0115             )
0116         ).second; insert_check.push_back(isInserted);
0117
0118         isInserted = h1D_map.insert( std::make_pair(
0119                 Form("h1D_centrality_MBDVtxZ%d", vtxz_bin),
0120                 new TH1D(Form("h1D_centrality_MBDVtxZ%d", vtxz_bin), Form("h1D_centrality_MBDVtxZ%d;Centrality [%];Entries",vtxz_bin), nCentralityFineBin, CentralityFineEdge_min, CentralityFineEdge_max)
0121             )
0122         ).second; insert_check.push_back(isInserted);
0123     }
0124
0125     for (int eta_bin = 0; eta_bin < nEtaBin; eta_bin++)
0126     {
0127         for (int vtxz_bin = 0; vtxz_bin < nVtxZBin; vtxz_bin++)
0128         {
0129             for (int Mbin = 0; Mbin < nCentrality_bin; Mbin++)
0130             {
0131                 isInserted = h1D_map.insert( std::make_pair(
0132                         Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, vtxz_bin),
0133                         new TH1D(Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d;Pair #Delta#phi [radian];Entries", Mbin, eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0134                     )
0135                 ).second; insert_check.push_back(isInserted);
0136
0137                 if (isWithRotate){
0138                     isInserted = h1D_map.insert( std::make_pair(
0139                             Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, vtxz_bin),
0140                             new TH1D(Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated;Pair #Delta#phi [radian];Entries", Mbin, eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0141                         )
0142                     ).second; insert_check.push_back(isInserted);
0143                 }
0144
0145
0146                 isInserted = h1D_map.insert( std::make_pair(
0147                         Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, vtxz_bin),
0148                         new TH1D(Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, vtxz_bin), Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d;Pair #Delta#phi [radian] (type A);Entries", Mbin, eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0149                     )
0150                 ).second; insert_check.push_back(isInserted);
0151
0152                 if (isWithRotate){
0153                     isInserted = h1D_map.insert( std::make_pair(
0154                             Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, vtxz_bin),
0155                             new TH1D(Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, vtxz_bin), Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated;Pair #Delta#phi [radian] (type A);Entries", Mbin, eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0156                         )
0157                     ).second; insert_check.push_back(isInserted);
0158                 }
0159
0160
0161                 // note : for truth eta
0162                 if (runnumber == -1 && eta_bin == 0){
0163                     isInserted = h1D_map.insert( std::make_pair(
0164                             Form("h1D_TrueEta_Mbin%d_VtxZ%d", Mbin, vtxz_bin),
0165                             new TH1D(Form("h1D_TrueEta_Mbin%d_VtxZ%d", Mbin, vtxz_bin), Form("h1D_TrueEta_Mbin%d_VtxZ%d;PHG4Particle #eta;Entries", Mbin, vtxz_bin), nEtaBin, EtaEdge_min, EtaEdge_max)
0166                         )
0167                     ).second; insert_check.push_back(isInserted);
0168                 }
0169             } // note : end Mbin
0170
0171             // // note : for inclusive [0 - 100]
0172             // isInserted = h1D_map.insert( std::make_pair(
0173                 //     Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d", eta_bin, vtxz_bin),
0174                 //     new TH1D(Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d", eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d;Pair #Delta#phi [radian];Entries", eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0175                 // )
0176             // ).second; insert_check.push_back(isInserted);
0177
0178             // isInserted = h1D_map.insert( std::make_pair(
0179                 //     Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d_rotated", eta_bin, vtxz_bin),
0180                 //     new TH1D(Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d_rotated", eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Inclusive100_Eta%d_VtxZ%d_rotated;Pair #Delta#phi [radian];Entries", eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0181                 // )
0182             // ).second; insert_check.push_back(isInserted);
0183
0184
0185             // // note : for inclusive [0 - 70]
0186             // isInserted = h1D_map.insert( std::make_pair(
0187                 //     Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d", eta_bin, vtxz_bin),
0188                 //     new TH1D(Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d", eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d;Pair #Delta#phi [radian];Entries", eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0189                 // )
0190             // ).second; insert_check.push_back(isInserted);
0191
0192             // isInserted = h1D_map.insert( std::make_pair(
0193                 //     Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d_rotated", eta_bin, vtxz_bin),
0194                 //     new TH1D(Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d_rotated", eta_bin, vtxz_bin), Form("h1D_DeltaPhi_Inclusive70_Eta%d_VtxZ%d_rotated;Pair #Delta#phi [radian];Entries", eta_bin, vtxz_bin), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0195                 // )
0196             // ).second; insert_check.push_back(isInserted);
0197
0198
0199             // note : for truth eta
0200             // if (runnumber == -1 && eta_bin == 0){
0201             //     isInserted = h1D_map.insert( std::make_pair(
0202                 //         Form("h1D_TrueEta_Inclusive100_VtxZ%d", vtxz_bin),
0203                 //         new TH1D(Form("h1D_TrueEta_Inclusive100_VtxZ%d", vtxz_bin), Form("h1D_TrueEta_Inclusive100_VtxZ%d;PHG4Particle #eta;Entries", vtxz_bin), nEtaBin, EtaEdge_min, EtaEdge_max)
0204                 // )
0205             //     ).second; insert_check.push_back(isInserted);
0206
0207             //     isInserted = h1D_map.insert( std::make_pair(
0208                 //         Form("h1D_TrueEta_Inclusive70_VtxZ%d", vtxz_bin),
0209                 //         new TH1D(Form("h1D_TrueEta_Inclusive70_VtxZ%d", vtxz_bin), Form("h1D_TrueEta_Inclusive70_VtxZ%d;PHG4Particle #eta;Entries", vtxz_bin), nEtaBin, EtaEdge_min, EtaEdge_max)
0210                 // )
0211             //     ).second; insert_check.push_back(isInserted);
0212             // }
0213         }
0214     }
0215
0216     // note : for the best pairs
0217     isInserted = h1D_map.insert( std::make_pair(
0218             Form("h1D_BestPair_Inclusive70_DeltaEta"),
0219             new TH1D(Form("h1D_BestPair_Inclusive70_DeltaEta"), Form("h1D_BestPair_Inclusive70_DeltaEta;Pair #Delta#eta;Entries"), nDeltaEtaBin, DeltaEtaEdge_min, DeltaEtaEdge_max)
0220         )
0221     ).second; insert_check.push_back(isInserted);
0222
0223     isInserted = h1D_map.insert( std::make_pair(
0224             Form("h1D_BestPair_Inclusive70_DeltaPhi"),
0225             new TH1D(Form("h1D_BestPair_Inclusive70_DeltaPhi"), Form("h1D_BestPair_Inclusive70_DeltaPhi;Pair #Delta#phi [radian];Entries"), nDeltaPhiBin, DeltaPhiEdge_min, DeltaPhiEdge_max)
0226         )
0227     ).second; insert_check.push_back(isInserted);
0228
0229     isInserted = h1D_map.insert( std::make_pair(
0230             Form("h1D_BestPair_Inclusive70_ClusPhiSize"),
0231             new TH1D(Form("h1D_BestPair_Inclusive70_ClusPhiSize"), Form("h1D_BestPair_Inclusive70_ClusPhiSize;ClusPhiSize;Entries"), 128,0,128)
0232         )
0233     ).second; insert_check.push_back(isInserted);
0234
0235     isInserted = h1D_map.insert( std::make_pair(
0236             Form("h1D_BestPair_Inclusive70_ClusAdc"),
0237             new TH1D(Form("h1D_BestPair_Inclusive70_ClusAdc"), Form("h1D_BestPair_Inclusive70_ClusAdc;ClusAdc;Entries"), 200,0,18000)
0238         )
0239     ).second; insert_check.push_back(isInserted);
0240
0241
0242
0243     // note : for 2D
0244     h2D_map.clear();
0245
0246     // note : truth nPHG4Particle in Eta-VtxZ
0247     if (runnumber == -1){
0248         for (int Mbin = 0; Mbin < nCentrality_bin; Mbin++)
0249         {
0250             isInserted = h2D_map.insert( std::make_pair(
0251                     Form("h2D_TrueEtaVtxZ_Mbin%d", Mbin),
0252                     new TH2D(Form("h2D_TrueEtaVtxZ_Mbin%d", Mbin), Form("h2D_TrueEtaVtxZ_Mbin%d;PHG4Particle #eta;TruthPV_trig_z [cm]", Mbin), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0253                 )
0254             ).second; insert_check.push_back(isInserted);
0255
0256             isInserted = h2D_map.insert( std::make_pair(
0257                     Form("h2D_TrueEtaVtxZ_Mbin%d_FineBin", Mbin),
0258                     new TH2D(Form("h2D_TrueEtaVtxZ_Mbin%d_FineBin", Mbin), Form("h2D_TrueEtaVtxZ_Mbin%d_FineBin;PHG4Particle #eta;TruthPV_trig_z [cm]", Mbin), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0259                 )
0260             ).second; insert_check.push_back(isInserted);
0261         }
0262
0263         isInserted = h1D_map.insert( std::make_pair(
0264                 "debug_h1D_NHadron_Inclusive100",
0265                 new TH1D("debug_h1D_NHadron_Inclusive100", "debug_h1D_NHadron_Inclusive100;N Hadrons;Entries", 200, 0, 10000)
0266             )
0267         ).second; insert_check.push_back(isInserted);
0268
0269         isInserted = h1D_map.insert( std::make_pair(
0270                 "debug_h1D_NHadron_OneEtaBin_Inclusive100",
0271                 new TH1D("debug_h1D_NHadron_OneEtaBin_Inclusive100", "debug_h1D_NHadron_OneEtaBin_Inclusive100;N Hadrons;Entries", 200, 0, 400)
0272             )
0273         ).second; insert_check.push_back(isInserted);
0274
0275         // isInserted = h2D_map.insert( std::make_pair(
0276             //     Form("h2D_Inclusive100_TrueEtaVtxZ"),
0277             //     new TH2D(Form("h2D_Inclusive100_TrueEtaVtxZ"), Form("h2D_Inclusive100_TrueEtaVtxZ;PHG4Particle #eta;TruthPV_trig_z [cm]"), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0278             // )
0279         // ).second; insert_check.push_back(isInserted);
0280
0281         // isInserted = h2D_map.insert( std::make_pair(
0282             //     Form("h2D_Inclusive100_TrueEtaVtxZ_FineBin"),
0283             //     new TH2D(Form("h2D_Inclusive100_TrueEtaVtxZ_FineBin"), Form("h2D_Inclusive100_TrueEtaVtxZ_FineBin;PHG4Particle #eta;TruthPV_trig_z [cm]"), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0284             // )
0285         // ).second; insert_check.push_back(isInserted);
0286
0287         // isInserted = h2D_map.insert( std::make_pair(
0288             //     Form("h2D_Inclusive70_TrueEtaVtxZ"),
0289             //     new TH2D(Form("h2D_Inclusive70_TrueEtaVtxZ"), Form("h2D_Inclusive70_TrueEtaVtxZ;PHG4Particle #eta;TruthPV_trig_z [cm]"), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0290             // )
0291         // ).second; insert_check.push_back(isInserted);
0292
0293         // isInserted = h2D_map.insert( std::make_pair(
0294             //     Form("h2D_Inclusive70_TrueEtaVtxZ_FineBin"),
0295             //     new TH2D(Form("h2D_Inclusive70_TrueEtaVtxZ_FineBin"), Form("h2D_Inclusive70_TrueEtaVtxZ_FineBin;PHG4Particle #eta;TruthPV_trig_z [cm]"), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0296             // )
0297         // ).second; insert_check.push_back(isInserted);
0298
0299         isInserted = h2D_map.insert( std::make_pair(
0300                 Form("h2D_TrueEvtCount_vtxZCentrality"),
0301                 new TH2D(Form("h2D_TrueEvtCount_vtxZCentrality"), Form("h2D_TrueEvtCount_vtxZCentrality;TruthPV_trig_z [cm];Centrality [%]"), nVtxZBin, VtxZEdge_min, VtxZEdge_max, nCentrality_bin, &centrality_edges[0]) // note : 0 - 5%
0302             )
0303         ).second; insert_check.push_back(isInserted);
0304     }
0305
0306     // note : for good proto tracklet in Eta-VtxZ
0307     for (int Mbin = 0; Mbin < nCentrality_bin; Mbin++){
0308
0309         isInserted = h2D_map.insert( std::make_pair(
0310                 Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin),
0311                 new TH2D(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin), Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d;Pair #eta;INTT vtxZ [cm]", Mbin), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0312             )
0313         ).second; insert_check.push_back(isInserted);
0314
0315         isInserted = h2D_map.insert( std::make_pair(
0316                 Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin),
0317                 new TH2D(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin), Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d;Pair #eta;INTT vtxZ [cm]", Mbin), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0318             )
0319         ).second; insert_check.push_back(isInserted);
0320
0321         // note : normal fine bin
0322         isInserted = h2D_map.insert( std::make_pair(
0323                 Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin),
0324                 new TH2D(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin), Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin;Pair #eta;INTT vtxZ [cm]", Mbin), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0325             )
0326         ).second; insert_check.push_back(isInserted);
0327
0328         isInserted = h2D_map.insert( std::make_pair(
0329                 Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin),
0330                 new TH2D(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin), Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin;Pair #eta;INTT vtxZ [cm]", Mbin), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0331             )
0332         ).second; insert_check.push_back(isInserted);
0333
0334
0335         // note : rotated
0336         isInserted = h2D_map.insert( std::make_pair(
0337                 Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin),
0338                 new TH2D(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin), Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated;Pair #eta;INTT vtxZ [cm]", Mbin), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0339             )
0340         ).second; insert_check.push_back(isInserted);
0341
0342         isInserted = h2D_map.insert( std::make_pair(
0343                 Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin),
0344                 new TH2D(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin), Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated;Pair #eta;INTT vtxZ [cm]", Mbin), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0345             )
0346         ).second; insert_check.push_back(isInserted);
0347     }
0348
0349
0350     isInserted = h2D_map.insert( std::make_pair(
0351             Form("h2D_RecoEvtCount_vtxZCentrality"),
0352             new TH2D(Form("h2D_RecoEvtCount_vtxZCentrality"), Form("h2D_RecoEvtCount_vtxZCentrality;INTT vtxZ [cm];Centrality [%]"), nVtxZBin, VtxZEdge_min, VtxZEdge_max, nCentrality_bin, &centrality_edges[0]) // note : 0 - 5%
0353         )
0354     ).second; insert_check.push_back(isInserted);
0355
0356     // note : best pair in Eta-VtxZ
0357     isInserted = h2D_map.insert( std::make_pair(
0358             Form("h2D_Inclusive100_BestPairEtaVtxZ"),
0359             new TH2D(Form("h2D_Inclusive100_BestPairEtaVtxZ"), Form("h2D_Inclusive100_BestPairEtaVtxZ;Tracklet #eta;INTT vtxZ [cm]"), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0360         )
0361     ).second; insert_check.push_back(isInserted);
0362
0363     isInserted = h2D_map.insert( std::make_pair(
0364             Form("h2D_Inclusive100_BestPairEtaVtxZ_FineBin"),
0365             new TH2D(Form("h2D_Inclusive100_BestPairEtaVtxZ_FineBin"), Form("h2D_Inclusive100_BestPairEtaVtxZ_FineBin;Tracklet #eta;INTT vtxZ [cm]"), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0366         )
0367     ).second; insert_check.push_back(isInserted);
0368
0369     isInserted = h2D_map.insert( std::make_pair(
0370             Form("h2D_Inclusive70_BestPairEtaVtxZ"),
0371             new TH2D(Form("h2D_Inclusive70_BestPairEtaVtxZ"), Form("h2D_Inclusive70_BestPairEtaVtxZ;Tracklet #eta;INTT vtxZ [cm]"), nEtaBin, EtaEdge_min, EtaEdge_max, nVtxZBin, VtxZEdge_min, VtxZEdge_max)
0372         )
0373     ).second; insert_check.push_back(isInserted);
0374
0375     isInserted = h2D_map.insert( std::make_pair(
0376             Form("h2D_Inclusive70_BestPairEtaVtxZ_FineBin"),
0377             new TH2D(Form("h2D_Inclusive70_BestPairEtaVtxZ_FineBin"), Form("h2D_Inclusive70_BestPairEtaVtxZ_FineBin;Tracklet #eta;INTT vtxZ [cm]"), 400, EtaEdge_min, EtaEdge_max, 400, VtxZEdge_min, VtxZEdge_max)
0378         )
0379     ).second; insert_check.push_back(isInserted);
0380
0381
0382     // note : vtxZ - centrality
0383     isInserted = h2D_map.insert( std::make_pair(
0384             Form("h2D_InttVtxZ_Centrality"),
0385             new TH2D(Form("h2D_InttVtxZ_Centrality"), Form("h2D_InttVtxZ_Centrality;INTT vtxZ [cm];Centrality"), nVtxZBin, VtxZEdge_min, VtxZEdge_max, nCentralityFineBin, CentralityFineEdge_min, CentralityFineEdge_max)
0386         )
0387     ).second; insert_check.push_back(isInserted);
0388
0389     isInserted = h2D_map.insert( std::make_pair(
0390             Form("h2D_MBDVtxZ_Centrality"),
0391             new TH2D(Form("h2D_MBDVtxZ_Centrality"), Form("h2D_MBDVtxZ_Centrality;INTT vtxZ [cm];Centrality"), nVtxZBin, VtxZEdge_min, VtxZEdge_max, nCentralityFineBin, CentralityFineEdge_min, CentralityFineEdge_max)
0392         )
0393     ).second; insert_check.push_back(isInserted);
0394
0395
0396     for (int isInserted_tag : insert_check){
0397         if (isInserted_tag == 0){
0398             std::cout<<"Histogram insertion failed"<<std::endl;
0399             exit(1);
0400         }
0401     }
0402 }
0403
0404 void TrackletHistogramFill::EvtCleanUp()
0405 {
0406     Used_Clus_index_map.clear();
0407     Pair_DeltaPhi_vec.clear();
0408 }
0409
0410 void TrackletHistogramFill::FillHistogram(
0411     // note : MBD & centrality relevant
0412     float MBD_z_vtx,
0413     int is_min_bias,
0414     float MBD_centrality,
0415     float MBD_south_charge_sum,
0416     float MBD_north_charge_sum,
0417     // float MBD_charge_sum,
0418     // float MBD_charge_asymm,
0419     int InttBcoFullDiff_next,
0420
0421     // // note : trigger tag
0422     int MBDNSg2,
0423     // int MBDNSg2_vtxZ10cm,
0424     // int MBDNSg2_vtxZ30cm,
0425     // int MBDNSg2_vtxZ60cm,
0426
0427     // std::vector<float> *ClusX,
0428     // std::vector<float> *ClusY,
0429     // std::vector<float> *ClusZ,
0430     // std::vector<int> *ClusLayer,
0431     std::vector<unsigned char> *ClusLadderZId,
0432     // std::vector<unsigned char> *ClusLadderPhiId,
0433     std::vector<int> *ClusAdc,
0434     std::vector<float> *ClusPhiSize,
0435
0436     // note : INTT vertex Z
0437     double INTTvtxZ,
0438     double INTTvtxZError,
0439     // double NgroupTrapezoidal,
0440     // double NgroupCoarse,
0441     double TrapezoidalFitWidth,
0442     double TrapezoidalFWHM,
0443
0444     // note : the tracklet pair
0445     std::vector<pair_str> *evt_TrackletPair_vec,
0446     std::vector<pair_str> *evt_TrackletPairRotate_vec,
0447
0448     // note : MC
0449     // float TruthPV_trig_x,
0450     // float TruthPV_trig_y,
0451     float TruthPV_trig_z,
0452     int NTruthVtx,
0453     int NPrimaryG4P,
0454     // std::vector<float> *PrimaryG4P_Pt,
0455     std::vector<float> *PrimaryG4P_Eta,
0456     // std::vector<float> *PrimaryG4P_Phi,
0457     // std::vector<float> *PrimaryG4P_E,
0458     // std::vector<float> *PrimaryG4P_PID,
0459     std::vector<int> *PrimaryG4P_isChargeHadron
0460 )
0461 {
0462     EvtCleanUp();
0463
0464     // ==============================================================================================================================================================================================================
0465     // note : for data
0466     if (runnumber != -1 && BcoFullDiffCut && InttBcoFullDiff_next <= cut_InttBcoFullDIff_next) {return;}
0467     if (runnumber != -1 && MBDNSg2 != 1) {return;} // todo: assume MC no trigger
0468
0469     // note : for MC
0470     if (runnumber == -1 && NTruthVtx != 1) {return;}
0471
0472     // note : both data and MC
0473     if (MBD_z_vtx != MBD_z_vtx) {return;}
0474     if (MBD_centrality != MBD_centrality) {return;}
0475     if (MBD_centrality < 0 || MBD_centrality > 1) {return;}
0476     if (INTTvtxZ != INTTvtxZ) {return;}
0477     if (MBD_z_vtx < cut_GlobalMBDvtxZ.first || MBD_z_vtx > cut_GlobalMBDvtxZ.second) {return;} // todo: the hard cut 60 cm
0478     if (is_min_bias != 1) {return;}
0479     // ==============================================================================================================================================================================================================
0480     // todo : No INTT vtxZ reweighting for the h1D_centrality_bin, as the comparison was done in RestDist
0481     int Mbin = h1D_centrality_bin -> Fill(MBD_centrality);
0482     Mbin = (Mbin == -1) ? -1 : Mbin - 1;
0483     if (Mbin == -1) {
0484         std::cout << "Mbin == -1, MBD_centrality = " << MBD_centrality << std::endl;
0485         return;
0486     }
0487
0488     int INTTvtxz_bin = h1D_vtxz_template->Fill(INTTvtxZ);
0489     INTTvtxz_bin = (INTTvtxz_bin == -1) ? -1 : INTTvtxz_bin - 1;
0490
0491     int MBDVtxz_bin = h1D_vtxz_template->Fill(MBD_z_vtx);
0492     MBDVtxz_bin = (MBDVtxz_bin == -1) ? -1 : MBDVtxz_bin - 1;
0493
0494     int TruthVtxZ_bin = h1D_vtxz_template->Fill(TruthPV_trig_z);
0495     TruthVtxZ_bin = (TruthVtxZ_bin == -1) ? -1 : TruthVtxZ_bin - 1;
0496
0497     if (runnumber == -1){
0498
0499         int NHadrons = 0;
0500         int NHadrons_OneEtaBin = 0;
0501
0502         if (h2D_map.find("h2D_TrueEvtCount_vtxZCentrality") != h2D_map.end()){
0503             h2D_map["h2D_TrueEvtCount_vtxZCentrality"] -> Fill(TruthPV_trig_z, MBD_centrality);
0504         }
0505
0506         for (int true_i = 0; true_i < NPrimaryG4P; true_i++){
0507             if (PrimaryG4P_isChargeHadron->at(true_i) != 1) { continue; }
0508
0509             // todo : for debug
0510             if (TruthPV_trig_z >= -10 && TruthPV_trig_z < 10){
0511                 NHadrons++;
0512                 // todo : the selected bin
0513                 if (PrimaryG4P_Eta->at(true_i) >= -1.1 && PrimaryG4P_Eta->at(true_i) < -0.9){
0514                     NHadrons_OneEtaBin++;
0515                 }
0516             }
0517
0518             if (h1D_map.find(Form("h1D_TrueEta_Mbin%d_VtxZ%d", Mbin, TruthVtxZ_bin)) != h1D_map.end()){
0519                 h1D_map[Form("h1D_TrueEta_Mbin%d_VtxZ%d", Mbin, TruthVtxZ_bin)] -> Fill(PrimaryG4P_Eta->at(true_i));
0520             }
0521
0522             if (h2D_map.find(Form("h2D_TrueEtaVtxZ_Mbin%d", Mbin)) != h2D_map.end()){
0523                 h2D_map[Form("h2D_TrueEtaVtxZ_Mbin%d", Mbin)] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0524             }
0525
0526             if (h2D_map.find(Form("h2D_TrueEtaVtxZ_Mbin%d_FineBin", Mbin)) != h2D_map.end()){
0527                 h2D_map[Form("h2D_TrueEtaVtxZ_Mbin%d_FineBin", Mbin)] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0528             }
0529
0530             // if (h2D_map.find("h2D_Inclusive100_TrueEtaVtxZ") != h2D_map.end()){
0531             //     h2D_map["h2D_Inclusive100_TrueEtaVtxZ"] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0532             // }
0533
0534             // if (h2D_map.find("h2D_Inclusive100_TrueEtaVtxZ_FineBin") != h2D_map.end()){
0535             //     h2D_map["h2D_Inclusive100_TrueEtaVtxZ_FineBin"] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0536             // }
0537
0538             // if (h1D_map.find("h1D_TrueEta_Inclusive100_VtxZ%d", TruthVtxZ_bin) != h1D_map.end()){
0539             //     h1D_map["h1D_TrueEta_Inclusive100_VtxZ%d", TruthVtxZ_bin] -> Fill(PrimaryG4P_Eta->at(true_i));
0540             // }
0541
0542             // if (0 <= MBD_centrality && MBD_centrality < 0.7){
0543             //     if (h2D_map.find("h2D_Inclusive70_TrueEtaVtxZ") != h2D_map.end()){
0544             //         h2D_map["h2D_Inclusive70_TrueEtaVtxZ"] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0545             //     }
0546
0547             //     if (h2D_map.find("h2D_Inclusive70_TrueEtaVtxZ_FineBin") != h2D_map.end()){
0548             //         h2D_map["h2D_Inclusive70_TrueEtaVtxZ_FineBin"] -> Fill(PrimaryG4P_Eta->at(true_i), TruthPV_trig_z);
0549             //     }
0550
0551             //     if (h1D_map.find("h1D_TrueEta_Inclusive70_VtxZ%d", TruthVtxZ_bin) != h1D_map.end()){
0552             //         h1D_map["h1D_TrueEta_Inclusive70_VtxZ%d", TruthVtxZ_bin] -> Fill(PrimaryG4P_Eta->at(true_i));
0553             //     }
0554
0555             // }
0556
0557         } // note : end of G4Particle loop
0558
0559         // todo : for debug
0560         if (TruthPV_trig_z >= -10 && TruthPV_trig_z < 10){
0561
0562             h1D_map["debug_h1D_NHadron_Inclusive100"] -> Fill(NHadrons);
0563             h1D_map["debug_h1D_NHadron_OneEtaBin_Inclusive100"] -> Fill(NHadrons_OneEtaBin);
0564         }
0565
0566
0567     } // note : end of truth
0568
0569     // note : FineBin
0570     // todo : No INTT vtxZ reweighting as this is the MBD_vtxZ filled
0571     if(h1D_map.find(Form("h1D_centrality_MBDVtxZ%d", MBDVtxz_bin)) != h1D_map.end()){
0572         h1D_map[Form("h1D_centrality_MBDVtxZ%d", MBDVtxz_bin)] -> Fill(MBD_centrality);
0573     }
0574     h2D_map["h2D_MBDVtxZ_Centrality"] -> Fill(MBD_z_vtx, MBD_centrality);
0575
0576
0577     // ==============================================================================================================================================================================================================
0578     if (INTT_vtxZ_QA && (MBD_z_vtx - INTTvtxZ < cut_vtxZDiff.first || MBD_z_vtx - INTTvtxZ > cut_vtxZDiff.second) ) {return;}
0579     if (INTT_vtxZ_QA && (TrapezoidalFitWidth < cut_TrapezoidalFitWidth.first || TrapezoidalFitWidth > cut_TrapezoidalFitWidth.second)){return;}
0580     if (INTT_vtxZ_QA && (TrapezoidalFWHM < cut_TrapezoidalFWHM.first || TrapezoidalFWHM > cut_TrapezoidalFWHM.second)){return;}
0581     if (INTT_vtxZ_QA && (INTTvtxZError < cut_INTTvtxZError.first || INTTvtxZError > cut_INTTvtxZError.second)){return;}
0582     // ==============================================================================================================================================================================================================
0583
0584     if (vtxZReweight && runnumber != -1){
0585         std::cout<<"Should have no vtxZ reweighting for data"<<std::endl;
0586         exit(1);
0587     }
0588
0589     if (vtxZReweight && h1D_vtxZReweightFactor == nullptr) {
0590         std::cout<<"vtxZReweightFactor histogram is not found"<<std::endl;
0591         exit(1);
0592     }
0593     vtxZReweightFactor = (vtxZReweight) ? h1D_vtxZReweightFactor -> GetBinContent(h1D_vtxZReweightFactor -> FindBin(INTTvtxZ)) : 1;
0594     vtxZReweightFactor = (vtxZReweightFactor > 250) ? 0 : vtxZReweightFactor; // todo: the rejection of the large weight
0595
0596     h1D_map["h1D_centrality"] -> Fill(MBD_centrality, vtxZReweightFactor);
0597
0598     if (h2D_map.find("h2D_InttVtxZ_Centrality") != h2D_map.end()){
0599         h2D_map["h2D_InttVtxZ_Centrality"] -> Fill(INTTvtxZ, MBD_centrality, vtxZReweightFactor);
0600     }
0601
0602     if (h1D_map.find(Form("h1D_centrality_InttVtxZ%d",INTTvtxz_bin)) != h1D_map.end()){
0603         h1D_map[Form("h1D_centrality_InttVtxZ%d",INTTvtxz_bin)] -> Fill(MBD_centrality, vtxZReweightFactor);
0604     }
0605
0606     // note : the event counting
0607     if (h2D_map.find("h2D_RecoEvtCount_vtxZCentrality") != h2D_map.end()){
0608         h2D_map["h2D_RecoEvtCount_vtxZCentrality"] -> Fill(INTTvtxZ, MBD_centrality, vtxZReweightFactor);
0609     }
0610
0611     std::vector<pair_str> p_evt_TrackletPair_vec = (evt_TrackletPair_vec != nullptr) ? *(evt_TrackletPair_vec) : std::vector<pair_str>(0);
0612     std::vector<pair_str> p_evt_TrackletPairRotate_vec = (evt_TrackletPairRotate_vec != nullptr) ? *(evt_TrackletPairRotate_vec) : std::vector<pair_str>(0);
0613
0614     // std::cout<<"test, evt_TrackletPair_vec->size(): "<<evt_TrackletPair_vec->size()<<" p_evt_TrackletPair_vec.size(): "<<p_evt_TrackletPair_vec.size()<<std::endl;
0615
0616     for (pair_str pair : p_evt_TrackletPair_vec)
0617     {
0618         int eta_bin = h1D_eta_template -> Fill(pair.pair_eta_fit);
0619         eta_bin = (eta_bin == -1) ? -1 : eta_bin - 1;
0620
0621         int inner_adc = ClusAdc->at(pair.inner_index);
0622         int outer_adc = ClusAdc->at(pair.outer_index);
0623         int inner_phi_size = ClusPhiSize->at(pair.inner_index);
0624         int outer_phi_size = ClusPhiSize->at(pair.outer_index);
0625
0626         if (isClusQA.first && (inner_adc <= isClusQA.second.first || outer_adc <= isClusQA.second.first)) {continue;} // note : adc
0627         if (isClusQA.first && (inner_phi_size > isClusQA.second.second || outer_phi_size > isClusQA.second.second)) {continue;} // note : phi size
0628
0629         // Pair_DeltaR_vec.push_back(sqrt(pow(pair.delta_phi,2)+pow(pair.delta_eta,2)));
0630         Pair_DeltaPhi_vec.push_back(fabs(pair.delta_phi));
0631
0632         h1D_PairDeltaEta_inclusive -> Fill(pair.delta_eta, vtxZReweightFactor);
0633         h1D_PairDeltaPhi_inclusive -> Fill(pair.delta_phi, vtxZReweightFactor);
0634
0635         // std::cout<<"Mbin: "<<Mbin<<" eta_bin: "<<eta_bin<<" INTTvtxz_bin: "<<INTTvtxz_bin<<std::endl;
0636
0637         if(h1D_map.find(Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, INTTvtxz_bin)) != h1D_map.end()){
0638             h1D_map[Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, INTTvtxz_bin)] -> Fill(pair.delta_phi, vtxZReweightFactor);
0639             // std::cout<<"In the if(), "<<Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, INTTvtxz_bin)<<std::endl;
0640         }
0641
0642         if (pair.delta_phi >= cut_GoodProtoTracklet_DeltaPhi.first && pair.delta_phi <= cut_GoodProtoTracklet_DeltaPhi.second){
0643
0644             if(h2D_map.find(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin)) != h2D_map.end()){
0645                 h2D_map[Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0646             }
0647
0648             if(h2D_map.find(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin)) != h2D_map.end()){
0649                 h2D_map[Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0650             }
0651         }
0652
0653         if (ClusLadderZId->at(pair.inner_index) != typeA_sensorZID[0] && ClusLadderZId->at(pair.inner_index) != typeA_sensorZID[1]) {continue;}
0654         if (ClusLadderZId->at(pair.outer_index) != typeA_sensorZID[0] && ClusLadderZId->at(pair.outer_index) != typeA_sensorZID[1]) {continue;}
0655
0656         if(h1D_map.find(Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, INTTvtxz_bin)) != h1D_map.end()){
0657             h1D_map[Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d", Mbin, eta_bin, INTTvtxz_bin)] -> Fill(pair.delta_phi, vtxZReweightFactor);
0658         }
0659
0660         if (pair.delta_phi >= cut_GoodProtoTracklet_DeltaPhi.first && pair.delta_phi <= cut_GoodProtoTracklet_DeltaPhi.second){
0661
0662             if(h2D_map.find(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin)) != h2D_map.end()){
0663                 h2D_map[Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0664             }
0665
0666             if(h2D_map.find(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin)) != h2D_map.end()){
0667                 h2D_map[Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_FineBin", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0668             }
0669         }
0670     }
0671
0672     if (isWithRotate){
0673         for (pair_str pair : p_evt_TrackletPairRotate_vec)
0674         {
0675             int eta_bin = h1D_eta_template -> Fill(pair.pair_eta_fit);
0676             eta_bin = (eta_bin == -1) ? -1 : eta_bin - 1;
0677
0678             int inner_adc = ClusAdc->at(pair.inner_index);
0679             int outer_adc = ClusAdc->at(pair.outer_index);
0680             int inner_phi_size = ClusPhiSize->at(pair.inner_index);
0681             int outer_phi_size = ClusPhiSize->at(pair.outer_index);
0682
0683             if (isClusQA.first && (inner_adc <= isClusQA.second.first || outer_adc <= isClusQA.second.first)) {continue;} // note : adc
0684             if (isClusQA.first && (inner_phi_size > isClusQA.second.second || outer_phi_size > isClusQA.second.second)) {continue;} // note : phi size
0685
0686             if(h1D_map.find(Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, INTTvtxz_bin)) != h1D_map.end()){
0687                 h1D_map[Form("h1D_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, INTTvtxz_bin)] -> Fill(pair.delta_phi, vtxZReweightFactor);
0688             }
0689
0690             if (pair.delta_phi >= cut_GoodProtoTracklet_DeltaPhi.first && pair.delta_phi <= cut_GoodProtoTracklet_DeltaPhi.second){
0691
0692                 if(h2D_map.find(Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin)) != h2D_map.end()){
0693                     h2D_map[Form("h2D_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0694                 }
0695             }
0696
0697             if (ClusLadderZId->at(pair.inner_index) != typeA_sensorZID[0] && ClusLadderZId->at(pair.inner_index) != typeA_sensorZID[1]) {continue;}
0698             if (ClusLadderZId->at(pair.outer_index) != typeA_sensorZID[0] && ClusLadderZId->at(pair.outer_index) != typeA_sensorZID[1]) {continue;}
0699
0700             if(h1D_map.find(Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, INTTvtxz_bin)) != h1D_map.end()){
0701                 h1D_map[Form("h1D_typeA_DeltaPhi_Mbin%d_Eta%d_VtxZ%d_rotated", Mbin, eta_bin, INTTvtxz_bin)] -> Fill(pair.delta_phi, vtxZReweightFactor);
0702             }
0703
0704             if (pair.delta_phi >= cut_GoodProtoTracklet_DeltaPhi.first && pair.delta_phi <= cut_GoodProtoTracklet_DeltaPhi.second){
0705
0706                 if(h2D_map.find(Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin)) != h2D_map.end()){
0707                     h2D_map[Form("h2D_typeA_GoodProtoTracklet_EtaVtxZ_Mbin%d_rotated", Mbin)] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0708                 }
0709             }
0710         }
0711     }
0712
0713     // if (MBD_centrality < 0 || MBD_centrality > 0.7) { return; } // note : 0 - 70% centrality
0714
0715     long long vec_size = Pair_DeltaPhi_vec.size();
0716     long long ind[Pair_DeltaPhi_vec.size()];
0717     TMath::Sort(vec_size, &Pair_DeltaPhi_vec[0], ind, false);
0718     for (int pair_i = 0; pair_i < vec_size; pair_i++)
0719     {
0720         pair_str pair = p_evt_TrackletPair_vec[ind[pair_i]];
0721         if (Used_Clus_index_map.find(pair.inner_index) != Used_Clus_index_map.end() || Used_Clus_index_map.find(pair.outer_index) != Used_Clus_index_map.end()) { continue; }
0722
0723         if (Pair_DeltaPhi_vec[ind[pair_i]] > cut_bestPair_DeltaPhi.second) { break; }
0724
0725         h2D_map["h2D_Inclusive100_BestPairEtaVtxZ"] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0726         h2D_map["h2D_Inclusive100_BestPairEtaVtxZ_FineBin"] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0727
0728
0729         if (MBD_centrality >= 0 && MBD_centrality <= 0.7){
0730             h2D_map["h2D_Inclusive70_BestPairEtaVtxZ"] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0731             h2D_map["h2D_Inclusive70_BestPairEtaVtxZ_FineBin"] -> Fill(pair.pair_eta_fit, INTTvtxZ, vtxZReweightFactor);
0732
0733             h1D_map["h1D_BestPair_Inclusive70_DeltaEta"] -> Fill(pair.delta_eta, vtxZReweightFactor);
0734             h1D_map["h1D_BestPair_Inclusive70_DeltaPhi"] -> Fill(pair.delta_phi, vtxZReweightFactor);
0735             h1D_map["h1D_BestPair_Inclusive70_ClusPhiSize"] -> Fill(ClusPhiSize->at(pair.inner_index), vtxZReweightFactor);
0736             h1D_map["h1D_BestPair_Inclusive70_ClusPhiSize"] -> Fill(ClusPhiSize->at(pair.outer_index), vtxZReweightFactor);
0737             h1D_map["h1D_BestPair_Inclusive70_ClusAdc"] -> Fill(ClusAdc->at(pair.inner_index), vtxZReweightFactor);
0738             h1D_map["h1D_BestPair_Inclusive70_ClusAdc"] -> Fill(ClusAdc->at(pair.outer_index), vtxZReweightFactor);
0739         }
0740
0741         Used_Clus_index_map.insert(std::make_pair(pair.inner_index, 1));
0742         Used_Clus_index_map.insert(std::make_pair(pair.outer_index, 1));
0743     }
0744
0745 }
0746
0747 void TrackletHistogramFill::EndRun()
0748 {
0749     file_out -> cd();
0750
0751     tree_out_par -> Fill();
0752     tree_out_par -> Write();
0753
0754     if (h1D_vtxZReweightFactor != nullptr){
0755         h1D_vtxZReweightFactor -> Write("h1D_vtxZReweightFactor");
0756     }
0757     h1D_eta_template -> Write("h1D_eta_template");
0758     h1D_vtxz_template -> Write("h1D_vtxz_template");
0759
0760     h1D_PairDeltaEta_inclusive -> Write("h1D_PairDeltaEta_inclusive");
0761     h1D_PairDeltaPhi_inclusive -> Write("h1D_PairDeltaPhi_inclusive");
0762
0763     h1D_centrality_bin -> Write();
0764
0765     for (auto &pair : h2D_map){
0766         pair.second -> Write();
0767     }
0768
0769     for (auto &pair : h1D_map){
0770         pair.second -> Write();
0771     }
0772
0773
0774
0775     file_out -> Close();
0776 }