/** * \file TrkCore.cpp * \author Elena Vannuccini */ // ......................................................... // ROOT header files // ......................................................... #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include // ......................................................... // other general header files // ......................................................... #include #include // ......................................................... // files in the inc directory // ......................................................... #include #include #include #include #include #include //#include #include // ......................................................... // YODA header files // ......................................................... #include #include #include #include #include // ......................................................... // tracker analysis FORTRAN77 routines // ......................................................... extern "C" { extern struct cTrkCalib pedsigbad_; extern struct cTrkLevel0 level0event_; extern struct cTrkLevel1 level1event_; extern struct cTrkLevel2 level2event_; extern struct cPath path_; extern struct cBPath bpath_; int fillpedsigfromdefault_(); int readmipparam_(); int readchargeparam_(); int readvkmask_(); int readalignparam_(); int readetaparam_(); int reductionflight_(); int analysisflight_(); } // ......................................................... // namespaces // ......................................................... using namespace std; using namespace pamela; // ================================================================================ // // // ================================================================================ /** * \brief Tracker data reduction routine. * * It performs data reduction LEVEL0->LEVEL1->LEVEL2, producing LEVEL1 and or LEVEL2 output, in ROOT or HBOOK format. * Input parameters: * @param run id of the run to be processed (if run=0 a whole file is reprocessed) * @param dbc pointer to BD server * @param file1 output LEVEL1 file name (null if no LEVEL1 output is required) * @param file2 output LEVEL2 file name (null if no LEVEL2 output is required) * @param standalone (flag to run the program in standalone mode, that is without reading RunInfo) */ //short int TrkCore(Int_t run, TSQLServer *dbc, TString file1, TString file2, Bool_t standalone) int TrkCore(ULong64_t run, TFile *f2, TSQLServer *dbc, int ncustom, char*vcustom[]) { // --------------------------- // Define some basic variables // --------------------------- TString Frameworklev1="paw"; TString Frameworklev2="root"; TString filename = 0; TString trkcalibfile = 0; Int_t last_trk_calib_used = 0; Long64_t nevents = 0LL; // LEVEL0 input TFile *f0 = 0; TTree *physicsTree = 0; TBranch *b_registry = 0; TBranch *b_trk = 0; TBranch *b_header = 0; RegistryEvent *reg = 0; EventHeader *header = 0; PscuHeader *pscu = 0; TrkLevel0 *l0_event = 0; // RunInfo ItoRunInfo *runinfo = 0; TArrayL *runlist = 0; ULong64_t from_run; ULong64_t to_run; Bool_t reprocessing = false; //LEVEL1 output - ROOT // TFile *f1; // TTree *t_level1; // TrkLevel1 *l1_event; //LEVEL2 output - ROOT TTree *t_level2 = 0; TTree *t_clone = 0; TrkLevel2 *l2_event = new TrkLevel2(); TrkLevel2 *l2_clone = new TrkLevel2(); // ----------------------- // ----------------------- // Handle input parameters // (data) // // - open/create output files, determining the environment root/hbook from the estension // - create the level1/level2 tree+branch/nt-uple // ----------------------- // ----------------------- TrkProcess *p = new TrkProcess(run,f2); p->HandleCustomPar(ncustom,vcustom); path_.debug = p->DEBUG; // temporaneo.. vedro` di fare di meglio if( p->DEBUG )p->Dump(); // p->Dump(); // =================================== // Open/Create level1 output file // (NB output files other than LEVEL2 are // created in a separate folder) // =================================== if(p->get1){ if(p->ifroot1){ throw -299; }else{ throw -299; }; }; // =================================== // Open/Create level2 output file // =================================== if(p->get2){ // if(p->ifroot2){ // root //------------------------------------------- // read RunInfo //------------------------------------------- if(!(p->standalone)){ // Open "RunInfo" tree and get list of run runinfo = new ItoRunInfo(f2); char *trkversion = TrkInfo(false); if( runinfo->Update(run,"TRK",trkversion) ) throw -205; runlist = runinfo->GetRunList(); reprocessing = runinfo->IsReprocessing();//???? if(p->DEBUG){ cout << "#events "<< runinfo->GetFileEntries() << endl;// #events in the file cout << "#runs "<< runinfo->GetRunEntries() << endl;// #runs in the file cout << "1-st run "<< runlist->At(0) << endl; cout << "last run "<< runlist->At(runinfo->GetRunEntries()-1) << endl; cout << "1-st event "<< runinfo->GetFirstEntry() << endl;// first event of our run cout << "last event+1 "<< runinfo->GetLastEntry() << endl;// first event AFTER the last event of our run cout << "reprocessing "<< runinfo->IsReprocessing() << endl << endl; }; }; //------------------------------------------- // //------------------------------------------- // Take (if present) the old Tracker tree t_clone = (TTree*)f2->Get("Tracker"); if( t_clone != NULL ) t_clone->SetBranchAddress("TrkLevel2",&l2_clone); // Create NEW Tracker tree t_level2 = new TTree("Tracker","PAMELA tracker level2 data "); t_level2->Branch("TrkLevel2","TrkLevel2",&l2_event); // }else{ // hbook // throw -299; // } }; // ------------------------------------------- // define runs to be processed/reprocessed // ------------------------------------------- if(run == 0){ // reprocessing ALL runs if(p->standalone)throw -298; // reprocessing not implemented from_run = runlist->At(0); to_run = runlist->At(runinfo->GetRunEntries()-1); }else{ // processing/reprocessing ONE single run from_run = run; to_run = run; }; // // init "expiration date" of calibrations and parameters // ULong64_t pedsig_time =0ULL; ULong64_t B_time =0ULL; ULong64_t mip_time =0ULL; ULong64_t charge_time =0ULL; ULong64_t eta_time =0ULL; ULong64_t align_time =0ULL; ULong64_t mask_time =0ULL; // // init event counter // Int_t ev_count =0; // // create query-results objects // GL_RUN q1 = GL_RUN(); GL_TRK_CALIB q2 = GL_TRK_CALIB(); GL_ROOT q3 = GL_ROOT(); GL_PARAM q4 = GL_PARAM(); // ------------------------------------------------------------ // if reprocessing one run, copy all the events BEFORE the run // ------------------------------------------------------------ if( !(p->standalone) ){ for(UInt_t i=0; iGetFirstEntry(); i++){ t_clone->GetEntry(i); *l2_event = *l2_clone; t_level2->Fill(); l2_event->Clear(); // COPY COPY COPY }; }; // ------------------------------------------------------------ // ------------------------------------------------------------ // START LOOP OVER RUNS TO PROCESSED/REPROCESSED // ------------------------------------------------------------ // ------------------------------------------------------------ for(ULong64_t idRun=from_run; idRun <= to_run; idRun++){ if(p->DEBUG) cout << " ========================= Run: "<< idRun << endl; ULong64_t runheadtime = 0ULL; ULong64_t runtrailtime = 0ULL; UInt_t evfrom = 0; UInt_t evto = 0; Int_t id_reg_run =-1; Int_t trk_calib_used = 0; if(p->standalone){ // ============================== // first query: retrieve run info // ============================== if (q1.Query_GL_RUN(idRun,dbc) )throw -50; id_reg_run = q1.ID_REG_RUN; runheadtime = q1.RUNHEADER_TIME; runtrailtime = q1.RUNTRAILER_TIME; evfrom = q1.EV_REG_PHYS_FROM; evto = q1.EV_REG_PHYS_TO; trk_calib_used = q1.TRK_CALIB_USED; }else{ // ============================== // get run info from RunInfo tree // ============================== if( runinfo->GetRunInfo(idRun) ) throw -205; id_reg_run = runinfo->ID_REG_RUN; runheadtime = runinfo->RUNHEADER_TIME; runtrailtime = runinfo->RUNTRAILER_TIME; evfrom = runinfo->EV_REG_PHYS_FROM; evto = runinfo->EV_REG_PHYS_TO; trk_calib_used = runinfo->TRK_CALIB_USED; }; // if(p->DEBUG){ cout << "ROOT file ID "<< id_reg_run << endl; cout << "RunHeader time "<< runheadtime << endl; cout << "RunTrailer time "<< runtrailtime << endl; cout << " from event "<< evfrom << endl; cout << " to event "<< evto << endl; cout << "trk-calibration used "<< trk_calib_used << endl; }; // ======================================================== // second query: search the LEVEL0 file that contains idRun // ======================================================== TString lastfilename = filename; if( q3.Query_GL_ROOT(id_reg_run,dbc) )throw -51; filename = q3.PATH + q3.NAME; // ======================================================== // Open the input LEVEL0 data file // ======================================================== if(filename.CompareTo(lastfilename)){ if(!lastfilename.IsNull())f0->Close(); //if( debug ) cout << "Opening LEVEL0 file: "<< filename << endl; if(p->DEBUG) cout << "Opening LEVEL0 file: "<< filename << endl; f0 = new TFile(filename); /* if ( !f0 ){ printf("\n\n ERROR: no data file, exiting...\n\n"); return(103); };*/ if ( !f0 ) throw -6; physicsTree = (TTree*)f0->Get("Physics"); if(!physicsTree) throw -7; b_header = physicsTree ->GetBranch("Header"); if(!b_header ) throw -8; b_registry = physicsTree ->GetBranch("Registry"); if(!b_registry ) throw -9; b_trk = physicsTree ->GetBranch("Tracker"); if(!b_trk ) throw -203; physicsTree->SetBranchAddress("Tracker" ,&l0_event); physicsTree->SetBranchAddress("Registry",®); physicsTree->SetBranchAddress("Header",&header); nevents = physicsTree->GetEntries(); if ( nevents < 1 ) throw -11; }; // ============================================================= // retrieve calibration file needed to process the run // ============================================================= if(trk_calib_used==104 && last_trk_calib_used!=104){ if(p->DEBUG )cout << "** Using default calibration **" << endl; if (q4.Query_GL_PARAM(runheadtime,"default calibration",dbc) )throw -52; path_.FillWith(q4.PATH+q4.NAME); fillpedsigfromdefault_(); if(path_.error) throw -216; }else if( (trk_calib_used==1 | trk_calib_used==2) && runheadtime > pedsig_time){ if( q2.Query_GL_TRK_CALIB(runheadtime,dbc) )throw -53; pedsig_time = q2.TO_TIME; if(q2.EV_REG_CALIBTRK1 != q2.EV_REG_CALIBTRK2) printf("WARNING!! ---> EV_REG_CALIBTRK1=%d it's different from EV_REG_CALIBTRK2=%d \n\n",q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2); if( q3.Query_GL_ROOT(q2.ID_REG_CALIBTRK,dbc) )throw -51; trkcalibfile = q3.PATH + q3.NAME; // store calibration informations if(p->DEBUG) cout << "Online calibration: "<< trkcalibfile << endl; // pedsigbad_ = TrkFlightPede(trkcalibfile,q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2); TFile* f0_c= f0; if(trkcalibfile.CompareTo(filename))f0_c = new TFile(trkcalibfile); pedsigbad_.FillFrom(f0_c,q2.EV_REG_CALIBTRK1,q2.EV_REG_CALIBTRK2); if(trkcalibfile.CompareTo(filename))f0_c->Close(); }; last_trk_calib_used = trk_calib_used; // ============================================================= // retrieve information about parameters to process LEVEL2 // ============================================================= // // magnetic field // if(p->DEBUG) cout << "Loading parameter files : "<< endl; if(runheadtime > B_time){ if( q4.Query_GL_PARAM(runheadtime,"field",dbc) )throw -52; B_time = q4.TO_TIME; l2_event->LoadField(q4.PATH+q4.NAME); l2_event->LoadField(q4.PATH+q4.NAME); if(path_.error) throw -215; }; // // mip conversion parameters // if(runheadtime > mip_time){ if( q4.Query_GL_PARAM(runheadtime,"trk mip",dbc) )throw -52; mip_time = q4.TO_TIME; path_.FillWith(q4.PATH+q4.NAME); readmipparam_(); if(path_.error) throw -212; }; // // charge correlation parameters // if(runheadtime > charge_time){ if( q4.Query_GL_PARAM(runheadtime,"trk charge",dbc) )throw -52; charge_time = q4.TO_TIME; path_.FillWith(q4.PATH+q4.NAME); readchargeparam_(); if(path_.error) throw -213; }; // // eta p.f.a. parameters // if(runheadtime > eta_time){ if( q4.Query_GL_PARAM(runheadtime,"trk pfa",dbc) )throw -52; eta_time = q4.TO_TIME; path_.FillWith(q4.PATH+q4.NAME); readetaparam_(); if(path_.error) throw -214; }; // // alignment parameters // if(runheadtime > align_time){ if( q4.Query_GL_PARAM(runheadtime,"trk alignment",dbc) )throw -52; align_time = q4.TO_TIME; path_.FillWith(q4.PATH+q4.NAME); readalignparam_(); if(path_.error) throw -211; }; // // viking mask // if(runheadtime > mask_time){ if( q4.Query_GL_PARAM(runheadtime,"trk mask",dbc) )throw -52; mask_time = q4.TO_TIME; path_.FillWith(q4.PATH+q4.NAME); readvkmask_(); if(path_.error) throw -210; }; // // if( p->DEBUG ){ cout<<"Total number of events : "<DEBUG && re%100 == 0 && re > 0 ) printf(" %iK \n",re/1000); if ( p->DEBUG && re%100 == 0 && re > 0 ) cout << "."; b_registry->GetEntry(re); b_trk->GetEntry(reg->event); b_header->GetEntry(reg->event); pscu = header->GetPscuHeader(); // pscu->Print(); // cout <Counter << endl; l0_event->GetCommonVar(&level0event_); // ============================= if(p->get1|p->get2) reductionflight_(); if(p->get2) analysisflight_(); // ============================= // ---------------- // LEVEL1 output // ---------------- if(p->get1){ if(p->ifroot1){ // root }else{ // hbook }; }; // ---------------- // LEVEL2 output // ---------------- if(p->get2){ // if(p->ifroot2){ // root l2_event->FillCommonVar(&level2event_); // l2_event->Dump(); t_level2->Fill(); l2_event->Clear(); // }else{ // hbook // } }; }; // end loop on events if(p->DEBUG)cout << " >>> processed "<< ev_count <<" events"<< endl; }; // end loop on runs // ------------------------------------------------------------ // if reprocessing one run, copy all the events AFTER the run // ------------------------------------------------------------ if( !(p->standalone) ){ for(UInt_t i=runinfo->GetLastEntry()+1; iGetFileEntries(); i++){ t_clone->GetEntry(i); *l2_event = *l2_clone; t_level2->Fill(); l2_event->Clear(); // COPY COPY COPY }; }; // --------------- // close the files // --------------- if(p->get1){ if(p->ifroot1){ // root }else{ // hbook } }; if(p->get2){ // if(p->ifroot2){ // root if( t_clone )t_clone->Delete("all");//delete old tree from file if( !(p->standalone) )runinfo->Close(); f2->Write("Tracker"); if( t_level2 )t_level2->Delete(); //delete new tree from memory // }else{ // hbook // } }; f0->Close(); return(0); }