Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:11:21 

0001 #include <TCanvas.h>
0002 #include <TCut.h>
0003 #include <TF1.h>
0004 #include <TFile.h>
0005 #include <TH1F.h>
0006 #include <TH2F.h>
0007 #include <TMath.h>
0008 #include <TObjString.h>
0009 #include <TRandom3.h>
0010 #include <TTree.h>
0011 #include <TTreeIndex.h>
0012 
0013 #include <algorithm>
0014 #include <fstream>
0015 #include <iostream>
0016 #include <iterator>
0017 #include <sstream>
0018 #include <string>
0019 #include <vector>
0020 
0021 #include "GenHadron.h"
0022 #include "Tracklet.h"
0023 #include "Vertex.h"
0024 
0025 // float theta2pseudorapidity(float theta) { return -1. * TMath::Log(TMath::Tan(theta / 2)); }
0026 
0027 // TF1 *ClusADCCut(int constscale, float etascale)
0028 // {
0029 //     TF1 *f = new TF1("f", Form("%d*TMath::CosH(%f*x)", constscale, etascale), -10, 10);
0030 //     return f;
0031 // }
0032 
0033 // TH1F *ClusADCCut_StepFunc_INTTPrivate()
0034 // {
0035 //     std::vector<float> thetastep = {0.001, 15, 20, 25, 30, 35, 45, 55, 125, 135, 145, 150, 155, 160, 165, 179.999};
0036 //     // reverse the thetastep
0037 //     std::reverse(thetastep.begin(), thetastep.end());
0038 //     std::vector<float> adccut_theta = {225, 165, 135, 120, 105, 90, 75, 60, 75, 90, 105, 120, 135, 165, 225};
0039 //     float etastep_array[thetastep.size()];
0040 //     cout << "etastep (INTT private): ";
0041 //     for (int i = 0; i < thetastep.size(); i++)
0042 //     {
0043 //         etastep_array[i] = theta2pseudorapidity(thetastep[i] * TMath::Pi() / 180);
0044 //         cout << etastep_array[i] << " ";
0045 //     }
0046 //     cout << endl;
0047 
0048 //     TH1F *hm_cut_inttprivate = new TH1F("hm_cut_inttprivate", "hm_cut_inttprivate", thetastep.size() - 1, etastep_array);
0049 //     cout << "ADC cut value: ";
0050 //     for (int j = 0; j < hm_cut_inttprivate->GetNbinsX(); j++)
0051 //     {
0052 //         hm_cut_inttprivate->SetBinContent(j + 1, adccut_theta[j]);
0053 //         cout << hm_cut_inttprivate->GetBinContent(j + 1) << " ";
0054 //     }
0055 //     cout << endl;
0056 
0057 //     return hm_cut_inttprivate;
0058 // }
0059 
0060 // TH1F *ClusADCCut_StepFunc(int constscale, float etascale)
0061 // {
0062 //     TF1 *f_cut = ClusADCCut(constscale, etascale);
0063 
0064 //     std::vector<float> adcstep;
0065 //     for (int i = 0; i < 20; i++)
0066 //     {
0067 //         adcstep.push_back(20 + i * 30);
0068 //     }
0069 //     std::vector<float> etastep;
0070 //     for (int i = 0; i < 20; i++)
0071 //     {
0072 //         etastep.insert(etastep.begin(), f_cut->GetX(adcstep[i], -10, 0));
0073 //         etastep.push_back(f_cut->GetX(adcstep[i], 0, 10));
0074 //     }
0075 //     // Remove nan values
0076 //     etastep.erase(std::remove_if(etastep.begin(), etastep.end(), [](float x) { return std::isnan(x); }), etastep.end());
0077 
0078 //     // print out cut information
0079 //     std::cout << "constscale = " << constscale << std::endl;
0080 //     // array of float from etastep
0081 //     float etastep_array[etastep.size()];
0082 //     cout << "etastep: ";
0083 //     for (int i = 0; i < etastep.size(); i++)
0084 //     {
0085 //         cout << etastep[i] << " ";
0086 //         etastep_array[i] = etastep[i];
0087 //     }
0088 //     cout << endl;
0089 //     TH1F *hm_cut = new TH1F("hm_cut", "hm_cut", etastep.size() - 1, etastep_array);
0090 //     cout << "ADC cut value: ";
0091 //     for (int j = 0; j < hm_cut->GetNbinsX(); j++)
0092 //     {
0093 //         if (hm_cut->GetBinLowEdge(j + 1) < 0)
0094 //         {
0095 //             hm_cut->SetBinContent(j + 1, f_cut->Eval(hm_cut->GetBinLowEdge(j + 1)));
0096 //         }
0097 //         else
0098 //         {
0099 //             hm_cut->SetBinContent(j + 1, f_cut->Eval(hm_cut->GetBinCenter(j + 1) + hm_cut->GetBinWidth(j + 1) / 2));
0100 //         }
0101 //         cout << hm_cut->GetBinContent(j + 1) << " ";
0102 //     }
0103 //     cout << endl;
0104 
0105 //     return hm_cut;
0106 
0107 //     // // construct the histogram
0108 //     // int nbins = 40, hmin = -10, hmax = 10;
0109 //     // TH1F *h = new TH1F("h", "h", nbins, hmin, hmax);
0110 //     // for (int i = 0; i < nbins; i++)
0111 //     // {
0112 //     //     h->SetBinContent(i + 1, f_cut->Eval(h->GetBinCenter(i + 1)));
0113 //     // }
0114 
0115 //     // // print the histogram
0116 //     // for (int i = 0; i < nbins; i++)
0117 //     // {
0118 //     //     std::cout << "h->GetBinContent(" << i + 1 << ") = " << h->GetBinContent(i + 1) << std::endl;
0119 //     // }
0120 
0121 //     // return h;
0122 
0123 //     // adcstep = [20+i*30 for i in range(0, 20)] # quantized step
0124 //     // etastep = []
0125 //     // f = TF1('f', '{}*TMath::CosH(x)'.format(scale), -10, 10)
0126 //     // for i in range(0, 20):
0127 //     //     etastep.insert(0, f.GetX(adcstep[i], -10, 0))
0128 //     //     etastep.append(f.GetX(adcstep[i], 0, 10))
0129 
0130 //     // # remove nan values
0131 //     // etastep = [x for x in etastep if not math.isnan(x)]
0132 //     // print(etastep)
0133 
0134 //     // hm_cut_v2 = TH1F('hm_cut_v2', 'hm_cut_v2', len(etastep)-1, array('d', etastep))
0135 //     // for j in range(hm_cut_v2.GetNbinsX()):
0136 //     //     # Set bin content at the quantized step
0137 //     //     if hm_cut_v2.GetBinLowEdge(j+1) < 0:
0138 //     //         hm_cut_v2.SetBinContent(j+1, f.Eval(hm_cut_v2.GetBinLowEdge(j+1)))
0139 //     //     else:
0140 //     //         hm_cut_v2.SetBinContent(j+1, f.Eval(hm_cut_v2.GetBinCenter(j+1) + hm_cut_v2.GetBinWidth(j+1)/2))
0141 //     //     # hm_cut_v2.SetBinContent(j+1, f.Eval(hm_cut_v2.GetBinCenter(j+1)))
0142 // }
0143 
0144 double cluszbin[56] = {-25.,
0145                        -22.57245 - 1.1, //
0146                        -22.57245 - 0.9,
0147                        -22.57245 + 0.9, //
0148                        -20.57245 - 0.9,
0149                        -20.57245 + 0.9, //
0150                        -18.57245 - 0.9,
0151                        -18.57245 + 0.9, //
0152                        -16.57245 - 0.9,
0153                        -16.57245 + 0.9, //
0154                        -14.57245 - 0.9,
0155                        -14.57245 + 0.9, //
0156                        -12.57245 - 0.7,
0157                        -12.57245 + 0.7, //
0158                        -10.97245 - 0.7,
0159                        -10.97245 + 0.7, //
0160                        -9.372450 - 0.7,
0161                        -9.372450 + 0.7, //
0162                        -7.772450 - 0.7,
0163                        -7.772450 + 0.7, //
0164                        -6.172450 - 0.7,
0165                        -6.172450 + 0.7, //
0166                        -4.572450 - 0.7,
0167                        -4.572450 + 0.7, //
0168                        -2.972450 - 0.7,
0169                        -2.972450 + 0.7, //
0170                        -1.372450 - 0.7,
0171                        -1.372450 + 0.7, //
0172                        0.4275496 - 0.7,
0173                        0.4275496 + 0.7, //
0174                        2.0275495 - 0.7,
0175                        2.0275495 + 0.7, //
0176                        3.6275494 - 0.7,
0177                        3.6275494 + 0.7, //
0178                        5.2275495 - 0.7,
0179                        5.2275495 + 0.7, //
0180                        6.8275494 - 0.7,
0181                        6.8275494 + 0.7, //
0182                        8.4275493 - 0.7,
0183                        8.4275493 + 0.7, //
0184                        10.027549 - 0.7,
0185                        10.027549 + 0.7, //
0186                        11.627549 - 0.7,
0187                        11.627549 + 0.7, //
0188                        13.627549 - 0.9,
0189                        13.627549 + 0.9, //
0190                        15.627549 - 0.9,
0191                        15.627549 + 0.9, //
0192                        17.627550 - 0.9,
0193                        17.627550 + 0.9, //
0194                        19.627550 - 0.9,
0195                        19.627550 + 0.9, //
0196                        21.627550 - 0.9,
0197                        21.627550 + 0.9, //
0198                        21.627550 + 1.1,
0199                        25.};
0200 
0201 int main(int argc, char *argv[])
0202 {
0203     if (argc != 5)
0204     {
0205         std::cout << "Usage: ./plotRecoCluster [isdata] [EvtVtx_map_filename] [infilename (ntuple)] [outfilename (histogram)]" << std::endl;
0206         return 0;
0207     }
0208 
0209     for (int i = 0; i < argc; i++)
0210     {
0211         std::cout << "argv[" << i << "] = " << argv[i] << std::endl;
0212     }
0213 
0214     bool IsData = (TString(argv[1]).Atoi() == 1) ? true : false;
0215     TString EvtVtx_map_filename = TString(argv[2]);
0216     TString infilename = TString(argv[3]);
0217     TString outfilename = TString(argv[4]);
0218 
0219     TString idxstr = (IsData) ? "INTT_BCO" : "event";
0220 
0221     // std::map<int, vector<float>> EvtVtx_map = EvtVtx_map_tklcluster(EvtVtx_map_filename.Data());
0222     std::map<int, vector<float>> EvtVtxMap_event = EvtVtx_map_event(EvtVtx_map_filename.Data());          // use with simulation
0223     std::map<uint64_t, vector<float>> EvtVtxMap_inttbco = EvtVtx_map_inttbco(EvtVtx_map_filename.Data()); // use with data
0224 
0225     // TH1F *hM_ClusX_all = new TH1F("hM_ClusX_all", "hM_ClusX_all", 200, -10, 10);
0226     // TH1F *hM_ClusX_layer1 = new TH1F("hM_ClusX_layer1", "hM_ClusX_layer1", 200, -10, 10);
0227     // TH1F *hM_ClusX_layer2 = new TH1F("hM_ClusX_layer2", "hM_ClusX_layer2", 200, -10, 10);
0228     // TH1F *hM_ClusY_all = new TH1F("hM_ClusY_all", "hM_ClusY_all", 200, -10, 10);
0229     // TH1F *hM_ClusY_layer1 = new TH1F("hM_ClusY_layer1", "hM_ClusY_layer1", 200, -10, 10);
0230     // TH1F *hM_ClusY_layer2 = new TH1F("hM_ClusY_layer2", "hM_ClusY_layer2", 200, -10, 10);
0231     TH1F *hM_ClusZ_all = new TH1F("hM_ClusZ_all", "hM_ClusZ_all", 55, cluszbin);
0232     TH1F *hM_ClusZ_layer1 = new TH1F("hM_ClusZ_layer1", "hM_ClusZ_layer1", 55, cluszbin);
0233     TH1F *hM_ClusZ_layer2 = new TH1F("hM_ClusZ_layer2", "hM_ClusZ_layer2", 55, cluszbin);
0234     // TH1F *hM_ClusR_all = new TH1F("hM_ClusR_all", "hM_ClusR_all", 120, 4, 7);
0235     // TH1F *hM_ClusR_layer1 = new TH1F("hM_ClusR_layer1", "hM_ClusR_layer1", 120, 4, 7);
0236     // TH1F *hM_ClusR_layer2 = new TH1F("hM_ClusR_layer2", "hM_ClusR_layer2", 120, 4, 7);
0237     // TH1F *hM_ClusEtaOri_all = new TH1F("hM_ClusEtaOri_all", "hM_ClusEtaOri_all", 160, -4, 4);
0238     // TH1F *hM_ClusEtaOri_layer1 = new TH1F("hM_ClusEtaOri_layer1", "hM_ClusEtaOri_layer1", 160, -4, 4);
0239     // TH1F *hM_ClusEtaOri_layer2 = new TH1F("hM_ClusEtaOri_layer2", "hM_ClusEtaOri_layer2", 160, -4, 4);
0240     TH1F *hM_ClusEtaPV_all = new TH1F("hM_ClusEtaPV_all", "hM_ClusEtaPV_all", 160, -4, 4);
0241     TH1F *hM_ClusEtaPV_layer1 = new TH1F("hM_ClusEtaPV_layer1", "hM_ClusEtaPV_layer1", 160, -4, 4);
0242     TH1F *hM_ClusEtaPV_layer2 = new TH1F("hM_ClusEtaPV_layer2", "hM_ClusEtaPV_layer2", 160, -4, 4);
0243     TH1F *hM_ClusEtaPV_all_ClusADCg35 = new TH1F("hM_ClusEtaPV_all_ClusADCg35", "hM_ClusEtaPV_all_ClusADCg35", 160, -4, 4);
0244     TH1F *hM_ClusEtaPV_layer1_ClusADCg35 = new TH1F("hM_ClusEtaPV_layer1_ClusADCg35", "hM_ClusEtaPV_layer1_ClusADCg35", 160, -4, 4);
0245     TH1F *hM_ClusEtaPV_layer2_ClusADCg35 = new TH1F("hM_ClusEtaPV_layer2_ClusADCg35", "hM_ClusEtaPV_layer2_ClusADCg35", 160, -4, 4);
0246     // remove clusters with cluster ADC/cluster phi size between [79,82] and [111,114]
0247     TH1F *hM_ClusEtaPV_all_ClusADCg35_ClusADCoverPhisizeCut = new TH1F("hM_ClusEtaPV_all_ClusADCg35_ClusADCoverPhisizeCut", "hM_ClusEtaPV_all_ClusADCg35_ClusADCoverPhisizeCut", 160, -4, 4);
0248     TH1F *hM_ClusEtaPV_layer1_ClusADCg35_ClusADCoverPhisizeCut = new TH1F("hM_ClusEtaPV_layer1_ClusADCg35_ClusADCoverPhisizeCut", "hM_ClusEtaPV_layer1_ClusADCg35_ClusADCoverPhisizeCut", 160, -4, 4);
0249     TH1F *hM_ClusEtaPV_layer2_ClusADCg35_ClusADCoverPhisizeCut = new TH1F("hM_ClusEtaPV_layer2_ClusADCg35_ClusADCoverPhisizeCut", "hM_ClusEtaPV_layer2_ClusADCg35_ClusADCoverPhisizeCut", 160, -4, 4);
0250     // cluster eta distribution for two strip sizes
0251     TH1F *hM_ClusEtaPV_1p6cmstrip_all = new TH1F("hM_ClusEtaPV_1p6cmstrip_all", "hM_ClusEtaPV_1p6cmstrip_all", 160, -4, 4); // 1.6 cm strip size; ClusLadderZId = 0 or 2
0252     TH1F *hM_ClusEtaPV_1p6cmstrip_layer1 = new TH1F("hM_ClusEtaPV_1p6cmstrip_layer1", "hM_ClusEtaPV_1p6cmstrip_layer1", 160, -4, 4);
0253     TH1F *hM_ClusEtaPV_1p6cmstrip_layer2 = new TH1F("hM_ClusEtaPV_1p6cmstrip_layer2", "hM_ClusEtaPV_1p6cmstrip_layer2", 160, -4, 4);
0254     TH1F *hM_ClusEtaPV_2cmstrip_all = new TH1F("hM_ClusEtaPV_2cmstrip_all", "hM_ClusEtaPV_2cmstrip_all", 160, -4, 4); // 2 cm strip size; ClusLadderZId = 1 or 3
0255     TH1F *hM_ClusEtaPV_2cmstrip_layer1 = new TH1F("hM_ClusEtaPV_2cmstrip_layer1", "hM_ClusEtaPV_2cmstrip_layer1", 160, -4, 4);
0256     TH1F *hM_ClusEtaPV_2cmstrip_layer2 = new TH1F("hM_ClusEtaPV_2cmstrip_layer2", "hM_ClusEtaPV_2cmstrip_layer2", 160, -4, 4);
0257     TH1F *hM_ClusEtaPV_ladderZID0_all = new TH1F("hM_ClusEtaPV_ladderZID0_all", "hM_ClusEtaPV_ladderZID0_all", 160, -4, 4); // ClusLadderZId = 0
0258     TH1F *hM_ClusEtaPV_ladderZID0_layer1 = new TH1F("hM_ClusEtaPV_ladderZID0_layer1", "hM_ClusEtaPV_ladderZID0_layer1", 160, -4, 4);
0259     TH1F *hM_ClusEtaPV_ladderZID0_layer2 = new TH1F("hM_ClusEtaPV_ladderZID0_layer2", "hM_ClusEtaPV_ladderZID0_layer2", 160, -4, 4);
0260     TH1F *hM_ClusEtaPV_ladderZID1_all = new TH1F("hM_ClusEtaPV_ladderZID1_all", "hM_ClusEtaPV_ladderZID1_all", 160, -4, 4); // ClusLadderZId = 1
0261     TH1F *hM_ClusEtaPV_ladderZID1_layer1 = new TH1F("hM_ClusEtaPV_ladderZID1_layer1", "hM_ClusEtaPV_ladderZID1_layer1", 160, -4, 4);
0262     TH1F *hM_ClusEtaPV_ladderZID1_layer2 = new TH1F("hM_ClusEtaPV_ladderZID1_layer2", "hM_ClusEtaPV_ladderZID1_layer2", 160, -4, 4);
0263     TH1F *hM_ClusEtaPV_ladderZID2_all = new TH1F("hM_ClusEtaPV_ladderZID2_all", "hM_ClusEtaPV_ladderZID2_all", 160, -4, 4); // ClusLadderZId = 2
0264     TH1F *hM_ClusEtaPV_ladderZID2_layer1 = new TH1F("hM_ClusEtaPV_ladderZID2_layer1", "hM_ClusEtaPV_ladderZID2_layer1", 160, -4, 4);
0265     TH1F *hM_ClusEtaPV_ladderZID2_layer2 = new TH1F("hM_ClusEtaPV_ladderZID2_layer2", "hM_ClusEtaPV_ladderZID2_layer2", 160, -4, 4);
0266     TH1F *hM_ClusEtaPV_ladderZID3_all = new TH1F("hM_ClusEtaPV_ladderZID3_all", "hM_ClusEtaPV_ladderZID3_all", 160, -4, 4); // ClusLadderZId = 3
0267     TH1F *hM_ClusEtaPV_ladderZID3_layer1 = new TH1F("hM_ClusEtaPV_ladderZID3_layer1", "hM_ClusEtaPV_ladderZID3_layer1", 160, -4, 4);
0268     TH1F *hM_ClusEtaPV_ladderZID3_layer2 = new TH1F("hM_ClusEtaPV_ladderZID3_layer2", "hM_ClusEtaPV_ladderZID3_layer2", 160, -4, 4);
0269 
0270     // TH1F *hM_ClusPhiOri_all = new TH1F("hM_ClusPhiOri_all", "hM_ClusPhiOri_all", 140, -3.5, 3.5);
0271     // TH1F *hM_ClusPhiOri_layer1 = new TH1F("hM_ClusPhiOri_layer1", "hM_ClusPhiOri_layer1", 140, -3.5, 3.5);
0272     // TH1F *hM_ClusPhiOri_layer2 = new TH1F("hM_ClusPhiOri_layer2", "hM_ClusPhiOri_layer2", 140, -3.5, 3.5);
0273     TH1F *hM_ClusPhiPV_all = new TH1F("hM_ClusPhiPV_all", "hM_ClusPhiPV_all", 140, -3.5, 3.5);
0274     TH1F *hM_ClusPhiPV_layer1 = new TH1F("hM_ClusPhiPV_layer1", "hM_ClusPhiPV_layer1", 140, -3.5, 3.5);
0275     TH1F *hM_ClusPhiPV_layer2 = new TH1F("hM_ClusPhiPV_layer2", "hM_ClusPhiPV_layer2", 140, -3.5, 3.5);
0276     TH1F *hM_ClusPhiPV_ClusPhiSize43or46 = new TH1F("hM_ClusPhiPV_ClusPhiSize43or46", "hM_ClusPhiPV_ClusPhiSize43or46", 330, -3.3, 3.3);
0277     // Very fine binning for ladder overlap study
0278     // INTT pitch size in phi = 78 um (Reference: https://indico.bnl.gov/event/15547/contributions/62868/attachments/41171/68968/2022_sPHENIX_School_RN.pdf)
0279     // INTT outer radius 102.62 mm
0280     // Azimuthal angle per pitch = (78*1E-6) / (102.62*1E-3) ~= 7.601E-4 rad ~= 0.0435 degree
0281     // Number of bins = (3.2*2) / (7.601E-4) = 8421
0282     TH1F *hM_ClusPhiPV_all_fine = new TH1F("hM_ClusPhiPV_all_fine", "hM_ClusPhiPV_all_fine", 8421, -3.2, 3.2);
0283     TH1F *hM_ClusPhiPV_layer1_fine = new TH1F("hM_ClusPhiPV_layer1_fine", "hM_ClusPhiPV_layer1_fine", 8421, -3.2, 3.2);
0284     TH1F *hM_ClusPhiPV_layer2_fine = new TH1F("hM_ClusPhiPV_layer2_fine", "hM_ClusPhiPV_layer2_fine", 8421, -3.2, 3.2);
0285     // cluster eta, weighted by the ratio of cluster ADC to the cluster phi size
0286     TH1F *hM_ClusEtaPV_all_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_all_weiClusADCoverPhiSize", "hM_ClusEtaPV_all_weiClusADCoverPhiSize", 160, -4, 4);
0287     TH1F *hM_ClusEtaPV_layer1_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_layer1_weiClusADCoverPhiSize", "hM_ClusEtaPV_layer1_weiClusADCoverPhiSize", 160, -4, 4);
0288     TH1F *hM_ClusEtaPV_layer2_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_layer2_weiClusADCoverPhiSize", "hM_ClusEtaPV_layer2_weiClusADCoverPhiSize", 160, -4, 4);
0289     TH1F *hM_ClusEtaPV_all_ClusADCg35_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_all_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_all_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4);
0290     TH1F *hM_ClusEtaPV_layer1_ClusADCg35_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_layer1_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_layer1_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4);
0291     TH1F *hM_ClusEtaPV_layer2_ClusADCg35_weiClusADCoverPhiSize = new TH1F("hM_ClusEtaPV_layer2_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_layer2_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4);
0292 
0293     TH2F *hM_ClusEtaPV_PVz_all_ClusADCg35 = new TH2F("hM_ClusEtaPV_PVz_all_ClusADCg35", "hM_ClusEtaPV_PVz_all_ClusADCg35", 160, -4, 4, 160, -40, 0);
0294     TH2F *hM_ClusEtaPV_PVz_layer1_ClusADCg35 = new TH2F("hM_ClusEtaPV_PVz_layer1_ClusADCg35", "hM_ClusEtaPV_PVz_layer1_ClusADCg35", 160, -4, 4, 160, -40, 0);
0295     TH2F *hM_ClusEtaPV_PVz_layer2_ClusADCg35 = new TH2F("hM_ClusEtaPV_PVz_layer2_ClusADCg35", "hM_ClusEtaPV_PVz_layer2_ClusADCg35", 160, -4, 4, 160, -40, 0);
0296     TH2F *hM_ClusEtaPV_PVz_all_LadderZID0 = new TH2F("hM_ClusEtaPV_PVz_all_LadderZID0", "hM_ClusEtaPV_PVz_all_LadderZID0", 160, -4, 4, 140, -35, 35);
0297     TH2F *hM_ClusEtaPV_PVz_all_LadderZID1 = new TH2F("hM_ClusEtaPV_PVz_all_LadderZID1", "hM_ClusEtaPV_PVz_all_LadderZID1", 160, -4, 4, 140, -35, 35);
0298     TH2F *hM_ClusEtaPV_PVz_all_LadderZID2 = new TH2F("hM_ClusEtaPV_PVz_all_LadderZID2", "hM_ClusEtaPV_PVz_all_LadderZID2", 160, -4, 4, 140, -35, 35);
0299     TH2F *hM_ClusEtaPV_PVz_all_LadderZID3 = new TH2F("hM_ClusEtaPV_PVz_all_LadderZID3", "hM_ClusEtaPV_PVz_all_LadderZID3", 160, -4, 4, 140, -35, 35);
0300     TH2F *hM_ClusEtaPV_PVz_all_ClusADCg35_weiClusADCoverPhiSize = new TH2F("hM_ClusEtaPV_PVz_all_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_PVz_all_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4, 160, -40, 0);
0301     TH2F *hM_ClusEtaPV_PVz_layer1_ClusADCg35_weiClusADCoverPhiSize = new TH2F("hM_ClusEtaPV_PVz_layer1_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_PVz_layer1_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4, 160, -40, 0);
0302     TH2F *hM_ClusEtaPV_PVz_layer2_ClusADCg35_weiClusADCoverPhiSize = new TH2F("hM_ClusEtaPV_PVz_layer2_ClusADCg35_weiClusADCoverPhiSize", "hM_ClusEtaPV_PVz_layer2_ClusADCg35_weiClusADCoverPhiSize", 160, -4, 4, 160, -40, 0);
0303 
0304     TH1F *hM_ClusADC_all = new TH1F("hM_ClusADC_all", "hM_ClusADC_all", 1800, 0, 18000);
0305     TH1F *hM_ClusADC_layer1 = new TH1F("hM_ClusADC_layer1", "hM_ClusADC_layer1", 1800, 0, 18000);
0306     TH1F *hM_ClusADC_layer2 = new TH1F("hM_ClusADC_layer2", "hM_ClusADC_layer2", 1800, 0, 18000);
0307     TH1F *hM_ClusADC_all_0to500 = new TH1F("hM_ClusADC_all_0to500", "hM_ClusADC_all_0to500", 500, 0, 500);
0308     TH1F *hM_ClusADC_layer1_0to500 = new TH1F("hM_ClusADC_layer1_0to500", "hM_ClusADC_layer1_0to500", 500, 0, 500);
0309     TH1F *hM_ClusADC_layer2_0to500 = new TH1F("hM_ClusADC_layer2_0to500", "hM_ClusADC_layer2_0to500", 500, 0, 500);
0310     // Select perpendicularly incident clusters, i.e cluster |eta|<=0.1
0311     TH1F *hM_ClusADC_all_etale0p1 = new TH1F("hM_ClusADC_all_etale0p1", "hM_ClusADC_all_etale0p1", 1800, 0, 18000);
0312     TH1F *hM_ClusADC_layer1_etale0p1 = new TH1F("hM_ClusADC_layer1_etale0p1", "hM_ClusADC_layer1_etale0p1", 1800, 0, 18000);
0313     TH1F *hM_ClusADC_layer2_etale0p1 = new TH1F("hM_ClusADC_layer2_etale0p1", "hM_ClusADC_layer2_etale0p1", 1800, 0, 18000);
0314     TH1F *hM_ClusADC_all_etale0p1_ClusZSize1 = new TH1F("hM_ClusADC_all_etale0p1_ClusZSize1", "hM_ClusADC_all_etale0p1_ClusZSize1", 1800, 0, 18000);
0315     TH1F *hM_ClusADC_layer1_etale0p1_ClusZSize1 = new TH1F("hM_ClusADC_layer1_etale0p1_ClusZSize1", "hM_ClusADC_layer1_etale0p1_ClusZSize1", 1800, 0, 18000);
0316     TH1F *hM_ClusADC_layer2_etale0p1_ClusZSize1 = new TH1F("hM_ClusADC_layer2_etale0p1_ClusZSize1", "hM_ClusADC_layer2_etale0p1_ClusZSize1", 1800, 0, 18000);
0317     TH1F *hM_ClusZSize_all = new TH1F("hM_ClusZSize_all", "hM_ClusZSize_all", 10, 0, 10);
0318     TH1F *hM_ClusZSize_layer1 = new TH1F("hM_ClusZSize_layer1", "hM_ClusZSize_layer1", 10, 0, 10);
0319     TH1F *hM_ClusZSize_layer2 = new TH1F("hM_ClusZSize_layer2", "hM_ClusZSize_layer2", 10, 0, 10);
0320     TH1F *hM_ClusZSize_all_etale0p1 = new TH1F("hM_ClusZSize_all_etale0p1", "hM_ClusZSize_all_etale0p1", 10, 0, 10);
0321     TH1F *hM_ClusZSize_layer1_etale0p1 = new TH1F("hM_ClusZSize_layer1_etale0p1", "hM_ClusZSize_layer1_etale0p1", 10, 0, 10);
0322     TH1F *hM_ClusZSize_layer2_etale0p1 = new TH1F("hM_ClusZSize_layer2_etale0p1", "hM_ClusZSize_layer2_etale0p1", 10, 0, 10);
0323     TH1F *hM_ClusPhiSize_all = new TH1F("hM_ClusPhiSize_all", "hM_ClusPhiSize_all", 130, 0, 130);
0324     TH1F *hM_ClusPhiSize_layer1 = new TH1F("hM_ClusPhiSize_layer1", "hM_ClusPhiSize_layer1", 130, 0, 130);
0325     TH1F *hM_ClusPhiSize_layer2 = new TH1F("hM_ClusPhiSize_layer2", "hM_ClusPhiSize_layer2", 130, 0, 130);
0326     // Select perpendicularly incident clusters, i.e cluster |eta|<=0.1
0327     TH1F *hM_ClusPhiSize_all_etale0p1 = new TH1F("hM_ClusPhiSize_all_etale0p1", "hM_ClusPhiSize_all_etale0p1", 130, 0, 130);
0328     TH1F *hM_ClusPhiSize_layer1_etale0p1 = new TH1F("hM_ClusPhiSize_layer1_etale0p1", "hM_ClusPhiSize_layer1_etale0p1", 130, 0, 130);
0329     TH1F *hM_ClusPhiSize_layer2_etale0p1 = new TH1F("hM_ClusPhiSize_layer2_etale0p1", "hM_ClusPhiSize_layer2_etale0p1", 130, 0, 130);
0330     // Select perpendicularly incident clusters with cluster Z size = 1
0331     TH1F *hM_ClusPhiSize_all_etale0p1_ClusZSize1 = new TH1F("hM_ClusPhiSize_all_etale0p1_ClusZSize1", "hM_ClusPhiSize_all_etale0p1_ClusZSize1", 130, 0, 130);
0332     TH1F *hM_ClusPhiSize_layer1_etale0p1_ClusZSize1 = new TH1F("hM_ClusPhiSize_layer1_etale0p1_ClusZSize1", "hM_ClusPhiSize_layer1_etale0p1_ClusZSize1", 130, 0, 130);
0333     TH1F *hM_ClusPhiSize_layer2_etale0p1_ClusZSize1 = new TH1F("hM_ClusPhiSize_layer2_etale0p1_ClusZSize1", "hM_ClusPhiSize_layer2_etale0p1_ClusZSize1", 130, 0, 130);
0334     TH2F *hM_ClusX_ClusY_all = new TH2F("hM_ClusX_ClusY_all", "hM_ClusX_ClusY_all", 260, -13, 13, 260, -13, 13);
0335     TH2F *hM_ClusX_ClusY_all_weiphisize = new TH2F("hM_ClusX_ClusY_all_weiphisize", "hM_ClusX_ClusY_all_weiphisize", 260, -13, 13, 260, -13, 13);
0336     TH2F *hM_ClusX_ClusY_all_weiclusadc = new TH2F("hM_ClusX_ClusY_all_weiclusadc", "hM_ClusX_ClusY_all_weiclusadc", 260, -13, 13, 260, -13, 13);
0337     TH2F *hM_ClusX_ClusY_ClusPhiSize43or46 = new TH2F("hM_ClusX_ClusY_ClusPhiSize43or46", "hM_ClusX_ClusY_ClusPhiSize43or46", 260, -13, 13, 260, -13, 13);
0338     TH2F *hM_ClusX_ClusY_ClusPhiSize43 = new TH2F("hM_ClusX_ClusY_ClusPhiSize43", "hM_ClusX_ClusY_ClusPhiSize43", 260, -13, 13, 260, -13, 13);
0339     TH2F *hM_ClusX_ClusY_ClusPhiSize46 = new TH2F("hM_ClusX_ClusY_ClusPhiSize46", "hM_ClusX_ClusY_ClusPhiSize46", 260, -13, 13, 260, -13, 13);
0340 
0341     TH2F *hM_ClusZ_ClusPhiPV_all = new TH2F("hM_ClusZ_ClusPhiPV_all", "hM_ClusZ_ClusPhiPV_all", 1000, -25, 25, 350, -3.5, 3.5);
0342     TH2F *hM_ClusZ_ClusPhiPV_layer1 = new TH2F("hM_ClusZ_ClusPhiPV_layer1", "hM_ClusZ_ClusPhiPV_layer1", 1000, -25, 25, 350, -3.5, 3.5);
0343     TH2F *hM_ClusZ_ClusPhiPV_layer2 = new TH2F("hM_ClusZ_ClusPhiPV_layer2", "hM_ClusZ_ClusPhiPV_layer2", 1000, -25, 25, 350, -3.5, 3.5);
0344     TH2F *hM_ClusZ_ClusPhiPV_all_coarse = new TH2F("hM_ClusZ_ClusPhiPV_all_coarse", "hM_ClusZ_ClusPhiPV_all_coarse", 55, cluszbin, 350, -3.5, 3.5);
0345     TH2F *hM_ClusZ_ClusPhiPV_all_coarse_weiphisize = new TH2F("hM_ClusZ_ClusPhiPV_all_coarse_weiphisize", "hM_ClusZ_ClusPhiPV_all_coarse_weiphisize", 55, cluszbin, 350, -3.5, 3.5);
0346     TH2F *hM_ClusZ_ClusPhiPV_all_coarse_weiclusadc = new TH2F("hM_ClusZ_ClusPhiPV_all_coarse_weiclusadc", "hM_ClusZ_ClusPhiPV_all_coarse_weiclusadc", 55, cluszbin, 350, -3.5, 3.5);
0347     // TH2F *hM_ClusEta_ClusZSize_all = new TH2F("hM_ClusEta_ClusZSize_all", "hM_ClusEta_ClusZSize_all", 160, -4, 4, 20, 0, 20);
0348     // TH2F *hM_ClusEta_ClusZSize_layer1 = new TH2F("hM_ClusEta_ClusZSize_layer1", "hM_ClusEta_ClusZSize_layer1", 160, -4, 4, 20, 0, 20);
0349     // TH2F *hM_ClusEta_ClusZSize_layer2 = new TH2F("hM_ClusEta_ClusZSize_layer2", "hM_ClusEta_ClusZSize_layer2", 160, -4, 4, 20, 0, 20);
0350     TH2F *hM_ClusPhiPV_ClusPhiSize_all = new TH2F("hM_ClusPhiPV_ClusPhiSize_all", "hM_ClusPhiPV_ClusPhiSize_all", 140, -3.5, 3.5, 100, 0, 100);
0351     TH2F *hM_ClusPhiPV_ClusPhiSize_layer1 = new TH2F("hM_ClusPhiPV_ClusPhiSize_layer1", "hM_ClusPhiPV_ClusPhiSize_layer1", 140, -3.5, 3.5, 100, 0, 100);
0352     TH2F *hM_ClusPhiPV_ClusPhiSize_layer2 = new TH2F("hM_ClusPhiPV_ClusPhiSize_layer2", "hM_ClusPhiPV_ClusPhiSize_layer2", 140, -3.5, 3.5, 100, 0, 100);
0353     TH2F *hM_ClusPhiPV_ClusADC_all = new TH2F("hM_ClusPhiPV_ClusADC_all", "hM_ClusPhiPV_ClusADC_all", 140, -3.5, 3.5, 180, 0, 18000);
0354     TH2F *hM_ClusPhiPV_ClusADC_layer1 = new TH2F("hM_ClusPhiPV_ClusADC_layer1", "hM_ClusPhiPV_ClusADC_layer1", 140, -3.5, 3.5, 180, 0, 18000);
0355     TH2F *hM_ClusPhiPV_ClusADC_layer2 = new TH2F("hM_ClusPhiPV_ClusADC_layer2", "hM_ClusPhiPV_ClusADC_layer2", 140, -3.5, 3.5, 180, 0, 18000);
0356     // TH2F *hM_ClusZSize_ClusPhiSize_all = new TH2F("hM_ClusZSize_ClusPhiSize_all", "hM_ClusZSize_ClusPhiSize_all", 20, 0, 20, 100, 0, 100);
0357     // TH2F *hM_ClusZSize_ClusPhiSize_layer1 = new TH2F("hM_ClusZSize_ClusPhiSize_layer1", "hM_ClusZSize_ClusPhiSize_layer1", 20, 0, 20, 100, 0, 100);
0358     // TH2F *hM_ClusZSize_ClusPhiSize_layer2 = new TH2F("hM_ClusZSize_ClusPhiSize_layer2", "hM_ClusZSize_ClusPhiSize_layer2", 20, 0, 20, 100, 0, 100);
0359     TH2F *hM_ClusEtaPV_ClusADC_all = new TH2F("hM_ClusEtaPV_ClusADC_all", "hM_ClusEtaPV_ClusADC_all", 160, -4, 4, 180, 0, 18000);
0360     TH2F *hM_ClusEtaPV_ClusADC_layer1 = new TH2F("hM_ClusEtaPV_ClusADC_layer1", "hM_ClusEtaPV_ClusADC_layer1", 160, -4, 4, 180, 0, 18000);
0361     TH2F *hM_ClusEtaPV_ClusADC_layer2 = new TH2F("hM_ClusEtaPV_ClusADC_layer2", "hM_ClusEtaPV_ClusADC_layer2", 160, -4, 4, 180, 0, 18000);
0362     TH2F *hM_ClusEtaPV_ClusADC_all_zoomin = new TH2F("hM_ClusEtaPV_ClusADC_all_zoomin", "hM_ClusEtaPV_ClusADC_all_zoomin", 120, -3, 3, 51, -0.5, 250.5);
0363     TH2F *hM_ClusEtaPV_ClusADC_layer1_zoomin = new TH2F("hM_ClusEtaPV_ClusADC_layer1_zoomin", "hM_ClusEtaPV_ClusADC_layer1_zoomin", 120, -3, 3, 50, -0.5, 249.5);
0364     TH2F *hM_ClusEtaPV_ClusADC_layer2_zoomin = new TH2F("hM_ClusEtaPV_ClusADC_layer2_zoomin", "hM_ClusEtaPV_ClusADC_layer2_zoomin", 120, -3, 3, 50, -0.5, 249.5);
0365     TH2F *hM_ClusEtaPV_ClusPhiSize_all = new TH2F("hM_ClusEtaPV_ClusPhiSize_all", "hM_ClusEtaPV_ClusPhiSize_all", 160, -4, 4, 100, 0, 100);
0366     TH2F *hM_ClusEtaPV_ClusPhiSize_layer1 = new TH2F("hM_ClusEtaPV_ClusPhiSize_layer1", "hM_ClusEtaPV_ClusPhiSize_layer1", 160, -4, 4, 100, 0, 100);
0367     TH2F *hM_ClusEtaPV_ClusPhiSize_layer2 = new TH2F("hM_ClusEtaPV_ClusPhiSize_layer2", "hM_ClusEtaPV_ClusPhiSize_layer2", 160, -4, 4, 100, 0, 100);
0368     TH2F *hM_ClusEtaPV_ClusPhiSize_all_zoomin = new TH2F("hM_ClusEtaPV_ClusPhiSize_all_zoomin", "hM_ClusEtaPV_ClusPhiSize_all_zoomin", 120, -3, 3, 20, 0, 20);
0369     TH2F *hM_ClusEtaPV_ClusPhiSize_layer1_zoomin = new TH2F("hM_ClusEtaPV_ClusPhiSize_layer1_zoomin", "hM_ClusEtaPV_ClusPhiSize_layer1_zoomin", 120, -3, 3, 20, 0, 20);
0370     TH2F *hM_ClusEtaPV_ClusPhiSize_layer2_zoomin = new TH2F("hM_ClusEtaPV_ClusPhiSize_layer2_zoomin", "hM_ClusEtaPV_ClusPhiSize_layer2_zoomin", 120, -3, 3, 20, 0, 20);
0371     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_all_zoomin = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_all_zoomin", "hM_ClusEtaPV_ClusADCoverClusPhiSize_all_zoomin", 120, -3, 3, 250, 0, 250);
0372     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_zoomin = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_zoomin", "hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_zoomin", 120, -3, 3, 250, 0, 250);
0373     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_zoomin = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_zoomin", "hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_zoomin", 120, -3, 3, 250, 0, 250);
0374     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_all_ClusADCg35 = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_all_ClusADCg35", "hM_ClusEtaPV_ClusADCoverClusPhiSize_all_ClusADCg35", 120, -3, 3, 250, 0, 250);
0375     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_ClusADCg35 = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_ClusADCg35", "hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_ClusADCg35", 120, -3, 3, 250, 0, 250);
0376     TH2F *hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_ClusADCg35 = new TH2F("hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_ClusADCg35", "hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_ClusADCg35", 120, -3, 3, 250, 0, 250);
0377 
0378     TH1F *hM_mutualdRcluster_all = new TH1F("hM_mutualdRcluster_all", "hM_mutualdRcluster_all", 200, 0, 0.1);
0379     TH1F *hM_mutualdRcluster_layer1 = new TH1F("hM_mutualdRcluster_layer1", "hM_mutualdRcluster_layer1", 200, 0, 0.1);
0380     TH1F *hM_mutualdRcluster_layer2 = new TH1F("hM_mutualdRcluster_layer2", "hM_mutualdRcluster_layer2", 200, 0, 0.1);
0381 
0382     // For cluster ADC cut optimization
0383     int constscale_low = 5, constscale_high = 55, constscale_step = 1;
0384     vector<TH1F *> v_hM_ClusEtaPV_EtaDepADCCut;
0385     vector<TH2F *> v_hM_ClusEtaPV_ClusADC_all_zoomin;
0386     // For cluster ADC cut optimization, histogram for cut values for each constscale
0387     vector<TH1F *> v_hM_ClusADCCut;
0388     vector<TF1 *> v_f_ClusADCCut;
0389     for (int constscale = constscale_low; constscale <= constscale_high; constscale += constscale_step)
0390     {
0391         v_hM_ClusEtaPV_EtaDepADCCut.push_back(new TH1F(Form("hM_ClusEtaPV_EtaDepADCCut_constscale%d_etascale1p0", constscale), Form("hM_ClusEtaPV_EtaDepADCCut_constscale%d_etascale1p0", constscale), 140, -3.5, 3.5));
0392         v_hM_ClusEtaPV_ClusADC_all_zoomin.push_back(new TH2F(Form("hM_ClusEtaPV_ClusADC_all_zoomin_constscale%d_etascale1p0", constscale), Form("hM_ClusEtaPV_ClusADC_all_zoomin_constscale%d_etascale1p0", constscale), 120, -3, 3, 50, 0, 500));
0393         v_hM_ClusADCCut.push_back(ClusADCCut_StepFunc(constscale, 1));
0394         v_f_ClusADCCut.push_back(ClusADCCut(constscale, 1));
0395     }
0396 
0397     vector<TH1F *> v_hM_ClusEtaPV_EtaDepADCCut_etascale;
0398     vector<TH2F *> v_hM_ClusEtaPV_ClusADC_all_zoomin_etascale;
0399     vector<TH1F *> v_hM_ClusADCCut_etascale;
0400     vector<TF1 *> v_f_ClusADCCut_etascale;
0401     float etascale_low = 0.6, etascale_high = 1.2, etascale_step = 0.02;
0402     for (float etascale = etascale_low; etascale <= etascale_high; etascale += etascale_step)
0403     {
0404         TString etascale_str = Form("%.2f", etascale);
0405         etascale_str.ReplaceAll(".", "p");
0406         v_hM_ClusEtaPV_EtaDepADCCut_etascale.push_back(new TH1F(Form("hM_ClusEtaPV_EtaDepADCCut_constscale30_etascale%s", etascale_str.Data()), Form("hM_ClusEtaPV_EtaDepADCCut_constscale30_etascale%s", etascale_str.Data()), 140, -3.5, 3.5));
0407         v_hM_ClusEtaPV_ClusADC_all_zoomin_etascale.push_back(new TH2F(Form("hM_ClusEtaPV_ClusADC_all_zoomin_constscale30_etascale%s", etascale_str.Data()), Form("hM_ClusEtaPV_ClusADC_all_zoomin_constscale30_etascale%s", etascale_str.Data()), 120, -3, 3, 50, 0, 500));
0408         v_hM_ClusADCCut_etascale.push_back(ClusADCCut_StepFunc(30, etascale));
0409         v_f_ClusADCCut_etascale.push_back(ClusADCCut(30, etascale));
0410     }
0411     // for the cluster adc cut from INTT private study
0412     TH1F *hM_ClusADCCut_inttprivate = ClusADCCut_StepFunc_INTTPrivate();
0413     TH1F *hM_ClusEtaPV_EtaDepADCCut_inttprivate = new TH1F("hM_ClusEtaPV_EtaDepADCCut_inttprivate", "hM_ClusEtaPV_EtaDepADCCut_inttprivate", 140, -3.5, 3.5);
0414     TH2F *hM_ClusEtaPV_ClusADC_all_zoomin_inttprivate = new TH2F("hM_ClusEtaPV_ClusADC_all_zoomin_inttprivate", "hM_ClusEtaPV_ClusADC_all_zoomin_inttprivate", 120, -3, 3, 50, 0, 500);
0415 
0416     TH2F *hM_ClusPhiSize_ClusADC_all = new TH2F("hM_ClusPhiSize_ClusADC_all", "hM_ClusPhiSize_ClusADC_all", 100, 0, 100, 180, 0, 18000);
0417     TH2F *hM_ClusPhiSize_ClusADC_layer1 = new TH2F("hM_ClusPhiSize_ClusADC_layer1", "hM_ClusPhiSize_ClusADC_layer1", 100, 0, 100, 180, 0, 18000);
0418     TH2F *hM_ClusPhiSize_ClusADC_layer2 = new TH2F("hM_ClusPhiSize_ClusADC_layer2", "hM_ClusPhiSize_ClusADC_layer2", 100, 0, 100, 180, 0, 18000);
0419     TH2F *hM_ClusPhiSize_ClusADC_all_zoomin = new TH2F("hM_ClusPhiSize_ClusADC_all_zoomin", "hM_ClusPhiSize_ClusADC_all_zoomin", 20, 0, 20, 50, 0, 500);
0420     TH2F *hM_ClusPhiSize_ClusADC_layer1_zoomin = new TH2F("hM_ClusPhiSize_ClusADC_layer1_zoomin", "hM_ClusPhiSize_ClusADC_layer1_zoomin", 20, 0, 20, 50, 0, 500);
0421     TH2F *hM_ClusPhiSize_ClusADC_layer2_zoomin = new TH2F("hM_ClusPhiSize_ClusADC_layer2_zoomin", "hM_ClusPhiSize_ClusADC_layer2_zoomin", 20, 0, 20, 50, 0, 500);
0422 
0423     // Vertex Z reweighting
0424     TH1F *hM_vtxzweight = VtxZ_ReweiHist("/sphenix/user/hjheng/sPHENIXRepo/analysis/dNdEta_Run2023/analysis_INTT/plot/RecoPV_ana/HIJING_ananew_20250106/VtxZ_reweight_HIJING_ananew_20250106.root", "VtxZ_reweight_VtxZm10to10");
0425 
0426     TFile *f = new TFile(infilename, "READ");
0427     TTree *t = (TTree *)f->Get("EventTree");
0428     t->BuildIndex(idxstr); // Reference: https://root-forum.cern.ch/t/sort-ttree-entries/13138
0429     TTreeIndex *index = (TTreeIndex *)t->GetTreeIndex();
0430     int event;
0431     bool is_min_bias, InttBco_IsToBeRemoved;
0432     float MBD_z_vtx, MBD_south_charge_sum, MBD_north_charge_sum;
0433     vector<int> *firedTriggers = 0;
0434     uint64_t INTT_BCO;
0435     vector<int> *ClusLayer = 0;
0436     vector<int> *ClusLadderZId = 0;
0437     vector<float> *ClusX = 0, *ClusY = 0, *ClusZ = 0, *ClusR = 0, *ClusPhi = 0, *ClusEta = 0, *ClusPhiSize = 0, *ClusZSize = 0;
0438     vector<unsigned int> *ClusAdc = 0;
0439     t->SetBranchAddress("event", &event);
0440     t->SetBranchAddress("is_min_bias", &is_min_bias);
0441     if (t->GetListOfBranches()->FindObject("InttBco_IsToBeRemoved"))
0442     {
0443         t->SetBranchAddress("InttBco_IsToBeRemoved", &InttBco_IsToBeRemoved);
0444     }
0445     else
0446     {
0447         InttBco_IsToBeRemoved = false;
0448     }
0449     t->SetBranchAddress("MBD_z_vtx", &MBD_z_vtx);
0450     t->SetBranchAddress("MBD_south_charge_sum", &MBD_south_charge_sum);
0451     t->SetBranchAddress("MBD_north_charge_sum", &MBD_north_charge_sum);
0452     if (IsData)
0453     {
0454         t->SetBranchAddress("firedTriggers", &firedTriggers);
0455     }
0456 
0457     t->SetBranchAddress("INTT_BCO", &INTT_BCO);
0458     t->SetBranchAddress("ClusLayer", &ClusLayer);
0459     t->SetBranchAddress("ClusLadderZId", &ClusLadderZId);
0460     t->SetBranchAddress("ClusX", &ClusX);
0461     t->SetBranchAddress("ClusY", &ClusY);
0462     t->SetBranchAddress("ClusZ", &ClusZ);
0463     t->SetBranchAddress("ClusR", &ClusR);
0464     t->SetBranchAddress("ClusPhi", &ClusPhi);
0465     t->SetBranchAddress("ClusEta", &ClusEta);
0466     t->SetBranchAddress("ClusAdc", &ClusAdc);
0467     t->SetBranchAddress("ClusPhiSize", &ClusPhiSize);
0468     t->SetBranchAddress("ClusZSize", &ClusZSize);
0469     for (int ev = 0; ev < t->GetEntriesFast(); ev++)
0470     {
0471         Long64_t local = t->LoadTree(index->GetIndex()[ev]);
0472         t->GetEntry(local);
0473 
0474         if (InttBco_IsToBeRemoved)
0475             continue;
0476 
0477         vector<float> PV = (IsData) ? EvtVtxMap_inttbco[INTT_BCO] : EvtVtxMap_event[event];
0478 
0479         if (PV.size() != 3)
0480         {
0481             cout << "[ERROR] No PV found for event " << event << endl;
0482             exit(1);
0483         }
0484 
0485         bool validMbdVtx = (MBD_z_vtx == MBD_z_vtx);
0486         bool MbdNScharnge = (MBD_north_charge_sum > 0 || MBD_south_charge_sum > 0);
0487         bool firedTrig = (IsData) ? (find(firedTriggers->begin(), firedTriggers->end(), 10) != firedTriggers->end()) : true;
0488         bool MbdZvtxRange = (MBD_z_vtx > -10 && MBD_z_vtx < 10);
0489         bool InttZvtxRange = (PV[2] > -10 && PV[2] < 10);
0490         bool EvtSel = validMbdVtx && MbdNScharnge && firedTrig && MbdZvtxRange && InttZvtxRange;
0491         cout << "Event " << event << " validMbdVtx: " << validMbdVtx << " MbdNScharnge: " << MbdNScharnge << " firedTrig: " << firedTrig << " MbdZvtxRange: " << MbdZvtxRange << " EvtSel: " << EvtSel << endl;
0492 
0493         if (!EvtSel)
0494             continue;
0495 
0496         float vtxzwei = (IsData) ? 1. : hM_vtxzweight->GetBinContent(hM_vtxzweight->FindBin(PV[2]));
0497 
0498         cout << "Number of clusters: " << ClusLayer->size() << endl;
0499 
0500         std::vector<std::vector<Hit *>> hits_layer =  {{}, {}};
0501         for (size_t i = 0; i < hits_layer.size(); i++)
0502         {
0503             for (auto &hit : hits_layer[i])
0504             {
0505                 delete hit;
0506             }
0507             CleanVec(hits_layer[i]);
0508         }
0509 
0510         
0511 
0512         for (size_t i = 0; i < ClusLayer->size(); i++)
0513         {
0514             if (ClusLayer->at(i) < 3 || ClusLayer->at(i) > 6)
0515             {
0516                 cout << "[ERROR] Unknown layer: " << ClusLayer->at(i) << endl; // Should not happen
0517                 exit(1);
0518             }
0519 
0520             int layer = (ClusLayer->at(i) == 3 || ClusLayer->at(i) == 4) ? 0 : 1;
0521 
0522             Hit *hit = new Hit(ClusX->at(i), ClusY->at(i), ClusZ->at(i), PV[0], PV[1], PV[2], layer);
0523             hits_layer[layer].push_back(hit);
0524 
0525             // hM_ClusX_all->Fill(ClusX->at(i));
0526             // hM_ClusY_all->Fill(ClusY->at(i));
0527             hM_ClusZ_all->Fill(ClusZ->at(i), vtxzwei);
0528             // hM_ClusR_all->Fill(ClusR->at(i));
0529             // hM_ClusEtaOri_all->Fill(ClusEta->at(i));
0530             hM_ClusEtaPV_all->Fill(hit->Eta(), vtxzwei);
0531             hM_ClusEtaPV_all_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0532 
0533             // hM_ClusPhiOri_all->Fill(ClusPhi->at(i));
0534             hM_ClusPhiPV_all->Fill(hit->Phi(), vtxzwei);
0535             hM_ClusPhiPV_all_fine->Fill(hit->Phi(), vtxzwei);
0536             hM_ClusADC_all->Fill(ClusAdc->at(i), vtxzwei);
0537             hM_ClusPhiSize_all->Fill(ClusPhiSize->at(i), vtxzwei);
0538             hM_ClusZSize_all->Fill(ClusZSize->at(i), vtxzwei);
0539             hM_ClusX_ClusY_all->Fill(ClusX->at(i), ClusY->at(i), vtxzwei);
0540             if (ClusPhiSize->at(i) == 43 || ClusPhiSize->at(i) == 46)
0541             {
0542                 hM_ClusPhiPV_ClusPhiSize43or46->Fill(hit->Phi(), vtxzwei);
0543                 hM_ClusX_ClusY_ClusPhiSize43or46->Fill(ClusX->at(i), ClusY->at(i), vtxzwei);
0544 
0545                 if (ClusPhiSize->at(i) == 43)
0546                 {
0547                     hM_ClusX_ClusY_ClusPhiSize43->Fill(ClusX->at(i), ClusY->at(i), vtxzwei);
0548                 }
0549 
0550                 if (ClusPhiSize->at(i) == 46)
0551                 {
0552                     hM_ClusX_ClusY_ClusPhiSize46->Fill(ClusX->at(i), ClusY->at(i), vtxzwei);
0553                 }
0554             }
0555 
0556             // Select perpendicularly incident clusters, i.e cluster |eta|<=0.1
0557             if (fabs(hit->Eta()) <= 0.1)
0558             {
0559                 hM_ClusZSize_all_etale0p1->Fill(ClusZSize->at(i), vtxzwei);
0560                 hM_ClusPhiSize_all_etale0p1->Fill(ClusPhiSize->at(i), vtxzwei);
0561                 hM_ClusADC_all_etale0p1->Fill(ClusAdc->at(i), vtxzwei);
0562                 if (ClusZSize->at(i) == 1)
0563                 {
0564                     hM_ClusADC_all_etale0p1_ClusZSize1->Fill(ClusAdc->at(i), vtxzwei);
0565                     hM_ClusPhiSize_all_etale0p1_ClusZSize1->Fill(ClusPhiSize->at(i), vtxzwei);
0566                 }
0567             }
0568 
0569             hM_ClusX_ClusY_all_weiphisize->Fill(ClusX->at(i), ClusY->at(i), ClusPhiSize->at(i) * vtxzwei);
0570             hM_ClusX_ClusY_all_weiclusadc->Fill(ClusX->at(i), ClusY->at(i), ClusAdc->at(i) * vtxzwei);
0571             hM_ClusZ_ClusPhiPV_all_coarse->Fill(ClusZ->at(i), hit->Phi(), vtxzwei);
0572             hM_ClusZ_ClusPhiPV_all_coarse_weiphisize->Fill(ClusZ->at(i), hit->Phi(), ClusPhiSize->at(i) * vtxzwei);
0573             hM_ClusZ_ClusPhiPV_all_coarse_weiclusadc->Fill(ClusZ->at(i), hit->Phi(), ClusAdc->at(i) * vtxzwei);
0574             hM_ClusZ_ClusPhiPV_all->Fill(ClusZ->at(i), hit->Phi(), vtxzwei);
0575             hM_ClusPhiPV_ClusPhiSize_all->Fill(hit->Phi(), ClusPhiSize->at(i), vtxzwei);
0576             hM_ClusEtaPV_ClusADC_all->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0577             hM_ClusEtaPV_ClusADC_all_zoomin->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0578             hM_ClusEtaPV_ClusPhiSize_all->Fill(hit->Eta(), ClusPhiSize->at(i), vtxzwei);
0579             hM_ClusEtaPV_ClusPhiSize_all_zoomin->Fill(hit->Eta(), ClusPhiSize->at(i), vtxzwei);
0580             hM_ClusPhiPV_ClusADC_all->Fill(hit->Phi(), ClusAdc->at(i), vtxzwei);
0581             hM_ClusPhiSize_ClusADC_all->Fill(ClusPhiSize->at(i), ClusAdc->at(i), vtxzwei);
0582             hM_ClusPhiSize_ClusADC_all_zoomin->Fill(ClusPhiSize->at(i), ClusAdc->at(i), vtxzwei);
0583             hM_ClusEtaPV_ClusADCoverClusPhiSize_all_zoomin->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0584             if (ClusAdc->at(i) > 35)
0585             {
0586                 hM_ClusEtaPV_all_ClusADCg35->Fill(hit->Eta(), vtxzwei);
0587                 hM_ClusEtaPV_all_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0588                 hM_ClusEtaPV_ClusADCoverClusPhiSize_all_ClusADCg35->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0589                 hM_ClusEtaPV_PVz_all_ClusADCg35->Fill(hit->Eta(), PV[2]);
0590                 hM_ClusEtaPV_PVz_all_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), PV[2], ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0591 
0592                 if (ClusLadderZId->at(i) == 0 || ClusLadderZId->at(i) == 2)
0593                 {
0594                     hM_ClusEtaPV_1p6cmstrip_all->Fill(hit->Eta(), vtxzwei);
0595                     if (ClusLadderZId->at(i) == 0)
0596                     {
0597                         hM_ClusEtaPV_ladderZID0_all->Fill(hit->Eta(), vtxzwei);
0598                         hM_ClusEtaPV_PVz_all_LadderZID0->Fill(hit->Eta(), PV[2], vtxzwei);
0599                     }
0600                     if (ClusLadderZId->at(i) == 2)
0601                     {
0602                         hM_ClusEtaPV_ladderZID2_all->Fill(hit->Eta(), vtxzwei);
0603                         hM_ClusEtaPV_PVz_all_LadderZID2->Fill(hit->Eta(), PV[2], vtxzwei);
0604                     }
0605                 }
0606                 if (ClusLadderZId->at(i) == 1 || ClusLadderZId->at(i) == 3)
0607                 {
0608                     hM_ClusEtaPV_2cmstrip_all->Fill(hit->Eta(), vtxzwei);
0609                     if (ClusLadderZId->at(i) == 1)
0610                     {
0611                         hM_ClusEtaPV_ladderZID1_all->Fill(hit->Eta(), vtxzwei);
0612                         hM_ClusEtaPV_PVz_all_LadderZID1->Fill(hit->Eta(), PV[2], vtxzwei);
0613                     }
0614                     if (ClusLadderZId->at(i) == 3)
0615                     {
0616                         hM_ClusEtaPV_ladderZID3_all->Fill(hit->Eta(), vtxzwei);
0617                         hM_ClusEtaPV_PVz_all_LadderZID3->Fill(hit->Eta(), PV[2], vtxzwei);
0618                     }
0619                 }
0620 
0621                 // cluster adc / phi size cut
0622                 // remove clusters with cluster ADC/cluster phi size between [79,82] and [111,114]
0623                 float clusadcoverphisize = ClusAdc->at(i) / ClusPhiSize->at(i);
0624                 bool Pass_adcoverphisize = !(clusadcoverphisize >= 79 && clusadcoverphisize <= 82) && !(clusadcoverphisize >= 111 && clusadcoverphisize <= 114);
0625                 if (Pass_adcoverphisize)
0626                 {
0627                     hM_ClusEtaPV_all_ClusADCg35_ClusADCoverPhisizeCut->Fill(hit->Eta(), vtxzwei);
0628                 }
0629             }
0630             // hM_ClusZSize_all->Fill(ClusZSize->at(i));
0631             // hM_ClusZSize_ClusPhiSize_all->Fill(ClusZSize->at(i), ClusPhiSize->at(i));
0632             // hM_ClusZSize_ClusADC_all->Fill(ClusZSize->at(i), ClusAdc->at(i));
0633 
0634             // apply the eta-dependent cluster ADC
0635             for (int constscale = constscale_low; constscale <= constscale_high; constscale += constscale_step)
0636             {
0637                 // Using TH1F (i.e step function)
0638                 if (ClusAdc->at(i) > v_hM_ClusADCCut[(constscale - constscale_low) / constscale_step]->GetBinContent(v_hM_ClusADCCut[(constscale - constscale_low) / constscale_step]->FindBin(hit->Eta())))
0639                 // Using TF1
0640                 // if (ClusAdc->at(i) > v_f_ClusADCCut[(constscale - constscale_low) / constscale_step]->Eval(hit->Eta()))
0641                 {
0642                     v_hM_ClusEtaPV_ClusADC_all_zoomin[(constscale - constscale_low) / constscale_step]->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0643                     v_hM_ClusEtaPV_EtaDepADCCut[(constscale - constscale_low) / constscale_step]->Fill(hit->Eta(), vtxzwei);
0644                 }
0645             }
0646 
0647             for (float etascale = etascale_low; etascale <= etascale_high; etascale += etascale_step)
0648             {
0649                 // Using TH1F (i.e step function)
0650                 if (ClusAdc->at(i) > v_hM_ClusADCCut_etascale[(int)((etascale - etascale_low) / etascale_step)]->GetBinContent(v_hM_ClusADCCut_etascale[(int)((etascale - etascale_low) / etascale_step)]->FindBin(hit->Eta())))
0651                 // if (ClusAdc->at(i) > v_f_ClusADCCut_etascale[(int)((etascale - etascale_low) / etascale_step)]->Eval(hit->Eta()))
0652                 {
0653                     v_hM_ClusEtaPV_ClusADC_all_zoomin_etascale[(int)((etascale - etascale_low) / etascale_step)]->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0654                     v_hM_ClusEtaPV_EtaDepADCCut_etascale[(int)((etascale - etascale_low) / etascale_step)]->Fill(hit->Eta(), vtxzwei);
0655                 }
0656             }
0657 
0658             // for intt pricate cluster adc cut
0659             if (ClusAdc->at(i) > hM_ClusADCCut_inttprivate->GetBinContent(hM_ClusADCCut_inttprivate->FindBin(hit->Eta())))
0660             {
0661                 hM_ClusEtaPV_ClusADC_all_zoomin_inttprivate->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0662                 hM_ClusEtaPV_EtaDepADCCut_inttprivate->Fill(hit->Eta(), vtxzwei);
0663             }
0664 
0665             if (layer == 0)
0666             {
0667 
0668                 // hM_ClusX_layer1->Fill(ClusX->at(i));
0669                 // hM_ClusY_layer1->Fill(ClusY->at(i));
0670                 hM_ClusZ_layer1->Fill(ClusZ->at(i), vtxzwei);
0671                 // hM_ClusR_layer1->Fill(ClusR->at(i));
0672                 // hM_ClusEtaOri_layer1->Fill(ClusEta->at(i));
0673                 hM_ClusEtaPV_layer1->Fill(hit->Eta(), vtxzwei);
0674                 hM_ClusEtaPV_layer1_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0675 
0676                 // hM_ClusPhiOri_layer1->Fill(ClusPhi->at(i));
0677                 hM_ClusPhiPV_layer1->Fill(hit->Phi(), vtxzwei);
0678                 hM_ClusPhiPV_layer1_fine->Fill(hit->Phi(), vtxzwei);
0679                 hM_ClusADC_layer1->Fill(ClusAdc->at(i), vtxzwei);
0680                 hM_ClusZSize_layer1->Fill(ClusZSize->at(i), vtxzwei);
0681                 hM_ClusPhiSize_layer1->Fill(ClusPhiSize->at(i), vtxzwei);
0682                 hM_ClusZ_ClusPhiPV_layer1->Fill(ClusZ->at(i), hit->Phi(), vtxzwei);
0683                 hM_ClusPhiPV_ClusPhiSize_layer1->Fill(hit->Phi(), ClusPhiSize->at(i), vtxzwei);
0684                 hM_ClusEtaPV_ClusADC_layer1->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0685                 hM_ClusEtaPV_ClusADC_layer1_zoomin->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0686                 hM_ClusEtaPV_ClusPhiSize_layer1->Fill(hit->Eta(), ClusPhiSize->at(i)), vtxzwei;
0687                 hM_ClusEtaPV_ClusPhiSize_layer1_zoomin->Fill(hit->Eta(), ClusPhiSize->at(i), vtxzwei);
0688                 hM_ClusPhiPV_ClusADC_layer1->Fill(hit->Phi(), ClusAdc->at(i), vtxzwei);
0689                 hM_ClusPhiSize_ClusADC_layer1->Fill(ClusPhiSize->at(i), ClusAdc->at(i), vtxzwei);
0690                 hM_ClusPhiSize_ClusADC_layer1_zoomin->Fill(ClusPhiSize->at(i), ClusAdc->at(i)), vtxzwei;
0691                 hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_zoomin->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0692 
0693                 // Select perpendicularly incident clusters, i.e cluster |eta|<=0.1
0694                 if (fabs(hit->Eta()) <= 0.1)
0695                 {
0696                     hM_ClusZSize_layer1_etale0p1->Fill(ClusZSize->at(i), vtxzwei);
0697                     hM_ClusPhiSize_layer1_etale0p1->Fill(ClusPhiSize->at(i), vtxzwei);
0698                     hM_ClusADC_layer1_etale0p1->Fill(ClusAdc->at(i), vtxzwei);
0699                     if (ClusZSize->at(i) == 1)
0700                     {
0701                         hM_ClusADC_layer1_etale0p1_ClusZSize1->Fill(ClusAdc->at(i), vtxzwei);
0702                         hM_ClusPhiSize_layer1_etale0p1_ClusZSize1->Fill(ClusPhiSize->at(i), vtxzwei);
0703                     }
0704                 }
0705 
0706                 if (ClusAdc->at(i) > 35)
0707                 {
0708                     hM_ClusEtaPV_layer1_ClusADCg35->Fill(hit->Eta(), vtxzwei);
0709                     hM_ClusEtaPV_layer1_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0710                     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_ClusADCg35->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0711                     hM_ClusEtaPV_PVz_layer1_ClusADCg35->Fill(hit->Eta(), PV[2]);
0712                     hM_ClusEtaPV_PVz_layer1_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), PV[2], ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0713 
0714                     if (ClusLadderZId->at(i) == 0 || ClusLadderZId->at(i) == 2)
0715                     {
0716                         hM_ClusEtaPV_1p6cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0717                         if (ClusLadderZId->at(i) == 0)
0718                         {
0719                             hM_ClusEtaPV_ladderZID0_layer1->Fill(hit->Eta(), vtxzwei);
0720                             hM_ClusEtaPV_1p6cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0721                         }
0722                         if (ClusLadderZId->at(i) == 2)
0723                         {
0724                             hM_ClusEtaPV_ladderZID2_layer1->Fill(hit->Eta(), vtxzwei);
0725                             hM_ClusEtaPV_1p6cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0726                         }
0727                     }
0728                     if (ClusLadderZId->at(i) == 1 || ClusLadderZId->at(i) == 3)
0729                     {
0730                         hM_ClusEtaPV_2cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0731                         if (ClusLadderZId->at(i) == 1)
0732                         {
0733                             hM_ClusEtaPV_ladderZID1_layer1->Fill(hit->Eta(), vtxzwei);
0734                             hM_ClusEtaPV_2cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0735                         }
0736                         if (ClusLadderZId->at(i) == 3)
0737                         {
0738                             hM_ClusEtaPV_ladderZID3_layer1->Fill(hit->Eta(), vtxzwei);
0739                             hM_ClusEtaPV_2cmstrip_layer1->Fill(hit->Eta(), vtxzwei);
0740                         }
0741                     }
0742 
0743                     // cluster adc / phi size cut
0744                     // remove clusters with cluster ADC/cluster phi size between [79,82] and [111,114]
0745                     float clusadcoverphisize = ClusAdc->at(i) / ClusPhiSize->at(i);
0746                     bool Pass_adcoverphisize = !(clusadcoverphisize >= 79 && clusadcoverphisize <= 82) && !(clusadcoverphisize >= 111 && clusadcoverphisize <= 114);
0747                     if (Pass_adcoverphisize)
0748                     {
0749                         hM_ClusEtaPV_layer1_ClusADCg35_ClusADCoverPhisizeCut->Fill(hit->Eta(), vtxzwei);
0750                     }
0751                 }
0752                 // hM_ClusZSize_layer1->Fill(ClusZSize->at(i));
0753                 // hM_ClusEta_ClusZSize_layer1->Fill(hit->Eta(), ClusZSize->at(i));
0754                 // hM_ClusZSize_ClusPhiSize_layer1->Fill(ClusZSize->at(i), ClusPhiSize->at(i));
0755                 // hM_ClusZSize_ClusADC_layer1->Fill(ClusZSize->at(i), ClusAdc->at(i));
0756             }
0757             else
0758             {
0759                 // hM_ClusX_layer2->Fill(ClusX->at(i));
0760                 // hM_ClusY_layer2->Fill(ClusY->at(i));
0761                 hM_ClusZ_layer2->Fill(ClusZ->at(i), vtxzwei);
0762                 // hM_ClusR_layer2->Fill(ClusR->at(i));
0763                 // hM_ClusEtaOri_layer2->Fill(ClusEta->at(i));
0764                 hM_ClusEtaPV_layer2->Fill(hit->Eta(), vtxzwei);
0765                 hM_ClusEtaPV_layer2_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0766                 // hM_ClusPhiOri_layer2->Fill(ClusPhi->at(i));
0767                 hM_ClusPhiPV_layer2->Fill(hit->Phi(), vtxzwei);
0768                 hM_ClusPhiPV_layer2_fine->Fill(hit->Phi(), vtxzwei);
0769                 hM_ClusADC_layer2->Fill(ClusAdc->at(i), vtxzwei);
0770                 hM_ClusZSize_layer2->Fill(ClusZSize->at(i), vtxzwei);
0771                 hM_ClusPhiSize_layer2->Fill(ClusPhiSize->at(i), vtxzwei);
0772                 hM_ClusZ_ClusPhiPV_layer2->Fill(ClusZ->at(i), hit->Phi(), vtxzwei);
0773                 hM_ClusPhiPV_ClusPhiSize_layer2->Fill(hit->Phi(), ClusPhiSize->at(i), vtxzwei);
0774                 hM_ClusEtaPV_ClusADC_layer2->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0775                 hM_ClusEtaPV_ClusADC_layer2_zoomin->Fill(hit->Eta(), ClusAdc->at(i), vtxzwei);
0776                 hM_ClusEtaPV_ClusPhiSize_layer2->Fill(hit->Eta(), ClusPhiSize->at(i), vtxzwei);
0777                 hM_ClusEtaPV_ClusPhiSize_layer2_zoomin->Fill(hit->Eta(), ClusPhiSize->at(i), vtxzwei);
0778                 hM_ClusPhiPV_ClusADC_layer2->Fill(hit->Phi(), ClusAdc->at(i), vtxzwei);
0779                 hM_ClusPhiSize_ClusADC_layer2->Fill(ClusPhiSize->at(i), ClusAdc->at(i), vtxzwei);
0780                 hM_ClusPhiSize_ClusADC_layer2_zoomin->Fill(ClusPhiSize->at(i), ClusAdc->at(i), vtxzwei);
0781                 hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_zoomin->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0782 
0783                 // Select perpendicularly incident clusters, i.e cluster |eta|<=0.1
0784                 if (fabs(hit->Eta()) <= 0.1)
0785                 {
0786                     hM_ClusZSize_layer2_etale0p1->Fill(ClusZSize->at(i), vtxzwei);
0787                     hM_ClusPhiSize_layer2_etale0p1->Fill(ClusPhiSize->at(i), vtxzwei);
0788                     hM_ClusADC_layer2_etale0p1->Fill(ClusAdc->at(i), vtxzwei);
0789                     if (ClusZSize->at(i) == 1)
0790                     {
0791                         hM_ClusADC_layer2_etale0p1_ClusZSize1->Fill(ClusAdc->at(i), vtxzwei);
0792                         hM_ClusPhiSize_layer2_etale0p1_ClusZSize1->Fill(ClusPhiSize->at(i), vtxzwei);
0793                     }
0794                 }
0795 
0796                 if (ClusAdc->at(i) > 35)
0797                 {
0798                     hM_ClusEtaPV_layer2_ClusADCg35->Fill(hit->Eta(), vtxzwei);
0799                     hM_ClusEtaPV_layer2_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0800                     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_ClusADCg35->Fill(hit->Eta(), ClusAdc->at(i) / ClusPhiSize->at(i), vtxzwei);
0801                     hM_ClusEtaPV_PVz_layer2_ClusADCg35->Fill(hit->Eta(), PV[2]);
0802                     hM_ClusEtaPV_PVz_layer2_ClusADCg35_weiClusADCoverPhiSize->Fill(hit->Eta(), PV[2], ClusAdc->at(i) / ClusPhiSize->at(i) * vtxzwei);
0803 
0804                     if (ClusLadderZId->at(i) == 0 || ClusLadderZId->at(i) == 2)
0805                     {
0806                         hM_ClusEtaPV_1p6cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0807                         if (ClusLadderZId->at(i) == 0)
0808                         {
0809                             hM_ClusEtaPV_ladderZID0_layer2->Fill(hit->Eta(), vtxzwei);
0810                             hM_ClusEtaPV_1p6cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0811                         }
0812                         if (ClusLadderZId->at(i) == 2)
0813                         {
0814                             hM_ClusEtaPV_ladderZID2_layer2->Fill(hit->Eta(), vtxzwei);
0815                             hM_ClusEtaPV_1p6cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0816                         }
0817                     }
0818                     if (ClusLadderZId->at(i) == 1 || ClusLadderZId->at(i) == 3)
0819                     {
0820                         hM_ClusEtaPV_2cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0821                         if (ClusLadderZId->at(i) == 1)
0822                         {
0823                             hM_ClusEtaPV_ladderZID1_layer2->Fill(hit->Eta(), vtxzwei);
0824                             hM_ClusEtaPV_2cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0825                         }
0826                         if (ClusLadderZId->at(i) == 3)
0827                         {
0828                             hM_ClusEtaPV_ladderZID3_layer2->Fill(hit->Eta(), vtxzwei);
0829                             hM_ClusEtaPV_2cmstrip_layer2->Fill(hit->Eta(), vtxzwei);
0830                         }
0831                     }
0832 
0833                     // cluster adc / phi size cut
0834                     // remove clusters with cluster ADC/cluster phi size between [79,82] and [111,114]
0835                     float clusadcoverphisize = ClusAdc->at(i) / ClusPhiSize->at(i);
0836                     bool Pass_adcoverphisize = !(clusadcoverphisize >= 79 && clusadcoverphisize <= 82) && !(clusadcoverphisize >= 111 && clusadcoverphisize <= 114);
0837                     if (Pass_adcoverphisize)
0838                     {
0839                         hM_ClusEtaPV_layer2_ClusADCg35_ClusADCoverPhisizeCut->Fill(hit->Eta(), vtxzwei);
0840                     }
0841                 }
0842                 // hM_ClusZSize_layer2->Fill(ClusZSize->at(i));
0843                 // hM_ClusEta_ClusZSize_layer2->Fill(ClusEta->at(i), ClusZSize->at(i));
0844                 // hM_ClusZSize_ClusPhiSize_layer2->Fill(ClusZSize->at(i), ClusPhiSize->at(i));
0845                 // hM_ClusZSize_ClusADC_layer2->Fill(ClusZSize->at(i), ClusAdc->at(i));
0846             }
0847         }
0848 
0849         
0850 
0851         // calculate the deltaR between two clusters
0852         for (size_t i = 0; i < hits_layer[0].size(); i++)
0853         {
0854             for (size_t j = i+1; j < hits_layer[0].size(); j++)
0855             {
0856                 float dR = deltaR(hits_layer[0][i]->Eta(), hits_layer[0][i]->Phi(), hits_layer[0][j]->Eta(), hits_layer[0][j]->Phi());
0857                 hM_mutualdRcluster_all->Fill(dR, vtxzwei);
0858                 hM_mutualdRcluster_layer1->Fill(dR, vtxzwei);
0859             }
0860         }
0861         for (size_t i = 0; i < hits_layer[1].size(); i++)
0862         {
0863             for (size_t j = i+1; j < hits_layer[1].size(); j++)
0864             {
0865                 float dR = deltaR(hits_layer[1][i]->Eta(), hits_layer[1][i]->Phi(), hits_layer[1][j]->Eta(), hits_layer[1][j]->Phi());
0866                 hM_mutualdRcluster_all->Fill(dR, vtxzwei);
0867                 hM_mutualdRcluster_layer2->Fill(dR, vtxzwei);
0868             }
0869         }
0870     }
0871 
0872     f->Close();
0873 
0874     TFile *fout = new TFile(outfilename, "RECREATE");
0875     fout->cd();
0876     // hM_ClusX_all->Write();
0877     // hM_ClusX_layer1->Write();
0878     // hM_ClusX_layer2->Write();
0879     // hM_ClusY_all->Write();
0880     // hM_ClusY_layer1->Write();
0881     // hM_ClusY_layer2->Write();
0882     hM_ClusZ_all->Write();
0883     hM_ClusZ_layer1->Write();
0884     hM_ClusZ_layer2->Write();
0885     // hM_ClusR_all->Write();
0886     // hM_ClusR_layer1->Write();
0887     // hM_ClusR_layer2->Write();
0888     hM_ClusEtaPV_all->Write();
0889     hM_ClusEtaPV_layer1->Write();
0890     hM_ClusEtaPV_layer2->Write();
0891     hM_ClusEtaPV_all_ClusADCg35->Write();
0892     hM_ClusEtaPV_layer1_ClusADCg35->Write();
0893     hM_ClusEtaPV_layer2_ClusADCg35->Write();
0894     hM_ClusEtaPV_all_ClusADCg35_weiClusADCoverPhiSize->Write();
0895     hM_ClusEtaPV_layer1_ClusADCg35_weiClusADCoverPhiSize->Write();
0896     hM_ClusEtaPV_layer2_ClusADCg35_weiClusADCoverPhiSize->Write();
0897     hM_ClusEtaPV_all_weiClusADCoverPhiSize->Write();
0898     hM_ClusEtaPV_layer1_weiClusADCoverPhiSize->Write();
0899     hM_ClusEtaPV_layer2_weiClusADCoverPhiSize->Write();
0900     hM_ClusEtaPV_all_ClusADCg35_ClusADCoverPhisizeCut->Write();
0901     hM_ClusEtaPV_layer1_ClusADCg35_ClusADCoverPhisizeCut->Write();
0902     hM_ClusEtaPV_layer2_ClusADCg35_ClusADCoverPhisizeCut->Write();
0903     hM_ClusEtaPV_PVz_all_ClusADCg35->Write();
0904     hM_ClusEtaPV_PVz_layer1_ClusADCg35->Write();
0905     hM_ClusEtaPV_PVz_layer2_ClusADCg35->Write();
0906     hM_ClusEtaPV_PVz_all_LadderZID0->Write();
0907     hM_ClusEtaPV_PVz_all_LadderZID1->Write();
0908     hM_ClusEtaPV_PVz_all_LadderZID2->Write();
0909     hM_ClusEtaPV_PVz_all_LadderZID3->Write();
0910     hM_ClusEtaPV_PVz_all_ClusADCg35_weiClusADCoverPhiSize->Write();
0911     hM_ClusEtaPV_PVz_layer1_ClusADCg35_weiClusADCoverPhiSize->Write();
0912     hM_ClusEtaPV_PVz_layer2_ClusADCg35_weiClusADCoverPhiSize->Write();
0913     hM_ClusEtaPV_1p6cmstrip_all->Write();
0914     hM_ClusEtaPV_1p6cmstrip_layer1->Write();
0915     hM_ClusEtaPV_1p6cmstrip_layer2->Write();
0916     hM_ClusEtaPV_2cmstrip_all->Write();
0917     hM_ClusEtaPV_2cmstrip_layer1->Write();
0918     hM_ClusEtaPV_2cmstrip_layer2->Write();
0919     hM_ClusEtaPV_ladderZID0_all->Write();
0920     hM_ClusEtaPV_ladderZID0_layer1->Write();
0921     hM_ClusEtaPV_ladderZID0_layer2->Write();
0922     hM_ClusEtaPV_ladderZID1_all->Write();
0923     hM_ClusEtaPV_ladderZID1_layer1->Write();
0924     hM_ClusEtaPV_ladderZID1_layer2->Write();
0925     hM_ClusEtaPV_ladderZID2_all->Write();
0926     hM_ClusEtaPV_ladderZID2_layer1->Write();
0927     hM_ClusEtaPV_ladderZID2_layer2->Write();
0928     hM_ClusEtaPV_ladderZID3_all->Write();
0929     hM_ClusEtaPV_ladderZID3_layer1->Write();
0930     hM_ClusEtaPV_ladderZID3_layer2->Write();
0931     // hM_ClusEtaOri_all->Write();
0932     // hM_ClusEtaOri_layer1->Write();
0933     // hM_ClusEtaOri_layer2->Write();
0934     // hM_ClusPhiOri_all->Write();
0935     // hM_ClusPhiOri_layer1->Write();
0936     // hM_ClusPhiOri_layer2->Write();
0937     hM_ClusPhiPV_all->Write();
0938     hM_ClusPhiPV_layer1->Write();
0939     hM_ClusPhiPV_layer2->Write();
0940     hM_ClusPhiPV_ClusPhiSize43or46->Write();
0941     hM_ClusPhiPV_all_fine->Write();
0942     hM_ClusPhiPV_layer1_fine->Write();
0943     hM_ClusPhiPV_layer2_fine->Write();
0944     hM_ClusADC_all->Write();
0945     hM_ClusADC_layer1->Write();
0946     hM_ClusADC_layer2->Write();
0947     hM_ClusADC_all_etale0p1->Write();
0948     hM_ClusADC_layer1_etale0p1->Write();
0949     hM_ClusADC_layer2_etale0p1->Write();
0950     hM_ClusADC_all_etale0p1_ClusZSize1->Write();
0951     hM_ClusADC_layer1_etale0p1_ClusZSize1->Write();
0952     hM_ClusADC_layer2_etale0p1_ClusZSize1->Write();
0953     hM_ClusZSize_all->Write();
0954     hM_ClusZSize_layer1->Write();
0955     hM_ClusZSize_layer2->Write();
0956     hM_ClusZSize_all_etale0p1->Write();
0957     hM_ClusZSize_layer1_etale0p1->Write();
0958     hM_ClusZSize_layer2_etale0p1->Write();
0959     hM_ClusPhiSize_all->Write();
0960     hM_ClusPhiSize_layer1->Write();
0961     hM_ClusPhiSize_layer2->Write();
0962     hM_ClusPhiSize_all_etale0p1->Write();
0963     hM_ClusPhiSize_layer1_etale0p1->Write();
0964     hM_ClusPhiSize_layer2_etale0p1->Write();
0965     hM_ClusPhiSize_all_etale0p1_ClusZSize1->Write();
0966     hM_ClusPhiSize_layer1_etale0p1_ClusZSize1->Write();
0967     hM_ClusPhiSize_layer2_etale0p1_ClusZSize1->Write();
0968     hM_ClusX_ClusY_all->Write();
0969     hM_ClusX_ClusY_ClusPhiSize43or46->Write();
0970     hM_ClusX_ClusY_ClusPhiSize43->Write();
0971     hM_ClusX_ClusY_ClusPhiSize46->Write();
0972     hM_ClusX_ClusY_all_weiphisize->Write();
0973     hM_ClusX_ClusY_all_weiclusadc->Write();
0974     hM_ClusZ_ClusPhiPV_all_coarse->Write();
0975     hM_ClusZ_ClusPhiPV_all_coarse_weiphisize->Write();
0976     hM_ClusZ_ClusPhiPV_all_coarse_weiclusadc->Write();
0977     hM_ClusZ_ClusPhiPV_all->Write();
0978     hM_ClusZ_ClusPhiPV_layer1->Write();
0979     hM_ClusZ_ClusPhiPV_layer2->Write();
0980     // hM_ClusEta_ClusZSize_all->Write();
0981     // hM_ClusEta_ClusZSize_layer1->Write();
0982     // hM_ClusEta_ClusZSize_layer2->Write();
0983     hM_ClusPhiPV_ClusPhiSize_all->Write();
0984     hM_ClusPhiPV_ClusPhiSize_layer1->Write();
0985     hM_ClusPhiPV_ClusPhiSize_layer2->Write();
0986     hM_ClusEtaPV_ClusADC_all->Write();
0987     hM_ClusEtaPV_ClusADC_layer1->Write();
0988     hM_ClusEtaPV_ClusADC_layer2->Write();
0989     hM_ClusEtaPV_ClusADC_all_zoomin->Write();
0990     hM_ClusEtaPV_ClusADC_layer1_zoomin->Write();
0991     hM_ClusEtaPV_ClusADC_layer2_zoomin->Write();
0992     hM_ClusEtaPV_ClusPhiSize_all->Write();
0993     hM_ClusEtaPV_ClusPhiSize_layer1->Write();
0994     hM_ClusEtaPV_ClusPhiSize_layer2->Write();
0995     hM_ClusEtaPV_ClusPhiSize_all_zoomin->Write();
0996     hM_ClusEtaPV_ClusPhiSize_layer1_zoomin->Write();
0997     hM_ClusEtaPV_ClusPhiSize_layer2_zoomin->Write();
0998     hM_ClusPhiPV_ClusADC_all->Write();
0999     hM_ClusPhiPV_ClusADC_layer1->Write();
1000     hM_ClusPhiPV_ClusADC_layer2->Write();
1001     // hM_ClusZSize_ClusPhiSize_all->Write();
1002     // hM_ClusZSize_ClusPhiSize_layer1->Write();
1003     // hM_ClusZSize_ClusPhiSize_layer2->Write();
1004     // hM_ClusZSize_ClusADC_all->Write();
1005     // hM_ClusZSize_ClusADC_layer1->Write();
1006     // hM_ClusZSize_ClusADC_layer2->Write();
1007     hM_ClusPhiSize_ClusADC_all->Write();
1008     hM_ClusPhiSize_ClusADC_layer1->Write();
1009     hM_ClusPhiSize_ClusADC_layer2->Write();
1010     hM_ClusPhiSize_ClusADC_all_zoomin->Write();
1011     hM_ClusPhiSize_ClusADC_layer1_zoomin->Write();
1012     hM_ClusPhiSize_ClusADC_layer2_zoomin->Write();
1013     hM_ClusEtaPV_ClusADCoverClusPhiSize_all_zoomin->Write();
1014     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_zoomin->Write();
1015     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_zoomin->Write();
1016     hM_ClusEtaPV_ClusADCoverClusPhiSize_all_ClusADCg35->Write();
1017     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer1_ClusADCg35->Write();
1018     hM_ClusEtaPV_ClusADCoverClusPhiSize_layer2_ClusADCg35->Write();
1019     // For cluster ADC cut optimization
1020     for (int constscale = constscale_low; constscale <= constscale_high; constscale += constscale_step)
1021     {
1022         v_hM_ClusEtaPV_EtaDepADCCut[(constscale - constscale_low) / constscale_step]->Write();
1023         v_hM_ClusEtaPV_ClusADC_all_zoomin[(constscale - constscale_low) / constscale_step]->Write();
1024     }
1025     for (int i = 0; i < v_hM_ClusEtaPV_EtaDepADCCut_etascale.size(); i++)
1026     {
1027         v_hM_ClusEtaPV_EtaDepADCCut_etascale[i]->Write();
1028         v_hM_ClusEtaPV_ClusADC_all_zoomin_etascale[i]->Write();
1029     }
1030     // for intt private cluster adc cut
1031     hM_ClusEtaPV_EtaDepADCCut_inttprivate->Write();
1032     hM_ClusEtaPV_ClusADC_all_zoomin_inttprivate->Write();
1033 
1034     hM_mutualdRcluster_all->Write();
1035     hM_mutualdRcluster_layer1->Write();
1036     hM_mutualdRcluster_layer2->Write();
1037     fout->Close();
1038 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!