File indexing completed on 2023-10-25 09:34:27
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038 #include <TROOT.h>
0039 #include <TChain.h>
0040 #include <TFile.h>
0041 #include <TH1D.h>
0042 #include <TH2F.h>
0043 #include <TProfile.h>
0044 #include <TFitResult.h>
0045 #include <TFitResultPtr.h>
0046 #include <TStyle.h>
0047 #include <TCanvas.h>
0048 #include <TLegend.h>
0049 #include <TPaveStats.h>
0050 #include <TPaveText.h>
0051
0052 #include <algorithm>
0053 #include <iomanip>
0054 #include <iostream>
0055 #include <fstream>
0056 #include <sstream>
0057 #include <map>
0058 #include <vector>
0059 #include <string>
0060
0061 void unpackDetId(unsigned int, int&, int&, int&, int&, int&);
0062
0063 #include "CalibMonitor.C"
0064 #include "CalibPlotProperties.C"
0065 #include "CalibTree.C"
0066
0067 int main(Int_t argc, Char_t* argv[]) {
0068 if (argc < 10) {
0069 std::cerr << "Please give N arguments \n"
0070 << "Mode (0 CalibMonitor; 1 CalibProperties; 2 CalibTree; 3 CalibSplit)\n"
0071 << "Input File Name\n"
0072 << "Output File Name(ROOT)\n"
0073 << "Flag\n"
0074 << "Directory Name\n"
0075 << "Prefix\n"
0076 << "PUcorr\n"
0077 << "Truncate\n"
0078 << "Nmax\n"
0079 << " .... Other parameters depending on mode\n"
0080 << std::endl;
0081 return -1;
0082 }
0083
0084 int mode = std::atoi(argv[1]);
0085 const char* infile = argv[2];
0086 std::string histfile(argv[3]);
0087 int flag = std::atoi(argv[4]);
0088 const char* dirname = argv[5];
0089 std::string prefix(argv[6]);
0090 int pucorr = std::atoi(argv[7]);
0091 int truncate = std::atoi(argv[8]);
0092 Long64_t nmax = static_cast<Long64_t>(std::atoi(argv[9]));
0093
0094 if (mode == 0) {
0095
0096 bool datamc = (argc > 10) ? (std::atoi(argv[10]) > 0) : true;
0097 int numb = (argc > 11) ? std::atoi(argv[11]) : 50;
0098 bool usegen = (argc > 12) ? (std::atoi(argv[12]) > 0) : false;
0099 double scale = (argc > 13) ? std::atof(argv[13]) : 1.0;
0100 int usescale = (argc > 14) ? std::atoi(argv[14]) : 0;
0101 int etalo = (argc > 15) ? std::atoi(argv[15]) : 0;
0102 int etahi = (argc > 16) ? std::atoi(argv[16]) : 30;
0103 int runlo = (argc > 17) ? std::atoi(argv[17]) : 0;
0104 int runhi = (argc > 18) ? std::atoi(argv[18]) : 99999999;
0105 int phimin = (argc > 19) ? std::atoi(argv[19]) : 1;
0106 int phimax = (argc > 20) ? std::atoi(argv[20]) : 72;
0107 int zside = (argc > 21) ? std::atoi(argv[21]) : 1;
0108 int nvxlo = (argc > 22) ? std::atoi(argv[22]) : 0;
0109 int nvxhi = (argc > 23) ? std::atoi(argv[23]) : 1000;
0110 bool exclude = (argc > 24) ? (std::atoi(argv[24]) > 0) : false;
0111 bool etamax = (argc > 25) ? (std::atoi(argv[25]) > 0) : false;
0112 bool debug = (argc > 26) ? (std::atoi(argv[26]) > 0) : false;
0113 bool append = (argc > 27) ? (std::atoi(argv[27]) > 0) : true;
0114 bool all = (argc > 28) ? (std::atoi(argv[28]) > 0) : true;
0115 const char* corrfile = (argc > 29) ? argv[29] : "";
0116 const char* rcorfile = (argc > 30) ? argv[30] : "";
0117 const char* dupfile = (argc > 31) ? argv[31] : "";
0118 const char* rbxfile = (argc > 32) ? argv[32] : "";
0119 const char* comfile = (argc > 33) ? argv[33] : "";
0120 const char* outfile = (argc > 34) ? argv[34] : "";
0121 if (strcmp(corrfile, "junk.txt") == 0)
0122 corrfile = "";
0123 if (strcmp(rcorfile, "junk.txt") == 0)
0124 rcorfile = "";
0125 if (strcmp(dupfile, "junk.txt") == 0)
0126 dupfile = "";
0127 if (strcmp(comfile, "junk.txt") == 0)
0128 comfile = "";
0129 if (strcmp(rbxfile, "junk.txt") == 0)
0130 rbxfile = "";
0131
0132 std::cout << "Execute CalibMonitor with infile:" << infile << " dirName: " << dirname << " dupFile: " << dupfile
0133 << " comFile:" << comfile << " outFile:" << outfile << " prefix: " << prefix << " corrFile: " << corrfile
0134 << " rcoFile:" << rcorfile << " puCorr:" << pucorr << " Flag:" << flag << " Numb:" << numb
0135 << " dataMC:" << datamc << " truncate:" << truncate << " useGen:" << usegen << " Scale:" << scale
0136 << " useScale:" << usescale << " etaRange:" << etalo << ":" << etahi << " runRange:" << runlo << ":"
0137 << runhi << " phiRange:" << phimin << ":" << phimax << " zside:" << zside << " nvxRange:" << nvxlo << ":"
0138 << nvxhi << " rbxFile:" << rbxfile << " exclude:" << exclude << " etaMax:" << etamax
0139 << " histFile:" << histfile << " append:" << append << " all:" << all << " nmax:" << nmax << std::endl;
0140
0141 CalibMonitor c1(infile,
0142 dirname,
0143 dupfile,
0144 comfile,
0145 outfile,
0146 prefix,
0147 corrfile,
0148 rcorfile,
0149 pucorr,
0150 flag,
0151 numb,
0152 datamc,
0153 truncate,
0154 usegen,
0155 scale,
0156 usescale,
0157 etalo,
0158 etahi,
0159 runlo,
0160 runhi,
0161 phimin,
0162 phimax,
0163 zside,
0164 nvxlo,
0165 nvxhi,
0166 rbxfile,
0167 exclude,
0168 etamax);
0169 c1.Loop(nmax, debug);
0170 c1.savePlot(histfile, append, all);
0171 } else if (mode == 1) {
0172
0173 bool datamc = (argc > 10) ? (std::atoi(argv[10]) > 0) : true;
0174 bool usegen = (argc > 11) ? (std::atoi(argv[11]) > 0) : false;
0175 double scale = (argc > 12) ? std::atof(argv[12]) : 1.0;
0176 int usescale = (argc > 13) ? std::atoi(argv[13]) : 0;
0177 int etalo = (argc > 14) ? std::atoi(argv[14]) : 0;
0178 int etahi = (argc > 15) ? std::atoi(argv[15]) : 30;
0179 int runlo = (argc > 16) ? std::atoi(argv[16]) : 0;
0180 int runhi = (argc > 17) ? std::atoi(argv[17]) : 99999999;
0181 int phimin = (argc > 18) ? std::atoi(argv[18]) : 1;
0182 int phimax = (argc > 19) ? std::atoi(argv[19]) : 72;
0183 int zside = (argc > 20) ? std::atoi(argv[20]) : 1;
0184 int nvxlo = (argc > 21) ? std::atoi(argv[21]) : 0;
0185 int nvxhi = (argc > 22) ? std::atoi(argv[22]) : 1000;
0186 bool exclude = (argc > 23) ? (std::atoi(argv[23]) > 0) : false;
0187 bool etamax = (argc > 24) ? (std::atoi(argv[24]) > 0) : false;
0188 bool append = (argc > 25) ? (std::atoi(argv[25]) > 0) : true;
0189 bool all = (argc > 26) ? (std::atoi(argv[26]) > 0) : true;
0190 const char* corrfile = (argc > 27) ? argv[27] : "";
0191 const char* rcorfile = (argc > 28) ? argv[28] : "";
0192 const char* dupfile = (argc > 29) ? argv[29] : "";
0193 const char* rbxfile = (argc > 30) ? argv[30] : "";
0194 if (strcmp(corrfile, "junk.txt") == 0)
0195 corrfile = "";
0196 if (strcmp(rcorfile, "junk.txt") == 0)
0197 rcorfile = "";
0198 if (strcmp(dupfile, "junk.txt") == 0)
0199 dupfile = "";
0200 if (strcmp(rbxfile, "junk.txt") == 0)
0201 rbxfile = "";
0202 bool debug(false);
0203
0204 CalibPlotProperties c1(infile,
0205 dirname,
0206 dupfile,
0207 prefix,
0208 corrfile,
0209 rcorfile,
0210 pucorr,
0211 flag,
0212 datamc,
0213 truncate,
0214 usegen,
0215 scale,
0216 usescale,
0217 etalo,
0218 etahi,
0219 runlo,
0220 runhi,
0221 phimin,
0222 phimax,
0223 zside,
0224 nvxlo,
0225 nvxhi,
0226 rbxfile,
0227 exclude,
0228 etamax);
0229 c1.Loop(nmax);
0230 c1.savePlot(histfile, append, all, debug);
0231 } else if (mode == 2) {
0232
0233 int maxIter = (argc > 10) ? std::atoi(argv[10]) : 30;
0234 const char* corrfile = (argc > 11) ? argv[11] : "";
0235 int applyl1 = (argc > 12) ? std::atoi(argv[12]) : 1;
0236 double l1cut = (argc > 13) ? std::atof(argv[13]) : 0.5;
0237 bool useiter = (argc > 14) ? (std::atoi(argv[14]) > 0) : true;
0238 bool useweight = (argc > 15) ? (std::atoi(argv[15]) > 0) : true;
0239 bool usemean = (argc > 16) ? (std::atoi(argv[16]) > 0) : false;
0240 int nmin = (argc > 17) ? std::atoi(argv[17]) : 0;
0241 bool inverse = (argc > 18) ? (std::atoi(argv[18]) > 0) : true;
0242 double ratmin = (argc > 19) ? std::atof(argv[19]) : 0.25;
0243 double ratmax = (argc > 20) ? std::atof(argv[20]) : 3.0;
0244 int ietamax = (argc > 21) ? std::atoi(argv[21]) : 25;
0245 int ietatrack = (argc > 22) ? std::atoi(argv[22]) : -1;
0246 int sysmode = (argc > 23) ? std::atoi(argv[23]) : -1;
0247 int rcorform = (argc > 24) ? std::atoi(argv[24]) : 0;
0248 bool usegen = (argc > 25) ? (std::atoi(argv[25]) > 0) : false;
0249 int runlo = (argc > 26) ? std::atoi(argv[26]) : 0;
0250 int runhi = (argc > 27) ? std::atoi(argv[27]) : 99999999;
0251 int phimin = (argc > 28) ? std::atoi(argv[28]) : 1;
0252 int phimax = (argc > 29) ? std::atoi(argv[29]) : 72;
0253 int zside = (argc > 30) ? std::atoi(argv[30]) : 0;
0254 int nvxlo = (argc > 31) ? std::atoi(argv[31]) : 0;
0255 int nvxhi = (argc > 32) ? std::atoi(argv[32]) : 1000;
0256 bool exclude = (argc > 33) ? (std::atoi(argv[33]) > 0) : false;
0257 int higheta = (argc > 34) ? std::atoi(argv[34]) : 1;
0258 double fraction = (argc > 35) ? std::atof(argv[35]) : 1.0;
0259 bool writehisto = (argc > 36) ? (std::atoi(argv[36]) > 0) : false;
0260 double pmin = (argc > 37) ? std::atof(argv[37]) : 40.0;
0261 double pmax = (argc > 38) ? std::atof(argv[38]) : 60.0;
0262 bool debug = (argc > 39) ? (std::atoi(argv[39]) > 0) : false;
0263 const char* rcorfile = (argc > 40) ? argv[40] : "";
0264 const char* dupfile = (argc > 41) ? argv[41] : "";
0265 const char* rbxfile = (argc > 42) ? argv[42] : "";
0266 const char* treename = (argc > 43) ? argv[43] : "CalibTree";
0267 if (strcmp(rcorfile, "junk.txt") == 0)
0268 rcorfile = "";
0269 if (strcmp(dupfile, "junk.txt") == 0)
0270 dupfile = "";
0271 if (strcmp(rbxfile, "junk.txt") == 0)
0272 rbxfile = "";
0273
0274 char name[500];
0275 sprintf(name, "%s/%s", dirname, treename);
0276 TChain* chain = new TChain(name);
0277 std::cout << "Create a chain for " << name << " from " << infile << std::endl;
0278
0279 if (!fillChain(chain, infile)) {
0280 std::cout << "*****No valid tree chain can be obtained*****" << std::endl;
0281 } else {
0282 std::cout << "Proceed with a tree chain with " << chain->GetEntries() << " entries" << std::endl;
0283 Long64_t nentryTot = chain->GetEntries();
0284 Long64_t nentries = (fraction > 0.01 && fraction < 0.99) ? (Long64_t)(fraction * nentryTot) : nentryTot;
0285 static const int maxIterMax = 100;
0286 if (maxIter > maxIterMax)
0287 maxIter = maxIterMax;
0288 std::cout << "Tree " << name << " " << chain << " in directory " << dirname << " from file " << infile
0289 << " with nentries (tracks): " << nentries << std::endl;
0290 unsigned int k(0), kmax(maxIter);
0291 std::cout << "Proceed using CalibTree with dupFile:" << dupfile << " rcorFile:" << rcorfile
0292 << " trunCate:" << truncate << " useIter:" << useiter << " useMean:" << usemean << " runRange:" << runlo
0293 << ":" << runhi << " phiRange:" << phimin << ":" << phimax << " zSide:" << zside
0294 << " nvxRange:" << nvxlo << ":" << nvxhi << " sysMode:" << sysmode << " rbxFile:" << rbxfile
0295 << " puCorr:" << pucorr << " rcorForm:" << rcorform << " useGen:" << usegen << " exclude:" << exclude
0296 << " highEta:" << higheta << " pRange:" << pmin << ":" << pmax << std::endl;
0297 CalibTree t(dupfile,
0298 rcorfile,
0299 truncate,
0300 useiter,
0301 usemean,
0302 runlo,
0303 runhi,
0304 phimin,
0305 phimax,
0306 zside,
0307 nvxlo,
0308 nvxhi,
0309 sysmode,
0310 rbxfile,
0311 pucorr,
0312 rcorform,
0313 usegen,
0314 exclude,
0315 higheta,
0316 pmin,
0317 pmax,
0318 chain);
0319 t.h_pbyE = new TH1D("pbyE", "pbyE", 100, -1.0, 9.0);
0320 t.h_Ebyp_bfr = new TProfile("Ebyp_bfr", "Ebyp_bfr", 60, -30, 30, 0, 10);
0321 t.h_Ebyp_aftr = new TProfile("Ebyp_aftr", "Ebyp_aftr", 60, -30, 30, 0, 10);
0322 t.h_cvg = new TH1D("Cvg0", "Convergence", kmax, 0, kmax);
0323 t.h_cvg->SetMarkerStyle(7);
0324 t.h_cvg->SetMarkerSize(5.0);
0325
0326 TFile* fout = new TFile(histfile.c_str(), "RECREATE");
0327 std::cout << "Output file: " << histfile << " opened in recreate mode" << std::endl;
0328 fout->cd();
0329
0330 double cvgs[maxIterMax], itrs[maxIterMax];
0331 t.getDetId(fraction, ietatrack, debug, nmax);
0332
0333 for (; k <= kmax; ++k) {
0334 std::cout << "Calling Loop() " << k << "th time" << std::endl;
0335 double cvg = t.Loop(k,
0336 fout,
0337 useweight,
0338 nmin,
0339 inverse,
0340 ratmin,
0341 ratmax,
0342 ietamax,
0343 ietatrack,
0344 applyl1,
0345 l1cut,
0346 k == kmax,
0347 fraction,
0348 writehisto,
0349 debug,
0350 nmax);
0351 itrs[k] = k;
0352 cvgs[k] = cvg;
0353 if (cvg < 0.00001)
0354 break;
0355 }
0356
0357 t.writeCorrFactor(corrfile, ietamax);
0358
0359 fout->cd();
0360 TGraph* g_cvg;
0361 g_cvg = new TGraph(k, itrs, cvgs);
0362 g_cvg->SetMarkerStyle(7);
0363 g_cvg->SetMarkerSize(5.0);
0364 g_cvg->Draw("AP");
0365 g_cvg->Write("Cvg");
0366 std::cout << "Finish looping after " << k << " iterations" << std::endl;
0367 t.makeplots(ratmin, ratmax, ietamax, useweight, fraction, debug, nmax);
0368 fout->Close();
0369 }
0370 } else {
0371
0372 double pmin = (argc > 10) ? std::atof(argv[10]) : 40.0;
0373 double pmax = (argc > 11) ? std::atof(argv[11]) : 60.0;
0374 bool debug = (argc > 12) ? (std::atoi(argv[12]) > 0) : false;
0375 CalibSplit c1(infile, dirname, histfile, pmin, pmax, debug);
0376 c1.Loop(nmax);
0377 }
0378 return 0;
0379 }
