calib_nov23/towerslope/tsdst_plus_tsfit.C

0001 #pragma once
0002 #if ROOT_VERSION_CODE >= ROOT_VERSION(6,00,0)
0003 #include <fun4all/SubsysReco.h>
0004 #include <fun4all/Fun4AllServer.h>
0005 #include <fun4all/Fun4AllInputManager.h>
0006 #include <fun4allraw/Fun4AllPrdfInputManager.h>
0007 #include <fun4all/Fun4AllDstInputManager.h>
0008 //#include <rawwaveformtowerbuilder/RawWaveformTowerBuilder.h>
0009 #include <fun4all/Fun4AllDstOutputManager.h>
0010 #include <fun4all/Fun4AllUtils.h>
0011
0012 #include <litecaloeval/LiteCaloEval.h>
0013 #include <caloreco/CaloTowerCalib.h>
0014 #include <caloreco/RawClusterBuilderTemplate.h>
0015 //#include <calib_emc_pi0/CaloCalibEmc_Pi0.h>
0016
0017 #include <phool/recoConsts.h>
0018
0019 #include "createLocalEMCalCalibFile.C"
0020 #include "mergeCDBTTrees.C"
0021 #include "TSCtoCDBTTree.C"
0022
0023 R__LOAD_LIBRARY(libfun4all.so)
0024 R__LOAD_LIBRARY(libfun4allraw.so)
0025 R__LOAD_LIBRARY(libcalo_reco.so)
0026 R__LOAD_LIBRARY(libLiteCaloEvalTowSlope.so)
0027 //R__LOAD_LIBRARY(libcalibCaloEmc_pi0.so)
0028
0029
0030
0031 #endif
0032
0033 // to get files from my local area
0034 void tsdst_plus_tsfit(int nevents = 1, const char *fname = "input.list", const char * outfile =  "defout.root", const char * fitoutfile = "deffitout.root")
0035 //,const int runNumber)
0036 {
0037
0038   Fun4AllServer *se = Fun4AllServer::instance();
0039
0040
0041   Fun4AllInputManager *in = new Fun4AllDstInputManager("in");
0042
0043   TString infiletstr(fname);
0044   if (infiletstr.Contains(".list"))
0045     in->AddListFile(fname);
0046   else
0047     in->fileopen(fname);
0048
0049   se->registerInputManager(in);
0050
0051
0052   // initalize CDB to pull what we "think" is the "base" calibration
0053   // into a local cdbttree so we can update it with the tower slope output
0054   recoConsts *rc = recoConsts::instance();
0055   ifstream file(fname);
0056   string first_file;
0057   getline(file, first_file);
0058
0059   pair<int, int> runseg = Fun4AllUtils::GetRunSegment(first_file);
0060   int runnumber = runseg.first;
0061   cout << "run number = " << runnumber << endl;
0062
0063   // global tag
0064   rc->set_StringFlag("CDB_GLOBALTAG","MDC2");
0065   // // 64 bit timestamp
0066   rc->set_uint64Flag("TIMESTAMP",runnumber);
0067
0068   string local_cdb_copy_name = "local_cdb_copy.root";
0069   createLocalEMCalCalibFile(local_cdb_copy_name.c_str(),runnumber);
0070
0071
0072   //JF to Blair Nov 23 -- leaving these next blocks, similar code
0073   //  would be needed if we run a second iteration of towerslope
0074   // which is fairly possible.
0075
0076   //CaloTowerBuilder *ca = new CaloTowerBuilder();
0077
0078   // // recoConsts *rc = recoConsts::instance();
0079   // // rc->set_StringFlag("CDB_GLOBALTAG","ProdA_2023"); // this points to the global  tag in the CDB
0080   // // The calibrations have a validity range set by the beam clock which is not read out of the prdfs as of now
0081   //rc->set_uint64Flag("TIMESTAMP",runNumber);
0082
0083
0084   // CaloTowerCalib *calib = new CaloTowerCalib("CEMCCALIB");
0085   // //  calib->setCalibName("cemc_abscalib_cosmic");// these two lines are needed to choose your own calibration
0086   // //  calib->setFieldName("cemc_abscalib_cosmic");
0087   // calib->set_detector_type(CaloTowerCalib::CEMC);
0088   // se->registerSubsystem(calib);
0089
0090   std::cout << "in tow dst " << std::endl;
0091
0092
0093   LiteCaloEval *eval7e = new LiteCaloEval("CEMCEVALUATOR2", "CEMC",outfile);
0094   //  LiteCaloEval *eval = new LiteCaloEval("HOCEMCEVALUATOR2", "HCALOUT", outputfile.c_str());
0095   //  eval->Verbosity(verbosity);
0096   eval7e->CaloType(LiteCaloEval::CEMC);
0097   eval7e->setInputTowerNodeName("TOWERINFO_CALIB_CEMC");
0098   //eval->CaloType(LiteCaloEval::HCALOUT);
0099   se->registerSubsystem(eval7e);
0100
0101
0102   se->run(nevents);
0103
0104   se->End();
0105   se->PrintTimer();
0106
0107
0108   LiteCaloEval modlce;
0109
0110   //reflce.CaloType(LiteCaloEval::CEMC);
0111   modlce.CaloType(LiteCaloEval::CEMC);
0112
0113
0114   modlce.Get_Histos(outfile,fitoutfile);
0115
0116   modlce.m_myminbin =  8;       // this really is min bin (no -1/+1 needed)
0117                                 //  but not sure if for summer 23 this should
0118                                 // be 8 or 9 I think 8 according to Joey Clement
0119                                 // on ~Nov 1, mattermost post
0120   modlce.m_mymaxbin =  95 + 1 ;  // this actually means maxbin is 1 less
0121   //  modlce.m_mymaxbin =  95+1;  // so that means maxbin is 95 which is
0122   //  really the max eta bin for emc, indexed 0 through 95
0123
0124
0125   modlce.setFitMin(0.12);
0126   modlce.setFitMax(0.7);
0127
0128   modlce.FitRelativeShifts(&modlce,110);
0129
0130
0131   // create the cdbttree from tsc output andd multiply the corrections
0132   // into the base calibration to pickup for pi0 first iteration
0133   TSCtoCDBTTree(fitoutfile,"tsc_output_cdb.root");
0134   mergeCDBTTrees("tsc_output_cdb.root",local_cdb_copy_name.c_str(),"calib_post_TSC.root");
0135
0136   gSystem->Exit(0);
0137 }