--- calo/ground/LEVEL2/macros/CaloLEVEL2.cra 2006/01/13 09:50:10 1.2 +++ calo/ground/LEVEL2/macros/CaloLEVEL2.cra 2006/03/09 15:55:29 1.3 @@ -1,12 +1,16 @@ // // Given a calibration and a data file this program create an ntuple with LEVEL2 calorimeter variables - Emiliano Mocchiutti // -// CaloLEVEL2.c version 4.01 (2006-01-11) +// CaloLEVEL2.c version 4.03 (2006-03-09) // // The only input needed is the path to the directory created by YODA for the data file you want to analyze. // // Changelog: // +// 4.02 - 4.03 (2006-03-09): Adapted to work with tracker software release r3v02. +// +// 4.01 - 4.02 (2006-01-26): Bug: wrong baselines calculation in some cases, fixed. +// // 4.00 - 4.01 (2006-01-11): Bugs: in Makefile (nothing to worry about) and here clevel1.trkchi2 never filled! fixed. Not really processing self trigger events! fixed. // Bug: not setting small values to zero after applying the cross-talk correction! fixed. // @@ -136,1364 +140,1373 @@ #endif #if defined (__CINT__) void unloadf77lib(TString Framework, Int_t MySQL){ - const char *franame = Framework; - // - char *pam_lib=gSystem->Getenv("PAM_LIB"); - if ( MySQL == 1 || MySQL == -1 ){ - stringstream libload; - libload.str(""); - libload << pam_lib << "/libreadmy"; - gSystem->Unload(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libreadmy_C"; - gSystem->Unload(libload.str().c_str()); - }; - if ( !strcmp(franame,"paw") ) { - stringstream libload; - libload.str(""); - libload << pam_lib << "/libclcalol2_C"; - gSystem->Unload(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libclcalol2"; - gSystem->Unload(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libopcalol2_C"; - gSystem->Unload(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libopcalol2"; - gSystem->Unload(libload.str().c_str()); - }; + const char *franame = Framework; + // + char *pam_lib=gSystem->Getenv("PAM_LIB"); + if ( MySQL == 1 || MySQL == -1 ){ stringstream libload; libload.str(""); - libload << pam_lib << "/libcrcalol2_C"; + libload << pam_lib << "/libreadmy"; gSystem->Unload(libload.str().c_str()); - char *pam_lib=gSystem->Getenv("PAM_LIB"); - stringstream libload; libload.str(""); - libload << pam_lib << "/libcrcalol2"; + libload << pam_lib << "/libreadmy_C"; gSystem->Unload(libload.str().c_str()); - char *pam_lib=gSystem->Getenv("PAM_LIB"); + }; + if ( !strcmp(franame,"paw") ) { stringstream libload; libload.str(""); - libload << pam_lib << "/libreadb2maps_C"; + libload << pam_lib << "/libclcalol2_C"; gSystem->Unload(libload.str().c_str()); - char *pam_lib=gSystem->Getenv("PAM_LIB"); - stringstream libload; libload.str(""); - libload << pam_lib << "/libreadb2maps"; + libload << pam_lib << "/libclcalol2"; gSystem->Unload(libload.str().c_str()); - char *pam_lib=gSystem->Getenv("PAM_LIB"); - stringstream libload; libload.str(""); - libload << pam_lib << "/libtrack"; + libload << pam_lib << "/libopcalol2_C"; gSystem->Unload(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libopcalol2"; + gSystem->Unload(libload.str().c_str()); + }; + stringstream libload; + libload.str(""); + libload << pam_lib << "/libcrcalol2_C"; + gSystem->Unload(libload.str().c_str()); + char *pam_lib=gSystem->Getenv("PAM_LIB"); + stringstream libload; + libload.str(""); + libload << pam_lib << "/libcrcalol2"; + gSystem->Unload(libload.str().c_str()); + char *pam_lib=gSystem->Getenv("PAM_LIB"); + stringstream libload; + libload.str(""); + libload << pam_lib << "/libreadb2maps_C"; + gSystem->Unload(libload.str().c_str()); + char *pam_lib=gSystem->Getenv("PAM_LIB"); + stringstream libload; + libload.str(""); + libload << pam_lib << "/libreadb2maps"; + gSystem->Unload(libload.str().c_str()); + char *pam_lib=gSystem->Getenv("PAM_LIB"); + stringstream libload; + libload.str(""); + libload << pam_lib << "/libtrack"; + gSystem->Unload(libload.str().c_str()); } #endif short int calolevel2core(Int_t ei, Int_t b[4] , TTree *otr, TTree *tree, CaLevel2 & clevel2, Evento & evento, Calib & calib, TString calcalibfile){ + // + // Define variables + // + struct CaLevel1 clevel1; + clevel1.xalig = 120.4; + clevel1.yalig = 118.6; + clevel1.zalig = -260.; + struct Mystruct mystruct; + clevel1.paw = calib.ispaw; + Int_t done = 0; + Int_t pre = -1; + Int_t etime; + Int_t S3 = 0; + Int_t S2 = 0; + Int_t S12 = 0; + Int_t S11 = 0; + Int_t se = 5; + bool isCOMP = 0; + bool isFULL = 0; + bool isRAW = 0; + Float_t ener; + Float_t estripnull[96][22][2]; + // + Int_t doneb = 0; + Int_t donec = 0; + Int_t ck = 0; + Int_t ipre = 0; + Int_t ip[3] = {0}; + Float_t base0, base1, base2; + base0 = 0.; + base1 = 0.; + base2 = 0.; + Float_t qpre[6] = {0}; + Float_t ene[96] = {0}; + // + pamela::calorimeter::CalorimeterEvent *de = 0; + pamela::PscuHeader *ph = 0; + pamela::EventHeader *eh = 0; + pamela::trigger::TriggerEvent *trig = 0; + otr->SetBranchAddress("Header", &eh); + otr->SetBranchAddress("Calorimeter.Event", &de); + otr->SetBranchAddress("Trigger.Event", &trig); + // + // other variables + // + char *yodala2; + Int_t chdone[4] = {0,0,0,0}; + Int_t pe = 0; + Float_t tmptrigty = -1.; + // + // start processing + // + // + // read from the header of the event the OBT and event number + // + otr->GetEntry(ei); + ph = eh->GetPscuHeader(); + clevel2.pkt_num = ph->GetCounter(); + etime = ph->GetOrbitalTime(); + clevel2.OBT = etime; + clevel2.pro_num = (ei+1); + // + // determine who generate the trigger for this event (TOF, S4/PULSER, CALO) + // + S3 = 0; + S2 = 0; + S12 = 0; + S11 = 0; + S3 = trig->patterntrig[2]; + S2 = trig->patterntrig[3]; + S12 = trig->patterntrig[4]; + S11 = trig->patterntrig[5]; + if ( trig->patterntrig[0] ) tmptrigty = 2.; + if ( S3 || S2 || S12 || S11 ) tmptrigty = 0.; + if ( trig->patterntrig[1] & (1<<0) || (!trig->patterntrig[0] && !S3 && !S2 && !S12 && !S11) ) tmptrigty = 1.; + clevel2.trigty = tmptrigty; + // + yodala2 = ""; + if ( calib.mysql ){ + if ( !calib.obtjump ){ + for (Int_t s = 0; s < 4; s++){ + if ( etime > calib.time[s][0] ){ + tryagain: + printf(" ** SECTION %i **\n",s); + stringstream qy; + qy.str(""); + Int_t chkmysql = 0; + stringstream myfile; + const char *tabula = calib.tablename; + qy << "select * from " << tabula; + qy << " where section = " << s; + qy << " and fromtime >= " << calib.time[s][0]; + qy << " limit 1"; + chkmysql = creadmysql(calib.db,(char *)qy.str().c_str(),mystruct); + if ( chkmysql ) { + printf("- ERROR: problems querying MySQL database! -\n- ERROR: Empty table? -\n\n"); + return(1); + } else { + if ( mystruct.totime == -1 ){ + printf("- ERROR: problems querying MySQL database! -\n- ERROR: Corrupted table? -\n\n"); + return(2); + }; + calib.time[s][0] = mystruct.totime; + calib.time[s][1] = mystruct.fromtime; + calib.ttime[s][0] = mystruct.calibtime; + const char *caliba = calib.basepath; + const char *yodala = calib.yodalev; + yodala2 = (char *)yodala; + const char *fila = mystruct.calcalibfile; + myfile.str(""); + myfile << caliba; + if ( !calib.DW ) { + myfile << fila; + } else { + myfile << "DW"; + string myfil = (const char*)fila; + Int_t myposiz = myfil.find("dw_"); + TString mytmp; + Int_t mspos = myposiz+2; + Int_t mepos = myposiz+13; + stringcopy(mytmp,fila,mspos,mepos); + const char *myf = mytmp; + myfile << mytmp; + }; + myfile << yodala << "/"; + calib.obtjump = mystruct.obtjump; + printf(" - event at time %i. From time %i to time %i \n use calibration at time %i, file %s\n",etime,calib.time[s][1],calib.time[s][0],calib.ttime[s][0],myfile.str().c_str()); + if ( calib.obtjump ) printf("\n WARNING: jump in the OBT! this calibration might not be the best! \n"); + if ( etime > calib.time[s][0] && !calib.obtjump ){ + printf("\n WARNING: event at time greater than the upper limit for this calibration!\nSearch for the correct calibration\n"); + goto tryagain; + }; + }; + Int_t pedeerr = 0; + pedeerr = CaloPede(myfile.str().c_str(),s,calib.ttime[s][0],calib); + if ( pedeerr ) { + printf("\n\nERROR: problems opening calibration file! \n\nERROR: YODA version of the calibration file is not %s? \n",yodala2); + return(1); + }; + printf("\n"); + }; + }; + }; + } else { // - // Define variables - // - struct CaLevel1 clevel1; - clevel1.xalig = 120.4; - clevel1.yalig = 118.6; - clevel1.zalig = -260.; - struct Mystruct mystruct; - clevel1.paw = calib.ispaw; - Int_t upnn ; - Int_t done = 0; - Int_t pre = -1; - Int_t etime; - Int_t S3 = 0; - Int_t S2 = 0; - Int_t S12 = 0; - Int_t S11 = 0; - Int_t se = 5; - bool isCOMP = 0; - bool isFULL = 0; - bool isRAW = 0; - Float_t ener, basel; - Float_t estripnull[96][22][2]; - // - Int_t doneb = 0; - Int_t nn; - Int_t ck = 0; - Int_t ipre = 0; - Int_t ip[3] = {0}; - Float_t base0, base1, base2; - base0 = 0.; - base1 = 0.; - base2 = 0.; - Float_t qpre[6] = {0}; - Float_t ene[96] = {0}; - // - pamela::calorimeter::CalorimeterEvent *de = 0; - pamela::PscuHeader *ph = 0; - pamela::EventHeader *eh = 0; - pamela::trigger::TriggerEvent *trig = 0; - otr->SetBranchAddress("Header", &eh); - otr->SetBranchAddress("Calorimeter.Event", &de); - otr->SetBranchAddress("Trigger.Event", &trig); - // - // other variables - // - char *yodala2; - Int_t chdone[4] = {0,0,0,0}; - Int_t pe = 0; - Float_t tmptrigty = -1.; - // - // start processing - // - // - // read from the header of the event the OBT and event number + // for each event check that the calibration we are using are still within calibration limits, if not call the next calibration // - otr->GetEntry(ei); - ph = eh->GetPscuHeader(); - clevel2.pkt_num = ph->GetCounter(); - etime = ph->GetOrbitalTime(); - clevel2.OBT = etime; - clevel2.pro_num = (ei+1); - // - // determine who generate the trigger for this event (TOF, S4/PULSER, CALO) - // - S3 = 0; - S2 = 0; - S12 = 0; - S11 = 0; - S3 = trig->patterntrig[2]; - S2 = trig->patterntrig[3]; - S12 = trig->patterntrig[4]; - S11 = trig->patterntrig[5]; - if ( trig->patterntrig[0] ) tmptrigty = 2.; - if ( S3 || S2 || S12 || S11 ) tmptrigty = 0.; - if ( trig->patterntrig[1] & (1<<0) || (!trig->patterntrig[0] && !S3 && !S2 && !S12 && !S11) ) tmptrigty = 1.; - clevel2.trigty = tmptrigty; - // - yodala2 = ""; - if ( calib.mysql ){ - if ( !calib.obtjump ){ - for (Int_t s = 0; s < 4; s++){ - if ( etime > calib.time[s][0] ){ - tryagain: - printf(" ** SECTION %i **\n",s); - stringstream qy; - qy.str(""); - Int_t chkmysql = 0; - stringstream myfile; - const char *tabula = calib.tablename; - qy << "select * from " << tabula; - qy << " where section = " << s; - qy << " and fromtime >= " << calib.time[s][0]; - qy << " limit 1"; - chkmysql = creadmysql(calib.db,(char *)qy.str().c_str(),mystruct); - if ( chkmysql ) { - printf("- ERROR: problems querying MySQL database! -\n- ERROR: Empty table? -\n\n"); - return(1); - } else { - if ( mystruct.totime == -1 ){ - printf("- ERROR: problems querying MySQL database! -\n- ERROR: Corrupted table? -\n\n"); - return(2); - }; - calib.time[s][0] = mystruct.totime; - calib.time[s][1] = mystruct.fromtime; - calib.ttime[s][0] = mystruct.calibtime; - const char *caliba = calib.basepath; - const char *yodala = calib.yodalev; - yodala2 = (char *)yodala; - const char *fila = mystruct.calcalibfile; - myfile.str(""); - myfile << caliba; - if ( !calib.DW ) { - myfile << fila; - } else { - myfile << "DW"; - string myfil = (const char*)fila; - Int_t myposiz = myfil.find("dw_"); - TString mytmp; - Int_t mspos = myposiz+2; - Int_t mepos = myposiz+13; - stringcopy(mytmp,fila,mspos,mepos); - const char *myf = mytmp; - myfile << mytmp; - }; - myfile << yodala << "/"; - calib.obtjump = mystruct.obtjump; - printf(" - event at time %i. From time %i to time %i \n use calibration at time %i, file %s\n",etime,calib.time[s][1],calib.time[s][0],calib.ttime[s][0],myfile.str().c_str()); - if ( calib.obtjump ) printf("\n WARNING: jump in the OBT! this calibration might not be the best! \n"); - if ( etime > calib.time[s][0] && !calib.obtjump ){ - printf("\n WARNING: event at time greater than the upper limit for this calibration!\nSearch for the correct calibration\n"); - goto tryagain; - }; - }; - Int_t pedeerr = 0; - pedeerr = CaloPede(myfile.str().c_str(),s,calib.ttime[s][0],calib); - if ( pedeerr ) { - printf("\n\nERROR: problems opening calibration file! \n\nERROR: YODA version of the calibration file is not %s? \n",yodala2); - return(1); - }; - printf("\n"); - }; + if ( !calib.obtjump ) { + for (Int_t s = 0; s < 4; s++){ + if ( calib.ttime[s][b[s]] ){ + while ( etime > calib.time[s][b[s]+1] ){ + printf(" CALORIMETER: \n" ); + printf(" - Section %i, event at time %i while old calibration time limit at %i. Use new calibration at time %i -\n",s,etime,calib.time[s][b[s]+1],calib.ttime[s][b[s]+1]); + printf(" END CALORIMETER. \n\n" ); + b[s]++; + TString pfile; + if ( calib.fcode[s][b[s]] != 10 ){ + TString file2f = ""; + stringcopy(file2f,calcalibfile); + TString pfile = whatnamewith(file2f,calib.fcode[s][b[s]]); + } else { + pfile = (TString)calcalibfile; }; + CaloPede(pfile,s,calib.ttime[s][b[s]],calib); + }; + }; + }; + }; + }; + // + // run over views and planes + // + for (Int_t l = 0; l < 2; l++){ + for (Int_t m = 0; m < 22; m++){ + // + // determine the section number + // + se = 5; + if (l == 0 && m%2 == 0) se = 3; + if (l == 0 && m%2 != 0) se = 2; + if (l == 1 && m%2 == 0) se = 1; + if (l == 1 && m%2 != 0) se = 0; + // + // determine what kind of event we are going to analyze + // + isCOMP = 0; + isFULL = 0; + isRAW = 0; + if ( de->stwerr[se] & (1 << 16) ) isCOMP = 1; + if ( de->stwerr[se] & (1 << 17) ) isFULL = 1; + if ( de->stwerr[se] & (1 << 3) ) isRAW = 1; + if ( !chdone[se] ){ + // + // check for any error in the event + // + clevel2.crc[se] = 0; + if ( de->perror[se] == 132 ){ + clevel2.crc[se] = 1; + pe++; + }; + clevel2.perr[se] = 0; + if ( de->perror[se] != 0 ){ + clevel2.perr[se] = 1; + pe++; + }; + clevel2.swerr[se] = 0; + for (Int_t j = 0; j < 7 ; j++){ + if ( (j != 3) && (de->stwerr[se] & (1 << j)) ){ + clevel2.swerr[se] = 1; + pe++; + }; + }; + chdone[se] = 1; + }; + if ( clevel2.crc[se] == 0 && (calib.good2 == 1 || clevel2.trigty == 2)){ + pre = -1; + if ( isRAW ){ + for (Int_t nn = 0; nn < 96; nn++){ + evento.dexy[l][m][nn] = de->dexy[l][m][nn] ; + evento.dexyc[l][m][nn] = de->dexyc[l][m][nn] ; + }; }; - } else { // - // for each event check that the calibration we are using are still within calibration limits, if not call the next calibration + // run over preamplifiers // - if ( !calib.obtjump ) { - for (Int_t s = 0; s < 4; s++){ - if ( calib.ttime[s][b[s]] ){ - while ( etime > calib.time[s][b[s]+1] ){ - printf(" CALORIMETER: \n" ); - printf(" - Section %i, event at time %i while old calibration time limit at %i. Use new calibration at time %i -\n",s,etime,calib.time[s][b[s]+1],calib.ttime[s][b[s]+1]); - printf(" END CALORIMETER. \n\n" ); - b[s]++; - TString pfile; - if ( calib.fcode[s][b[s]] != 10 ){ - TString file2f = ""; - stringcopy(file2f,calcalibfile); - TString pfile = whatnamewith(file2f,calib.fcode[s][b[s]]); - } else { - pfile = (TString)calcalibfile; - }; - CaloPede(pfile,s,calib.ttime[s][b[s]],calib); - }; - }; + pre = -1; + for (Int_t i = 0; i < 3; i++){ + for (Int_t j = 0; j < 2; j++){ + pre = j + i*2; + // + // baseline check and calculation + // + if ( !isRAW ) { + evento.base[l][m][pre] = de->base[l][m][pre] ; + } else { + // + // if it is a raw event and we haven't checked + // yet, calculate the baseline. + // + CaloFindBaseRaw(calib,evento,l,m,pre); }; + }; }; - }; - // - // run over views and planes - // - for (Int_t l = 0; l < 2; l++){ - for (Int_t m = 0; m < 22; m++){ - // - // determine the section number + // + // run over strips + // + pre = -1; + for (Int_t i = 0 ; i < 3 ; i++){ + // ck = 0; + ip[i] = 0; + for (Int_t n = i*32 ; n < (i+1)*32 ; n++){ + if (n%16 == 0) { + done = 0; + doneb = 0; + donec = 0; + ck = 0; + pre++; + qpre[pre] = 0; + }; + // + // baseline check and calculation // - se = 5; - if (l == 0 && m%2 == 0) se = 3; - if (l == 0 && m%2 != 0) se = 2; - if (l == 1 && m%2 == 0) se = 1; - if (l == 1 && m%2 != 0) se = 0; + if ( !isRAW ) { + evento.dexy[l][m][n] = de->dexy[l][m][n] ; + evento.dexyc[l][m][n] = de->dexyc[l][m][n] ; + }; // - // determine what kind of event we are going to analyze + // no suitable new baseline, use old ones and compress data! // - isCOMP = 0; - isFULL = 0; - isRAW = 0; - if ( de->stwerr[se] & (1 << 16) ) isCOMP = 1; - if ( de->stwerr[se] & (1 << 17) ) isFULL = 1; - if ( de->stwerr[se] & (1 << 3) ) isRAW = 1; - if ( !chdone[se] ){ - // - // check for any error in the event - // - clevel2.crc[se] = 0; - if ( de->perror[se] == 132 ){ - clevel2.crc[se] = 1; - pe++; + if ( !done ){ + if ( (evento.base[l][m][pre] == 31000. || evento.base[l][m][pre] == 0.) ){ + ck = 1; + if (pre%2 == 0) { + ip[i] = pre + 1; + } else { + ip[i] = pre - 1; }; - clevel2.perr[se] = 0; - if ( de->perror[se] != 0 ){ - clevel2.perr[se] = 1; - pe++; - }; - clevel2.swerr[se] = 0; - for (Int_t j = 0; j < 7 ; j++){ - if ( (j != 3) && (de->stwerr[se] & (1 << j)) ){ - clevel2.swerr[se] = 1; - pe++; - }; - }; - chdone[se] = 1; - }; - if ( clevel2.crc[se] == 0 && (calib.good2 == 1 || clevel2.trigty == 2)){ - pre = -1; - if ( isRAW ){ - for (Int_t nn = 0; nn < 96; nn++){ - evento.dexy[l][m][nn] = de->dexy[l][m][nn] ; - evento.dexyc[l][m][nn] = de->dexyc[l][m][nn] ; + if ( (evento.base[l][m][ip[i]] == 31000. || evento.base[l][m][ip[i]] == 0.) ){ + // + ck = 2; + if ( calib.sbase[l][m][pre] == 31000. || calib.sbase[l][m][pre] == 0. ) { + ck = 3; }; - }; - // - // run over preamplifiers - // - pre = -1; - for (Int_t i = 0; i < 3; i++){ - for (Int_t j = 0; j < 2; j++){ - pre = j + i*2; - // - // baseline check and calculation - // - if ( !isRAW ) { - evento.base[l][m][pre] = de->base[l][m][pre] ; - } else { - // - // if it is a raw event and we haven't checked - // yet, calculate the baseline. - // - CaloFindBaseRaw(calib,evento,l,m,pre); - }; - }; }; + // upnn = n+16; + //if ( upnn > 96 ) nn = 96; + //for ( Int_t nn = n; nndexyc[l][m][nn]; + //}; + //CaloCompressData(calib,evento,l,m,pre); + done = 1; + }; + }; + // + // CALIBRATION ALGORITHM + // + if ( !doneb ){ + switch (ck) { + case 0: + base0 = evento.base[l][m][pre]; + base2 = calib.calbase[l][m][pre]; + break; + case 1: + base0 = evento.base[l][m][ip[i]]; + base2 = calib.calbase[l][m][ip[i]]; + break; + case 2: + base0 = calib.sbase[l][m][pre]; + base2 = calib.calbase[l][m][pre]; + break; + case 3: + base0 = calib.calbase[l][m][pre]; + base2 = calib.calbase[l][m][pre]; + break; + }; + base1 = calib.calbase[l][m][pre]; + doneb = 1; + }; + ener = evento.dexyc[l][m][n]; + clevel1.estrip[n][m][l] = 0.; + if ( base0>0 && base0 < 30000. ){ + if ( !donec && (base0 - base1 + base2) != 0. ){ + calib.sbase[l][m][pre] = base0 - base1 + base2; + donec = 1; + }; + if ( ener > 0. ){ + clevel1.estrip[n][m][l] = (ener - calib.calped[l][m][n] - base0 - base1 + base2)/calib.mip[l][m][n] ; // - // run over strips - // - pre = -1; - for (Int_t i = 0 ; i < 3 ; i++){ - ck = 0; - ip[i] = 0; - for (Int_t n = i*32 ; n < (i+1)*32 ; n++){ - if (n%16 == 0) { - done = 0; - doneb = 0; - pre++; - qpre[pre] = 0; - }; - // - // baseline check and calculation - // - if ( !isRAW ) { - evento.dexy[l][m][n] = de->dexy[l][m][n] ; - evento.dexyc[l][m][n] = de->dexyc[l][m][n] ; - }; - // - // no suitable new baseline, use old ones and compress data! - // - if ( !done ){ - if ( (evento.base[l][m][pre] == 31000. || evento.base[l][m][pre] == 0.) ){ - ck = 1; - if (pre%2 == 0) - ip[i] = pre + 1; - else - ip[i] = pre - 1; - if ( (evento.base[l][m][ip[i]] == 31000. || evento.base[l][m][ip[i]] == 0.) ){ - // - evento.base[l][m][pre] = calib.sbase[l][m][pre]; - // - ck = 2; - upnn = n+16; - if ( upnn > 96 ) nn = 96; - for ( Int_t nn = n; nndexyc[l][m][nn]; - }; - CaloCompressData(calib,evento,l,m,pre); - done = 1; - }; - }; - }; - // - // CALIBRATION ALGORITHM - // - basel = evento.base[l][m][pre]; - if ( !doneb ){ - switch (ck) { - case 0: - base0 = evento.base[l][m][pre]; - base2 = calib.calbase[l][m][pre]; - break; - case 1: - base0 = evento.base[l][m][ip[i]]; - base2 = calib.calbase[l][m][ip[i]]; - break; - case 2: - base0 = calib.sbase[l][m][pre]; - base2 = calib.calbase[l][m][pre]; - break; - }; - base1 = calib.calbase[l][m][pre]; - doneb = 1; - }; - ener = evento.dexyc[l][m][n]; - clevel1.estrip[n][m][l] = 0.; - if ( base0>0 && base0 < 30000. && ener > 0. ){ - clevel1.estrip[n][m][l] = (ener - calib.calped[l][m][n] - base0 - base1 + base2)/calib.mip[l][m][n] ; - // - // OK, now in estrip we have the energy deposit in MIP of all the strips for this event (at the end of loops of course) - // - // if ( clevel1.estrip[n][m][l] < evento.emin || calib.calgood[l][m][n] != 0 ) { - // if ( clevel1.estrip[n][m][l] < evento.emin || calib.calrms[l][m][n] > 26 ) { - // clevel1.estrip[n][m][l] = 0.; - // }; - qpre[pre] += clevel1.estrip[n][m][l]; - }; - calib.sbase[l][m][pre] = evento.base[l][m][pre]; - }; - if (ck == 1){ - if (ip[i]%2 == 0) - ipre = ip[i] + 1; - else - ipre = ip[i] - 1; - for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){ - clevel1.estrip[j][m][l] = clevel1.estrip[j][m][l] + (qpre[ipre] - qpre[ip[i]]) * 0.00478; - }; - }; - if (ck == 2){ - for (Int_t j = i*32 ; j < (i+1)*32 ; j++){ - ipre = j/16 + 1; - clevel1.estrip[j][m][l] = clevel1.estrip[j][m][l] + qpre[ipre] * 0.00478; - }; - }; - }; + // OK, now in estrip we have the energy deposit in MIP of all the strips for this event (at the end of loops of course) // - Int_t j4 = -4; - Int_t jjj = -3; - Int_t jj = -2; - for (Int_t j = 0 ; j < 100 ; j++){ - jj++; - jjj++; - j4++; - if ( j < 96 ) ene[j] = clevel1.estrip[j][m][l]; - if ( jj >= 0 && jj < 96 ){ - if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] = ene[jj-1] - clevel1.estrip[jj][m][l] * 0.01581; - if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] = ene[jj+1] - clevel1.estrip[jj][m][l] * 0.01581; - }; - if ( jjj >= 0 && jjj < 96 ){ - if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1.estrip[jjj-1][m][l] = clevel1.estrip[jjj-1][m][l] - ene[jjj] * 0.01581; - if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1.estrip[jjj+1][m][l] = clevel1.estrip[jjj+1][m][l] - ene[jjj] * 0.01581; - }; - if ( j4 >= 0 && j4 < 96 ){ - if ( clevel1.estrip[j4][m][l]!=0. && ( clevel1.estrip[j4][m][l] < evento.emin || calib.calrms[l][m][j4] > 26 )){ - clevel1.estrip[j4][m][l] = 0.; - }; - }; - }; - // - } else { - for (Int_t nn = 0; nn < 96; nn++){ - clevel1.estrip[nn][m][l] = 0.; - }; + qpre[pre] += clevel1.estrip[n][m][l]; + }; + }; + }; + if (ck == 1){ + if (ip[i]%2 == 0) + ipre = ip[i] + 1; + else + ipre = ip[i] - 1; + for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){ + clevel1.estrip[j][m][l] = clevel1.estrip[j][m][l] + (qpre[ipre] - qpre[ip[i]]) * 0.00478; }; + }; + if (ck == 2){ + for (Int_t j = i*32 ; j < (i+1)*32 ; j++){ + ipre = j/16 + 1; + clevel1.estrip[j][m][l] = clevel1.estrip[j][m][l] + qpre[ipre] * 0.00478; + }; + }; + }; + // + Int_t j4 = -4; + Int_t jjj = -3; + Int_t jj = -2; + for (Int_t j = 0 ; j < 100 ; j++){ + jj++; + jjj++; + j4++; + if ( j < 96 ) ene[j] = clevel1.estrip[j][m][l]; + if ( jj >= 0 && jj < 96 ){ + if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] = ene[jj-1] - clevel1.estrip[jj][m][l] * 0.01581; + if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] = ene[jj+1] - clevel1.estrip[jj][m][l] * 0.01581; + }; + if ( jjj >= 0 && jjj < 96 ){ + if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1.estrip[jjj-1][m][l] = clevel1.estrip[jjj-1][m][l] - ene[jjj] * 0.01581; + if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1.estrip[jjj+1][m][l] = clevel1.estrip[jjj+1][m][l] - ene[jjj] * 0.01581; + }; + if ( j4 >= 0 && j4 < 96 ){ + if ( clevel1.estrip[j4][m][l]!=0. && ( clevel1.estrip[j4][m][l] < evento.emin || calib.calrms[l][m][j4] > 26 )){ + clevel1.estrip[j4][m][l] = 0.; + }; + }; + }; + // + } else { + for (Int_t nn = 0; nn < 96; nn++){ + clevel1.estrip[nn][m][l] = 0.; }; - }; - // - // - // - if ( !pe && calib.good2 ){ - clevel2.good = 1; - } else { - clevel2.good = 0; - }; - // - // all we need from the tracker is the state vector al_p and the boolean for good events: - // - clevel1.trkchi2 = calib.trkchi2; - if ( calib.good2 ){ - for (Int_t e = 0; e<5 ; e++){ - clevel1.al_p[e][0] = calib.al_p[e][0]; - clevel1.al_p[e][1] = calib.al_p[e][1]; - }; - clevel1.good2 = 1; - } else { - for (Int_t e = 0; e<5 ; e++){ - clevel1.al_p[e][0] = 0.; - clevel1.al_p[e][1] = 0.; }; - clevel1.good2 = 0; }; - // - // for each event generate level2 data (clevel2) calling the fortran routine and passing calibrated data (clevel1): - // - cfillcalol2(clevel1,clevel2); - // - // save in the class the level2 calorimeter informations - // - // - // fill the rootple if not in paw mode - // - if ( clevel1.paw == 0. ) tree->Fill(); - // - // clear estrip values. - // - memcpy(clevel1.estrip, estripnull, sizeof(estripnull)); - // - return(0); + }; + // + // + // + if ( !pe && calib.good2 ){ + clevel2.good = 1; + } else { + clevel2.good = 0; + }; + // + // all we need from the tracker is the state vector al_p and the boolean for good events: + // + clevel1.trkchi2 = calib.trkchi2; + if ( calib.good2 ){ + for (Int_t e = 0; e<5 ; e++){ + clevel1.al_p[e][0] = calib.al_p[e][0]; + clevel1.al_p[e][1] = calib.al_p[e][1]; + }; + clevel1.good2 = 1; + } else { + for (Int_t e = 0; e<5 ; e++){ + clevel1.al_p[e][0] = 0.; + clevel1.al_p[e][1] = 0.; + }; + clevel1.good2 = 0; + }; + // + // for each event generate level2 data (clevel2) calling the fortran routine and passing calibrated data (clevel1): + // + cfillcalol2(clevel1,clevel2); + // + // save in the class the level2 calorimeter informations + // + // + // fill the rootple if not in paw mode + // + if ( clevel1.paw == 0. ) tree->Fill(); + // + // clear estrip values. + // + memcpy(clevel1.estrip, estripnull, sizeof(estripnull)); + // + return(0); } // short int CaloLEVEL2(TString filename, TString TrackerDir="", TString outDir ="", TString Framework = "root", Int_t FORCE = 0){ - // - // define variables - // - Int_t swcode = 04; - Int_t swtrkcode = 2; - Int_t MySQL = 1; - Int_t debug = 0; - // + // + // define variables + // + Int_t swcode = 04; + Int_t swtrkcode = 0; + Int_t MySQL = 1; + Int_t debug = 0; + // #if defined(__CINT__) - if ( !gSystem->Getenv("PAM_INC") || !gSystem->Getenv("PAM_LIB") || !gSystem->Getenv("PAM_CALIB") ){ - printf("\n ERROR: you first must set the PAMELA enviromental variables\n read the README file \n"); - printf("\n Exiting... \n"); - return(11); - }; - emicheckLib(); + if ( !gSystem->Getenv("PAM_INC") || !gSystem->Getenv("PAM_LIB") || !gSystem->Getenv("PAM_CALIB") ){ + printf("\n ERROR: you first must set the PAMELA enviromental variables\n read the README file \n"); + printf("\n Exiting... \n"); + return(11); + }; + emicheckLib(); #endif - Int_t trkdir = 1; - TString trackerdir; - if ( TrackerDir=="" ){ - trackerdir = filename; - trkdir = 0; - } else { - trackerdir = TrackerDir; - }; - TString calcalibfile; - // - // other variables - // - TFile *headerFile; - TFile *caloFile; - TFile *triggerFile; - TString filety; - stringstream file2; - stringstream file3; - const char *file4; - Int_t okcalo = 0; - TTree *ctree = 0; - TFile *chfile; - Long64_t cnevents=0; - stringstream qy; - qy.str(""); - Int_t chkmysql = 0; - Int_t imtrack; - // - // load libraries to determine the track information. - // + Int_t trkdir = 1; + TString trackerdir; + if ( TrackerDir=="" ){ + trackerdir = filename; + trkdir = 0; + } else { + trackerdir = TrackerDir; + }; + TString calcalibfile; + // + // other variables + // + TFile *headerFile; + TFile *caloFile; + TFile *triggerFile; + TString filety; + stringstream file2; + stringstream file3; + const char *file4; + Int_t okcalo = 0; + TTree *ctree = 0; + TFile *chfile; + Long64_t cnevents=0; + stringstream qy; + qy.str(""); + Int_t chkmysql = 0; + Int_t imtrack; + // + // load libraries to determine the track information. + // #if defined (__CINT__) - const char *pam_lib=gSystem->Getenv("PAM_LIB"); - stringstream libload; - libload.str(""); - libload << pam_lib << "/libtrack"; - gSystem->Load(libload.str().c_str()); + const char *pam_lib=gSystem->Getenv("PAM_LIB"); + stringstream libload; + libload.str(""); + libload << pam_lib << "/libtrack"; + gSystem->Load(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libreadb2maps"; + gSystem->Load(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libreadb2maps_C"; + gSystem->Load(libload.str().c_str()); + // + // open libraries to extract calorimeter level2 data + // + libload.str(""); + libload << pam_lib << "/libcrcalol2"; + gSystem->Load(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libcrcalol2_C"; + gSystem->Load(libload.str().c_str()); +#endif + const char *franame = Framework; +#if defined (__CINT__) + if ( !strcmp(franame,"paw") ) { libload.str(""); - libload << pam_lib << "/libreadb2maps"; + libload << pam_lib << "/libopcalol2"; gSystem->Load(libload.str().c_str()); libload.str(""); - libload << pam_lib << "/libreadb2maps_C"; + libload << pam_lib << "/libopcalol2_C"; gSystem->Load(libload.str().c_str()); - // - // open libraries to extract calorimeter level2 data - // libload.str(""); - libload << pam_lib << "/libcrcalol2"; + libload << pam_lib << "/libclcalol2"; gSystem->Load(libload.str().c_str()); libload.str(""); - libload << pam_lib << "/libcrcalol2_C"; + libload << pam_lib << "/libclcalol2_C"; gSystem->Load(libload.str().c_str()); + }; #endif - const char *franame = Framework; -#if defined (__CINT__) - if ( !strcmp(franame,"paw") ) { - libload.str(""); - libload << pam_lib << "/libopcalol2"; - gSystem->Load(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libopcalol2_C"; - gSystem->Load(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libclcalol2"; - gSystem->Load(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libclcalol2_C"; - gSystem->Load(libload.str().c_str()); - }; -#endif - if ( MySQL != 1 ){ - printf("\nWARNING: you MUST run this script using the MySQL database. Reverting MySQL to 1 \n\n"); - MySQL = 1; - }; - //MySQL = 0; - gDirectory->GetList()->Delete(); - const char* startingdir = gSystem->WorkingDirectory(); - TString path; - stringcopy(path,startingdir); - // + if ( MySQL != 1 ){ + printf("\nWARNING: you MUST run this script using the MySQL database. Reverting MySQL to 1 \n\n"); + MySQL = 1; + }; + //MySQL = 0; + gDirectory->GetList()->Delete(); + const char* startingdir = gSystem->WorkingDirectory(); + TString path; + stringcopy(path,startingdir); + // #if defined (__CINT__) - Int_t MySQLload = 0; + Int_t MySQLload = 0; #endif - // - // Open files - // - headerFile = emigetFile(filename, "Physics", "Header"); - if ( !headerFile ){ - printf("\n\n ERROR: no header file, exiting...\n\n"); + // + // Open files + // + headerFile = emigetFile(filename, "Physics", "Header"); + if ( !headerFile ){ + printf("\n\n ERROR: no header file, exiting...\n\n"); #if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQLload); + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQLload); #endif - return(1); - }; - caloFile = emigetFile(filename, "Calorimeter"); - if ( !caloFile ){ - printf("\n\n ERROR: no calorimeter file! \n\n"); + return(1); + }; + caloFile = emigetFile(filename, "Calorimeter"); + if ( !caloFile ){ + printf("\n\n ERROR: no calorimeter file! \n\n"); #if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQLload); + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQLload); #endif - return(2); - }; - // - TTree *otr = (TTree*)headerFile->Get("Pscu"); - otr->AddFriend("Calorimeter",caloFile); - triggerFile = emigetFile(filename, "Trigger"); - if ( !triggerFile ){ - printf("\n\n ERROR: No trigger file! \n"); + return(2); + }; + // + TTree *otr = (TTree*)headerFile->Get("Pscu"); + otr->AddFriend("Calorimeter",caloFile); + triggerFile = emigetFile(filename, "Trigger"); + if ( !triggerFile ){ + printf("\n\n ERROR: No trigger file! \n"); #if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQLload); + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQLload); #endif - return(1); - } else { - otr->AddFriend("Trigger",triggerFile); - }; - // - // Define some basic variables - // - struct Evento evento; - struct Calib calib; - struct Mystruct mystruct; - evento.emin = 0.7; - Long64_t nevents = otr->GetEntries(); - Int_t snevents = (Int_t)nevents; - if ( nevents < 1 ) { - printf(" ERROR: the file is empty!\n\n"); + return(1); + } else { + otr->AddFriend("Trigger",triggerFile); + }; + // + // Define some basic variables + // + struct Evento evento; + struct Calib calib; + struct Mystruct mystruct; + evento.emin = 0.7; + Long64_t nevents = otr->GetEntries(); + Int_t snevents = (Int_t)nevents; + if ( nevents < 1 ) { + printf(" ERROR: the file is empty!\n\n"); +#if defined (__CINT__) + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQLload); +#endif + return(7); + }; + Int_t b[4]; + Int_t paw = 0; + // + const string myfil = (const char *)filename; + Int_t myposiz = myfil.find("dw_"); + calib.DW = 0; + if ( myposiz == -1 ) { + myposiz = myfil.find("DW_"); + calib.DW = 1; + }; + if ( myposiz == -1 ) return(6); + stringcopy(calib.basepath,filename,0,myposiz); + stringcopy(calib.yodalev,filename,myposiz+13,myposiz+15); + // + // calibration filename + // + calcalibfile = filename; + // + // create the ntuple level2 name depending on framework and starting from the YODA name + // + if ( !strcmp(franame,"paw") ) paw = 1; + // + if ( paw ){ + calib.ispaw = 1.; + filety = "rz"; + } else { + calib.ispaw = 0.; + filety = "root"; + }; + // + TString numb("2"); // level + TString detc("Calorimeter"); // detector + char *file = getLEVname(filename,detc,numb,filety); + stringstream truefile; + truefile.str(""); + truefile << file; + // + if ( outDir == "" ){ + file4 = filename; + file3.str(""); + file3 << file4 << "/Physics/Level2"; + } else { + file4 = outDir; + file3.str(""); + file3 << file4; + } + // + file2.str(""); + file2 << file3.str().c_str() << "/"; + file2 << truefile.str().c_str(); + printf("\n Filename will be: \n %s \n\n",file2.str().c_str()); + // + // check if Level2 directory exists, if not create it. + // + Int_t ERR = gSystem->MakeDirectory(file3.str().c_str()); + // + if ( !ERR ) { + printf(" LEVEL2 directory doesn't exist! Creating LEVEL2 directory \n\n"); + } else { + // + // directory exists, check if file exists (if we are not in the force mode) + // + if ( !FORCE ) { + printf(" Not in FORCE mode, check the existence of LEVEL2 data: \n\n"); + Int_t nofile = 0; + if ( paw ){ + ifstream mypawfile; + mypawfile.open(file2.str().c_str()); + if ( mypawfile ){ + nofile = 1; + } else { + printf("Error in opening file: file %s does not exist \n",file2.str().c_str()); + }; + } else { + TFile tfile(file2.str().c_str()); + if ( !tfile.IsZombie() ) nofile = 1; + }; + if ( nofile ){ + printf(" ERROR: file already exists! Use FORCE = 1 to override \n\n"); #if defined (__CINT__) gSystem->ChangeDirectory(path); unloadf77lib(Framework,MySQLload); #endif - return(7); - }; - Int_t b[4]; - Int_t paw = 0; - // - const string myfil = (const char *)filename; - Int_t myposiz = myfil.find("dw_"); - calib.DW = 0; - if ( myposiz == -1 ) { - myposiz = myfil.find("DW_"); - calib.DW = 1; - }; - if ( myposiz == -1 ) return(6); - stringcopy(calib.basepath,filename,0,myposiz); - stringcopy(calib.yodalev,filename,myposiz+13,myposiz+15); - // - // calibration filename - // - calcalibfile = filename; - // - // create the ntuple level2 name depending on framework and starting from the YODA name - // - if ( !strcmp(franame,"paw") ) paw = 1; - // - if ( paw ){ - calib.ispaw = 1.; - filety = "rz"; - } else { - calib.ispaw = 0.; - filety = "root"; + return(3); + } else { + printf("\n OK, I will create it!\n\n"); + }; }; - // - TString numb("2"); // level - TString detc("Calorimeter"); // detector - char *file = getLEVname(filename,detc,numb,filety); - stringstream truefile; - truefile.str(""); - truefile << file; - // - if ( outDir == "" ){ - file4 = filename; - file3.str(""); - file3 << file4 << "/Physics/Level2"; + }; + // + // clear calibration time intervals + // + for (Int_t s=0; s<4;s++){ + for (Int_t d = 0; d<51; d++){ + calib.ttime[s][d] = 0 ; + calib.time[s][d] = 0 ; + }; + }; + // + // open and read the rootple containing the ADC to MIP conversion values + // + stringstream calfile; + calfile.str(""); + const char *pam_calib = pathtocalibration(); + calfile << pam_calib << "/CaloADC2MIP.root"; + // + printf(" Using calorimeter calibration file: \n %s \n",calfile.str().c_str()); + chfile = new TFile(calfile.str().c_str(),"READ","Calorimeter CALIBRATION data"); + CalorimeterCalibration *ccalo = 0; + if ( chfile->IsZombie() ){ + if ( !FORCE ){ + printf("\n ERROR: no calorimeter calibration file! \n If you want to proceed using 26 as conversion factor for all strips use FORCE = 1. \n\n"); + return(10); } else { - file4 = outDir; - file3.str(""); - file3 << file4; - } - // - file2.str(""); - file2 << file3.str().c_str() << "/"; - file2 << truefile.str().c_str(); - printf("\n Filename will be: \n %s \n\n",file2.str().c_str()); - // - // check if Level2 directory exists, if not create it. - // - Int_t ERR = gSystem->MakeDirectory(file3.str().c_str()); - // - if ( !ERR ) { - printf(" LEVEL2 directory doesn't exist! Creating LEVEL2 directory \n\n"); - } else { - // - // directory exists, check if file exists (if we are not in the force mode) - // - if ( !FORCE ) { - printf(" Not in FORCE mode, check the existence of LEVEL2 data: \n\n"); - Int_t nofile = 0; - if ( paw ){ - ifstream mypawfile; - mypawfile.open(file2.str().c_str()); - if ( mypawfile ){ - nofile = 1; - } else { - printf("Error in opening file: file %s does not exist \n",file2.str().c_str()); - }; - } else { - TFile tfile(file2.str().c_str()); - if ( !tfile.IsZombie() ) nofile = 1; - }; - if ( nofile ){ - printf(" ERROR: file already exists! Use FORCE = 1 to override \n\n"); -#if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQLload); -#endif - return(3); - } else { - printf("\n OK, I will create it!\n\n"); - }; - }; + printf("\n WARNING: no calorimeter calibration file! \n WARNING: using 26 as conversion factor for all strips. \n"); + okcalo = 0; }; - // - // clear calibration time intervals - // - for (Int_t s=0; s<4;s++){ - for (Int_t d = 0; d<51; d++){ - calib.ttime[s][d] = 0 ; - calib.time[s][d] = 0 ; - }; - }; - // - // open and read the rootple containing the ADC to MIP conversion values - // - stringstream calfile; - calfile.str(""); - const char *pam_calib = pathtocalibration(); - calfile << pam_calib << "/CaloADC2MIP.root"; - // - printf(" Using calorimeter calibration file: \n %s \n",calfile.str().c_str()); - chfile = new TFile(calfile.str().c_str(),"READ","Calorimeter CALIBRATION data"); - CalorimeterCalibration *ccalo = 0; - if ( chfile->IsZombie() ){ - if ( !FORCE ){ - printf("\n ERROR: no calorimeter calibration file! \n If you want to proceed using 26 as conversion factor for all strips use FORCE = 1. \n\n"); - return(10); - } else { - printf("\n WARNING: no calorimeter calibration file! \n WARNING: using 26 as conversion factor for all strips. \n"); - okcalo = 0; + } else { + okcalo = 1; + ctree = (TTree*)chfile->Get("CaloADC"); + ctree->SetBranchAddress("Event", &ccalo); + // + cnevents = ctree->GetEntries(); + ctree->GetEntry(cnevents-1); + }; + // + if ( okcalo ) { + for (Int_t m = 0; m < 2 ; m++ ){ + for (Int_t k = 0; k < 22; k++ ){ + for (Int_t l = 0; l < 96; l++ ){ + if ( (ccalo->fp[1][m][k][l] > 20. && ccalo->fp[1][m][k][l] < 32.) || ccalo->mask[m][k][l] == 1. ) { + calib.mip[m][k][l] = ccalo->mip[m][k][l]; + if ( calib.mip[m][k][l] == 0 ) calib.mip[m][k][l] = 26. ; + } else { + calib.mip[m][k][l] = 26. ; + }; + calib.calped[m][k][l] = 0. ; + evento.dexy[m][k][l] = 0. ; + evento.dexyc[m][k][l] = 0. ; }; - } else { - okcalo = 1; - ctree = (TTree*)chfile->Get("CaloADC"); - ctree->SetBranchAddress("Event", &ccalo); - // - cnevents = ctree->GetEntries(); - ctree->GetEntry(cnevents-1); + }; }; - // - if ( okcalo ) { - for (Int_t m = 0; m < 2 ; m++ ){ - for (Int_t k = 0; k < 22; k++ ){ - for (Int_t l = 0; l < 96; l++ ){ - if ( (ccalo->fp[1][m][k][l] > 20. && ccalo->fp[1][m][k][l] < 32.) || ccalo->mask[m][k][l] == 1. ) { - calib.mip[m][k][l] = ccalo->mip[m][k][l]; - if ( calib.mip[m][k][l] == 0 ) calib.mip[m][k][l] = 26. ; - } else { - calib.mip[m][k][l] = 26. ; - }; - calib.calped[m][k][l] = 0. ; - evento.dexy[m][k][l] = 0. ; - evento.dexyc[m][k][l] = 0. ; - }; - }; - }; - } else { - for (Int_t m = 0; m < 2 ; m++ ){ - for (Int_t k = 0; k < 22; k++ ){ - for (Int_t l = 0; l < 96; l++ ){ - calib.mip[m][k][l] = 26. ; - calib.calped[m][k][l] = 0. ; - evento.dexy[m][k][l] = 0. ; - evento.dexyc[m][k][l] = 0. ; - }; - }; + } else { + for (Int_t m = 0; m < 2 ; m++ ){ + for (Int_t k = 0; k < 22; k++ ){ + for (Int_t l = 0; l < 96; l++ ){ + calib.mip[m][k][l] = 26. ; + calib.calped[m][k][l] = 0. ; + evento.dexy[m][k][l] = 0. ; + evento.dexyc[m][k][l] = 0. ; }; + }; }; - chfile->Close(); - // - // Are we able to connect to the MySQL server? - // - if ( MySQL ){ - calib.mysql = 1; - gSystem->ChangeDirectory(path); + }; + chfile->Close(); + // + // Are we able to connect to the MySQL server? + // + if ( MySQL ){ + calib.mysql = 1; + gSystem->ChangeDirectory(path); #if defined (__CINT__) - libload.str(""); - libload << pam_lib << "/libreadmy"; - gSystem->Load(libload.str().c_str()); - libload.str(""); - libload << pam_lib << "/libreadmy_C"; - gSystem->Load(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libreadmy"; + gSystem->Load(libload.str().c_str()); + libload.str(""); + libload << pam_lib << "/libreadmy_C"; + gSystem->Load(libload.str().c_str()); #endif - calib.db = "romemuons"; - string myfil = (const char*)filename; - Int_t myposiz = myfil.find("dw_"); - if ( myposiz == -1 ) myposiz = myfil.find("DW_"); - if ( myposiz == -1 ) return(6); - // - TString mytmp; - Int_t mspos = myposiz+2; - Int_t mepos = myposiz+13; - stringcopy(mytmp,filename,mspos,mepos); - const char *myfile = mytmp; - // - stringstream tablename; - tablename.str(""); - tablename << "calocalib_dw" << myfile; - TString tabula; - tabula = tablename.str().c_str(); - // - printf("\n Try the connection to the MySQL database in Trieste...\n"); - qy.str(""); - qy << "select * from " << tablename.str().c_str(); - chkmysql = creadmysql(calib.db,(char *)qy.str().c_str(),mystruct); - if ( chkmysql || mystruct.totime == -1 ) { - printf("\n- ERROR: problems querying MySQL database! -\n\n- ERROR: Cannot use the MySQL database! -\n\n"); - if ( mystruct.totime == -1 ) printf("- ERROR: it seems there are no data in table called %s -\n\n",tablename.str().c_str()); - if ( chkmysql == 1 ) printf("- ERROR: problems in MySQL initialization -\n\n"); - if ( chkmysql == 2 ) printf("- ERROR: problems in MySQL login, check username and password -\n\n"); - if ( chkmysql == 3 ) printf("- ERROR: it seems there is no table called %s -\n\n",tablename.str().c_str()); - // if ( !FORCE ){ - // printf(" Exiting! Use FORCE = 1 or MySQL = 0 to run in standalone mode \n\n"); - printf(" Exiting! Check your installation or contact calorimeter developers \n\n"); + calib.db = "romemuons"; + string myfil = (const char*)filename; + Int_t myposiz = myfil.find("dw_"); + if ( myposiz == -1 ) myposiz = myfil.find("DW_"); + if ( myposiz == -1 ) return(6); + // + TString mytmp; + Int_t mspos = myposiz+2; + Int_t mepos = myposiz+13; + stringcopy(mytmp,filename,mspos,mepos); + const char *myfile = mytmp; + // + stringstream tablename; + tablename.str(""); + tablename << "calocalib_dw" << myfile; + TString tabula; + tabula = tablename.str().c_str(); + // + printf("\n Try the connection to the MySQL database in Trieste...\n"); + qy.str(""); + qy << "select * from " << tablename.str().c_str(); + chkmysql = creadmysql(calib.db,(char *)qy.str().c_str(),mystruct); + if ( chkmysql || mystruct.totime == -1 ) { + printf("\n- ERROR: problems querying MySQL database! -\n\n- ERROR: Cannot use the MySQL database! -\n\n"); + if ( mystruct.totime == -1 ) printf("- ERROR: it seems there are no data in table called %s -\n\n",tablename.str().c_str()); + if ( chkmysql == 1 ) printf("- ERROR: problems in MySQL initialization -\n\n"); + if ( chkmysql == 2 ) printf("- ERROR: problems in MySQL login, check username and password -\n\n"); + if ( chkmysql == 3 ) printf("- ERROR: it seems there is no table called %s -\n\n",tablename.str().c_str()); + // if ( !FORCE ){ + // printf(" Exiting! Use FORCE = 1 or MySQL = 0 to run in standalone mode \n\n"); + printf(" Exiting! Check your installation or contact calorimeter developers \n\n"); #if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQL); + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQL); #endif - return(5); - //}; - // printf("WARNING: running in standalone mode! \n\n"); - // printf("WARNING: I will search for the best calibration\n\n"); - // calib.mysql = 0; - // MySQL = -1; - } else { - printf(" ...OK, the connection is fine! \n\n Using database \"%s\", table \"%s\"\n\n",calib.db,tablename.str().c_str()); - }; - calib.tablename = tablename.str().c_str(); + return(5); + //}; + // printf("WARNING: running in standalone mode! \n\n"); + // printf("WARNING: I will search for the best calibration\n\n"); + // calib.mysql = 0; + // MySQL = -1; } else { - calib.tablename = ""; - calib.mysql = 0; + printf(" ...OK, the connection is fine! \n\n Using database \"%s\", table \"%s\"\n\n",calib.db,tablename.str().c_str()); }; - if ( MySQL < 1 ){ - // - // first of all find the best calibrations for this file - // - Int_t wused = 0; - TString nfilen; - stringcopy(nfilen,file2.str().c_str()); - CaloFindCalibs(filename, calcalibfile, wused, calib); - if ( wused == 1 ) calcalibfile = filename; - // - // print on the screen the results: - // - const char *ffile = filename; - printf(" ------ %s ------- \n \n",ffile); - Int_t calibex = 0; - for (Int_t s=0; s<4;s++){ - printf(" ** SECTION %i **\n",s); - for (Int_t d = 0; d<51; d++){ - TString pfile; - if ( calib.ttime[s][d] != 0 ) { - calibex++; - if ( calib.fcode[s][d] != 10 ){ - TString file2f = ""; - stringcopy(file2f,calcalibfile); - pfile = whatnamewith(file2f,calib.fcode[s][d]); - } else { - pfile = (TString)calcalibfile; - }; - const char *ffile = pfile; - printf(" - from time %i to time %i use calibration at\n time %i, file: %s \n",calib.time[s][d],calib.time[s][d+1],calib.ttime[s][d],ffile); - if ( !strcmp(ffile,"wrong") ) calibex--; - }; - }; - printf("\n"); + calib.tablename = tablename.str().c_str(); + } else { + calib.tablename = ""; + calib.mysql = 0; + }; + if ( MySQL < 1 ){ + // + // first of all find the best calibrations for this file + // + Int_t wused = 0; + TString nfilen; + stringcopy(nfilen,file2.str().c_str()); + CaloFindCalibs(filename, calcalibfile, wused, calib); + if ( wused == 1 ) calcalibfile = filename; + // + // print on the screen the results: + // + const char *ffile = filename; + printf(" ------ %s ------- \n \n",ffile); + Int_t calibex = 0; + for (Int_t s=0; s<4;s++){ + printf(" ** SECTION %i **\n",s); + for (Int_t d = 0; d<51; d++){ + TString pfile; + if ( calib.ttime[s][d] != 0 ) { + calibex++; + if ( calib.fcode[s][d] != 10 ){ + TString file2f = ""; + stringcopy(file2f,calcalibfile); + pfile = whatnamewith(file2f,calib.fcode[s][d]); + } else { + pfile = (TString)calcalibfile; + }; + const char *ffile = pfile; + printf(" - from time %i to time %i use calibration at\n time %i, file: %s \n",calib.time[s][d],calib.time[s][d+1],calib.ttime[s][d],ffile); + if ( !strcmp(ffile,"wrong") ) calibex--; }; - printf(" ----------------------------------------------------------------------- \n \n"); + }; + printf("\n"); + }; + printf(" ----------------------------------------------------------------------- \n \n"); - // - // const char *file2 = nfilen; // the full path and name of the level1 ntuple - file2.str(""); - file2 << nfilen; - if ( calibex < 4 ) { - printf("No full calibration data in this file, sorry!\n\n "); - // - // remove the empty file! - // - stringstream remfile; - remfile.str(""); - remfile << "rm -f " << file2.str().c_str(); - gSystem->Exec(remfile.str().c_str()); + // + // const char *file2 = nfilen; // the full path and name of the level1 ntuple + file2.str(""); + file2 << nfilen; + if ( calibex < 4 ) { + printf("No full calibration data in this file, sorry!\n\n "); + // + // remove the empty file! + // + stringstream remfile; + remfile.str(""); + remfile << "rm -f " << file2.str().c_str(); + gSystem->Exec(remfile.str().c_str()); #if defined (__CINT__) - gSystem->ChangeDirectory(path); - unloadf77lib(Framework,MySQL); + gSystem->ChangeDirectory(path); + unloadf77lib(Framework,MySQL); #endif - return(4); - }; - - // - // fill with the first calibration values the common containing calibration informations - // - for (Int_t s = 0; s < 4; s++){ - b[s]=0; - TString pfile; - if ( calib.fcode[s][b[s]] != 10 ){ - TString file2f = ""; - stringcopy(file2f,calcalibfile); - pfile = whatnamewith(file2f,calib.fcode[s][b[s]]); - } else { - pfile = (TString)calcalibfile; - }; - CaloPede(pfile,s,calib.ttime[s][b[s]],calib); - }; - } else { - calib.obtjump = 0; - calib.time[0][0] = -1; - calib.time[1][0] = -1; - calib.time[2][0] = -1; - calib.time[3][0] = -1; + return(4); }; - - // - // Open file to write depending on the framework - // - // - // In any case we must create a tree - // - TTree *tree = 0; - // - // Open tracker ntuple or rootples and load magnetic field maps. - // + // - printf(" TRACKER: loading the magnetic field maps...\n\n\n"); - pam_calib = pathtocalibration(); - stringstream bdir; - bdir.str(""); - bdir << pam_calib << "."; - creadB(bdir.str().c_str()); + // fill with the first calibration values the common containing calibration informations // - printf(" ...done! \n"); + for (Int_t s = 0; s < 4; s++){ + b[s]=0; + TString pfile; + if ( calib.fcode[s][b[s]] != 10 ){ + TString file2f = ""; + stringcopy(file2f,calcalibfile); + pfile = whatnamewith(file2f,calib.fcode[s][b[s]]); + } else { + pfile = (TString)calcalibfile; + }; + CaloPede(pfile,s,calib.ttime[s][b[s]],calib); + }; + } else { + calib.obtjump = 0; + calib.time[0][0] = -1; + calib.time[1][0] = -1; + calib.time[2][0] = -1; + calib.time[3][0] = -1; + }; + + // + // Open file to write depending on the framework + // + // + // In any case we must create a tree + // + TTree *tree = 0; + // + // Open tracker ntuple or rootples and load magnetic field maps. + // + // + printf(" TRACKER: loading the magnetic field maps...\n\n\n"); + pam_calib = pathtocalibration(); + stringstream bdir; + bdir.str(""); + bdir << pam_calib << "."; + creadB(bdir.str().c_str()); + // + printf(" ...done! \n"); + // + printf("\n Check the existence of tracker data... \n"); + Int_t isrootple = 1; + struct Tracklev2 trk; + TFile *trfile = emigetFile(trackerdir,"Physics.Level2","Tracker"); + TTree *tr = 0; + Int_t trkpaw = 0; + Int_t trnevents = 0; + stringstream filetrk; + if ( !trfile ) { + const char *file0 = filename; + string fil = (const char *)filename; + Int_t posiz = fil.find("dw_"); + if ( posiz == -1 ) posiz = fil.find("DW_"); + if ( posiz == -1 ) return(5); + // + TString trntmp; + Int_t spos = posiz+3; + Int_t epos = posiz+13; + stringcopy(trntmp,filename,spos,epos); + const char *trname = trntmp; + if ( !trkdir ){ + filetrk.str(""); + filetrk << file0 << "/Physics/Level2/DW_"; + filetrk << trname << "_level2.rz"; + } else { + const char *trackdir = trackerdir; + filetrk.str(""); + filetrk << trackdir << "/DW_"; + filetrk << trname << "_level2.rz"; + }; // - printf("\n Check the existence of tracker data... \n"); - Int_t isrootple = 1; - struct Tracklev2 trk; - TFile *trfile = emigetFile(trackerdir,"Physics.Level2","Tracker"); - TTree *tr = 0; - Int_t trkpaw = 0; - Int_t trnevents = 0; - stringstream filetrk; - if ( !trfile ) { - const char *file0 = filename; - string fil = (const char *)filename; - Int_t posiz = fil.find("dw_"); - if ( posiz == -1 ) posiz = fil.find("DW_"); - if ( posiz == -1 ) return(5); - // - TString trntmp; - Int_t spos = posiz+3; - Int_t epos = posiz+13; - stringcopy(trntmp,filename,spos,epos); - const char *trname = trntmp; - if ( !trkdir ){ - filetrk.str(""); - filetrk << file0 << "/Physics/Level2/DW_"; - filetrk << trname << "_level2.rz"; - } else { - const char *trackdir = trackerdir; - filetrk.str(""); - filetrk << trackdir << "/DW_"; - filetrk << trname << "_level2.rz"; - }; - // - ifstream mypawfile; - mypawfile.open(filetrk.str().c_str()); - if ( mypawfile ){ - trkpaw = 1; - isrootple = 0; - printf(" ...found tracker level2 NTUPLE:\n %s \n\n",filetrk.str().c_str()); + ifstream mypawfile; + mypawfile.open(filetrk.str().c_str()); + if ( mypawfile ){ + trkpaw = 1; + isrootple = 0; + printf(" ...found tracker level2 NTUPLE:\n %s \n\n",filetrk.str().c_str()); #if defined (__CINT__) - const char *sdir=gSystem->Getenv("PAM_LIB"); - bdir.str(""); - bdir << sdir << "/liboptrklev2.so"; - gSystem->Load(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/liboptrklev2_C.so"; - gSystem->Load(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libretrklev2.so"; - gSystem->Load(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libretrklev2_C.so"; - gSystem->Load(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libcltrklev2.so"; - gSystem->Load(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libcltrklev2_C.so"; - gSystem->Load(bdir.str().c_str()); + const char *sdir=gSystem->Getenv("PAM_LIB"); + bdir.str(""); + bdir << sdir << "/liboptrklev2.so"; + gSystem->Load(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/liboptrklev2_C.so"; + gSystem->Load(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libretrklev2.so"; + gSystem->Load(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libretrklev2_C.so"; + gSystem->Load(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libcltrklev2.so"; + gSystem->Load(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libcltrklev2_C.so"; + gSystem->Load(bdir.str().c_str()); #endif - trnevents = 0; - coptrklev2((char *)filetrk.str().c_str(),trk,trnevents); - } else { - printf("Error in opening file: file %s does not exist \n",file2.str().c_str()); - printf("\n No tracker data! You must first process data with tracker programs\n and then you can run CaloLEVEL2. \n Exiting... \n\n"); - return(8); - }; + trnevents = 0; + coptrklev2((char *)filetrk.str().c_str(),trk,trnevents); } else { - // - const char *file0 = filename; - string fil = (const char *)filename; - Int_t posiz = fil.find("dw_"); - if ( posiz == -1 ) posiz = fil.find("DW_"); - if ( posiz == -1 ) return(5); - // - TString trntmp; - Int_t spos = posiz; - Int_t epos = posiz+15; - stringcopy(trntmp,filename,spos,epos); - const char *trname = trntmp; - // - if ( !trkdir ){ - filetrk.str(""); - filetrk << file0 << "/Physics/Level2/"; - filetrk << trname << ".Physics.Level2.Tracker.Event.root"; - } else { - const char *trkdirname = trackerdir; - filetrk.str(""); - filetrk << trkdirname << "/"; - filetrk << trname << ".Physics.Level2.Tracker.Event.root"; - }; - // - printf(" ...found tracker level2 ROOTPLE:\n %s \n",filetrk.str().c_str()); - tr = (TTree*) trfile->Get("TrkLevel2"); - settrklev2(tr,trk); - trnevents = tr->GetEntries(); - }; - if ( trnevents != nevents ){ - printf("\n WARNING: different length for tracker and calorimeter ntuples!\n"); - printf(" calorimeter: %i events \n",(int)nevents); - printf(" tracker : %i events \n\n",(int)trnevents); - }; - // - stringstream name; - TString tmpfile2; - stringcopy(tmpfile2,file2.str().c_str()); - const char *tmpname = tmpfile2; - TFile *hfile = 0; - struct CaLevel2 clevel2; - if ( paw ){ - name.str(""); - name << tmpname; - copencalol2((char *)name.str().c_str(),trkpaw,swcode,swtrkcode); + printf("Error in opening file: file %s does not exist \n",file2.str().c_str()); + printf("\n No tracker data! You must first process data with tracker programs\n and then you can run CaloLEVEL2. \n Exiting... \n\n"); + return(8); + }; + } else { + // + const char *file0 = filename; + string fil = (const char *)filename; + Int_t posiz = fil.find("dw_"); + if ( posiz == -1 ) posiz = fil.find("DW_"); + if ( posiz == -1 ) return(5); + // + TString trntmp; + Int_t spos = posiz; + Int_t epos = posiz+15; + stringcopy(trntmp,filename,spos,epos); + const char *trname = trntmp; + // + if ( !trkdir ){ + filetrk.str(""); + filetrk << file0 << "/Physics/Level2/"; + filetrk << trname << ".Physics.Level2.Tracker.Event.root"; } else { - char *type; - type = "NEW"; - if ( FORCE ) type = "RECREATE"; - hfile = new TFile(file2.str().c_str(),type,"Calorimeter LEVEL2 data"); - // - // hfile = new TFile(file2,type,"Calorimeter LEVEL2 data",0); - // - tree = new TTree("CaloLevel2","PAMELA Level2 calorimeter data"); - tree->Branch("OBT",&clevel2.OBT,"OBT/I"); - tree->Branch("pkt_num",&clevel2.pkt_num,"pkt_num/I"); - tree->Branch("pro_num",&clevel2.pro_num,"pro_num/I"); - tree->Branch("trigty",&clevel2.trigty,"trigty/F"); - tree->Branch("good",&clevel2.good,"good/I"); - tree->Branch("perr",clevel2.perr,"perr[4]/I"); - tree->Branch("swerr",clevel2.swerr,"swerr[4]/I"); - tree->Branch("crc",clevel2.crc,"crc[4]/I"); - tree->Branch("nstrip",&clevel2.nstrip,"nstrip/F"); - tree->Branch("qtot",&clevel2.qtot,"qtot/F"); - tree->Branch("ncore",&clevel2.ncore,"ncore/F"); - tree->Branch("qcore",&clevel2.qcore,"qcore/F"); - tree->Branch("impx",&clevel2.impx,"impx/F"); - tree->Branch("impy",&clevel2.impy,"impy/F"); - tree->Branch("tanx",&clevel2.tanx,"tanx/F"); - tree->Branch("tany",&clevel2.tany,"tany/F"); - tree->Branch("nint",&clevel2.nint,"nint/F"); - tree->Branch("ncyl",&clevel2.ncyl,"ncyl/F"); - tree->Branch("qcyl",&clevel2.qcyl,"qcyl/F"); - tree->Branch("qtrack",&clevel2.qtrack,"qtrack/F"); - tree->Branch("qmax",&clevel2.qmax,"qmax/F"); - tree->Branch("nx22",&clevel2.nx22,"nx22/F"); - tree->Branch("qx22",&clevel2.nx22,"qx22/F"); - tree->Branch("qq",clevel2.qq,"qq[4]/F"); - tree->Branch("qtrackx",&clevel2.qtrackx,"qtrackx/F"); - tree->Branch("qtracky",&clevel2.qtracky,"qtracky/F"); - tree->Branch("dxtrack",&clevel2.dxtrack,"dxtrack/F"); - tree->Branch("dytrack",&clevel2.dytrack,"dytrack/F"); - tree->Branch("qlast",&clevel2.qlast,"qlast/F"); - tree->Branch("nlast",&clevel2.nlast,"nlast/F"); - tree->Branch("qpre",&clevel2.qpre,"qpre/F"); - tree->Branch("npre",&clevel2.npre,"npre/F"); - tree->Branch("qpresh",&clevel2.qpresh,"qpresh/F"); - tree->Branch("npresh",&clevel2.npresh,"npresh/F"); - tree->Branch("qlow",&clevel2.qlow,"qlow/F"); - tree->Branch("nlow",&clevel2.nlow,"nlow/F"); - tree->Branch("qtr",&clevel2.qtr,"qtr/F"); - tree->Branch("ntr",&clevel2.ntr,"ntr/F"); - tree->Branch("planetot",&clevel2.planetot,"planetot/F"); - tree->Branch("qmean",&clevel2.qmean,"qmean/F"); - tree->Branch("varcfit",clevel2.varcfit,"varcfit[2]/F"); - tree->Branch("npcfit",clevel2.npcfit,"npcfit[2]/I"); - tree->Branch("thex",&clevel2.thex,"thex/F"); - tree->Branch("they",&clevel2.they,"they/F"); - tree->Branch("cibar",clevel2.cibar,"cibar[22][2]/I"); - tree->Branch("tibar",clevel2.tibar,"tibar[22][2]/I"); - tree->Branch("cbar",clevel2.cbar,"cbar[22][2]/F"); - tree->Branch("tbar",clevel2.tbar,"tbar[22][2]/F"); - // - TTree *software = 0; - software = new TTree("Software","Software used to generate data"); - software->Branch("swcode",&swcode,"swcode/I"); - software->Branch("swtrkcode",&swtrkcode,"swtrkcode/I"); - software->Fill(); - }; + const char *trkdirname = trackerdir; + filetrk.str(""); + filetrk << trkdirname << "/"; + filetrk << trname << ".Physics.Level2.Tracker.Event.root"; + }; + // + printf(" ...found tracker level2 ROOTPLE:\n %s \n",filetrk.str().c_str()); + tr = (TTree*) trfile->Get("TrkLevel2"); + settrklev2(tr,trk); + trnevents = tr->GetEntries(); + }; + if ( trnevents != nevents ){ + printf("\n WARNING: different length for tracker and calorimeter ntuples!\n"); + printf(" calorimeter: %i events \n",(int)nevents); + printf(" tracker : %i events \n\n",(int)trnevents); + }; + // + swtrkcode = trk.swcode; + // + stringstream name; + TString tmpfile2; + stringcopy(tmpfile2,file2.str().c_str()); + const char *tmpname = tmpfile2; + TFile *hfile = 0; + struct CaLevel2 clevel2; + if ( paw ){ + name.str(""); + name << tmpname; + copencalol2((char *)name.str().c_str(),trkpaw,swcode,swtrkcode); + } else { + char *type; + type = "NEW"; + if ( FORCE ) type = "RECREATE"; + hfile = new TFile(file2.str().c_str(),type,"Calorimeter LEVEL2 data"); + // + // hfile = new TFile(file2,type,"Calorimeter LEVEL2 data",0); + // + tree = new TTree("CaloLevel2","PAMELA Level2 calorimeter data"); + tree->Branch("OBT",&clevel2.OBT,"OBT/I"); + tree->Branch("pkt_num",&clevel2.pkt_num,"pkt_num/I"); + tree->Branch("pro_num",&clevel2.pro_num,"pro_num/I"); + tree->Branch("trigty",&clevel2.trigty,"trigty/F"); + tree->Branch("good",&clevel2.good,"good/I"); + tree->Branch("perr",clevel2.perr,"perr[4]/I"); + tree->Branch("swerr",clevel2.swerr,"swerr[4]/I"); + tree->Branch("crc",clevel2.crc,"crc[4]/I"); + tree->Branch("nstrip",&clevel2.nstrip,"nstrip/F"); + tree->Branch("qtot",&clevel2.qtot,"qtot/F"); + tree->Branch("ncore",&clevel2.ncore,"ncore/F"); + tree->Branch("qcore",&clevel2.qcore,"qcore/F"); + tree->Branch("impx",&clevel2.impx,"impx/F"); + tree->Branch("impy",&clevel2.impy,"impy/F"); + tree->Branch("tanx",&clevel2.tanx,"tanx/F"); + tree->Branch("tany",&clevel2.tany,"tany/F"); + tree->Branch("nint",&clevel2.nint,"nint/F"); + tree->Branch("ncyl",&clevel2.ncyl,"ncyl/F"); + tree->Branch("qcyl",&clevel2.qcyl,"qcyl/F"); + tree->Branch("qtrack",&clevel2.qtrack,"qtrack/F"); + tree->Branch("qmax",&clevel2.qmax,"qmax/F"); + tree->Branch("nx22",&clevel2.nx22,"nx22/F"); + tree->Branch("qx22",&clevel2.nx22,"qx22/F"); + tree->Branch("qq",clevel2.qq,"qq[4]/F"); + tree->Branch("qtrackx",&clevel2.qtrackx,"qtrackx/F"); + tree->Branch("qtracky",&clevel2.qtracky,"qtracky/F"); + tree->Branch("dxtrack",&clevel2.dxtrack,"dxtrack/F"); + tree->Branch("dytrack",&clevel2.dytrack,"dytrack/F"); + tree->Branch("qlast",&clevel2.qlast,"qlast/F"); + tree->Branch("nlast",&clevel2.nlast,"nlast/F"); + tree->Branch("qpre",&clevel2.qpre,"qpre/F"); + tree->Branch("npre",&clevel2.npre,"npre/F"); + tree->Branch("qpresh",&clevel2.qpresh,"qpresh/F"); + tree->Branch("npresh",&clevel2.npresh,"npresh/F"); + tree->Branch("qlow",&clevel2.qlow,"qlow/F"); + tree->Branch("nlow",&clevel2.nlow,"nlow/F"); + tree->Branch("qtr",&clevel2.qtr,"qtr/F"); + tree->Branch("ntr",&clevel2.ntr,"ntr/F"); + tree->Branch("planetot",&clevel2.planetot,"planetot/F"); + tree->Branch("qmean",&clevel2.qmean,"qmean/F"); + tree->Branch("varcfit",clevel2.varcfit,"varcfit[2]/F"); + tree->Branch("npcfit",clevel2.npcfit,"npcfit[2]/I"); + tree->Branch("thex",&clevel2.thex,"thex/F"); + tree->Branch("they",&clevel2.they,"they/F"); + tree->Branch("cibar",clevel2.cibar,"cibar[22][2]/I"); + tree->Branch("tibar",clevel2.tibar,"tibar[22][2]/I"); + tree->Branch("cbar",clevel2.cbar,"cbar[22][2]/F"); + tree->Branch("tbar",clevel2.tbar,"tbar[22][2]/F"); + // + TTree *software = 0; + software = new TTree("Software","Software used to generate data"); + software->Branch("swcode",&swcode,"swcode/I"); + software->Branch("swtrkcode",&swtrkcode,"swtrkcode/I"); + software->Fill(); + }; + // + // run over all the events + // + printf("\n Processed events: \n\n"); + // + pamela::PscuHeader *ph = 0; + pamela::EventHeader *eh = 0; + Int_t caloerr = 0; + Int_t i = -1; + Int_t itr; + if ( isrootple ){ + itr = -1; + } else { + itr = 0; + }; + Int_t retval = 0; + Int_t syncro = 1; + Int_t trklost = 0; + Int_t calolost = 0; + Int_t pktnum = 0; + Int_t obt = 0; + // + // + while ( i < (nevents-1) ){ + // + if ( i%1000 == 0 && i > 0 ) printf(" %iK \n",i/1000); + // printf(" %i \n",i); + // + // look for tracker data + // + itr++; + i++; + syncro = 1; +#if !defined (__CINT__) + trkcalosync: +#endif +#if defined (__CINT__) + trkcalosync: printf(""); +#endif // - // run over all the events + // check if we have tracker data + // + if ( i >= nevents ){ + printf(" WARNING: no more calorimeter data.\n"); + if ( nevents > trnevents ) calolost += (snevents-trnevents); + if ( nevents < trnevents ) trklost += (-snevents+trnevents); + retval = 9; + goto closeandexit; + }; + if ( itr > trnevents && isrootple ){ + printf(" WARNING: no more tracker data.\n"); + if ( nevents > trnevents ) calolost += (snevents-trnevents); + if ( nevents < trnevents ) trklost += (-snevents+trnevents); + retval = 9; + goto closeandexit; + }; + if ( itr >= trnevents+1 && !isrootple ){ + printf(" WARNING: no more tracker data.\n"); + if ( nevents > trnevents ) calolost += (snevents-trnevents); + if ( nevents < trnevents ) trklost += (-snevents+trnevents); + retval = 9; + goto closeandexit; + }; // - printf("\n Processed events: \n\n"); + // retrieve tracker informations // - pamela::PscuHeader *ph = 0; - pamela::EventHeader *eh = 0; - Int_t caloerr = 0; - Int_t i = -1; - Int_t itr; if ( isrootple ){ - itr = -1; + if ( debug ) printf(" itr %i i %i nevents %i trnevents %i \n",itr,i,(int)nevents,trnevents); + tr->GetEntry(itr); } else { - itr = 0; + if ( debug ) printf(" itr %i i %i nevents %i trnevents %i \n",itr,i,(int)nevents,trnevents); + cretrklev2(itr,trk); }; - Int_t retval = 0; - Int_t syncro = 1; - Int_t trklost = 0; - Int_t calolost = 0; - Int_t pktnum = 0; - Int_t obt = 0; // - while ( i < (nevents-1) ){ - // - if ( i%1000 == 0 && i > 0 ) printf(" %iK \n",i/1000); - // printf(" %i \n",i); - // - // look for tracker data - // + // check synchronization tracker and calorimeter informations: + // + otr->SetBranchAddress("Header", &eh); + otr->GetEntry(i); + ph = eh->GetPscuHeader(); + pktnum = ph->GetCounter(); + obt = ph->GetOrbitalTime(); + if ( pktnum != trk.pkt_num || obt != trk.obt ){ + if ( pktnum > trk.pkt_num || obt > trk.obt ){ + if ( debug ) printf("itr++ %i pktnum calo %i trk %i obt calo %i trk %i \n",itr,pktnum,trk.pkt_num,obt,trk.obt); itr++; + trklost++; + if ( syncro ) printf(" WARNING: lost sync! try to recover... \n"); + syncro = 0; + goto trkcalosync; + }; + if ( pktnum < trk.pkt_num || obt < trk.obt ){ + if ( debug ) printf("i++ %i pktnum calo %i trk %i obt calo %i trk %i \n",i,pktnum,trk.pkt_num,obt,trk.obt); i++; - syncro = 1; -#if !defined (__CINT__) - trkcalosync: -#endif -#if defined (__CINT__) - trkcalosync: printf(""); -#endif - // - // check if we have tracker data - // - if ( i >= nevents ){ - printf(" WARNING: no more calorimeter data.\n"); - if ( nevents > trnevents ) calolost += (snevents-trnevents); - if ( nevents < trnevents ) trklost += (-snevents+trnevents); - retval = 9; - goto closeandexit; - }; - if ( itr > trnevents && isrootple ){ - printf(" WARNING: no more tracker data.\n"); - if ( nevents > trnevents ) calolost += (snevents-trnevents); - if ( nevents < trnevents ) trklost += (-snevents+trnevents); - retval = 9; - goto closeandexit; - }; - if ( itr >= trnevents+1 && !isrootple ){ - printf(" WARNING: no more tracker data.\n"); - if ( nevents > trnevents ) calolost += (snevents-trnevents); - if ( nevents < trnevents ) trklost += (-snevents+trnevents); - retval = 9; - goto closeandexit; - }; - // - // retrieve tracker informations - // - if ( isrootple ){ - if ( debug ) printf(" itr %i i %i nevents %i trnevents %i \n",itr,i,(int)nevents,trnevents); - tr->GetEntry(itr); - } else { - if ( debug ) printf(" itr %i i %i nevents %i trnevents %i \n",itr,i,(int)nevents,trnevents); - cretrklev2(itr,trk); + calolost++; + if ( syncro ) printf(" WARNING: lost sync! try to recover... \n"); + syncro = 0; + goto trkcalosync; + }; + }; + // + // here we have synchronized data + // + if ( !syncro ) { + printf(" ...synchronization recovered! \n"); + printf(" Sync info: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n\n",trklost,calolost); + syncro = 1; + if ( debug ) printf("pktnum calo %i trk %i obt calo %i trk %i \n",pktnum,trk.pkt_num,obt,trk.obt); + }; + // + // store track information in the calib structure + // + if ( trk.ntrk > 0 ){ + imtrack = 0; + for (Int_t e = 0; eSetBranchAddress("Header", &eh); - otr->GetEntry(i); - ph = eh->GetPscuHeader(); - pktnum = ph->GetCounter(); - obt = ph->GetOrbitalTime(); - if ( pktnum != trk.pkt_num || obt != trk.obt ){ - if ( pktnum > trk.pkt_num || obt > trk.obt ){ - if ( debug ) printf("itr++ %i pktnum calo %i trk %i obt calo %i trk %i \n",itr,pktnum,trk.pkt_num,obt,trk.obt); - itr++; - trklost++; - if ( syncro ) printf(" WARNING: lost sync! try to recover... \n"); - syncro = 0; - goto trkcalosync; - }; - if ( pktnum < trk.pkt_num || obt < trk.obt ){ - if ( debug ) printf("i++ %i pktnum calo %i trk %i obt calo %i trk %i \n",i,pktnum,trk.pkt_num,obt,trk.obt); - i++; - calolost++; - if ( syncro ) printf(" WARNING: lost sync! try to recover... \n"); - syncro = 0; - goto trkcalosync; - }; - }; - // - // here we have synchronized data - // - if ( !syncro ) { - printf(" ...synchronization recovered! \n"); - printf(" Sync info: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n\n",trklost,calolost); - syncro = 1; - if ( debug ) printf("pktnum calo %i trk %i obt calo %i trk %i \n",pktnum,trk.pkt_num,obt,trk.obt); - }; - // - // store track information in the calib structure - // - if ( trk.ntrk > 0 ){ - imtrack = 0; - for (Int_t e = 0; etrk.chi2[1] ) { - calib.trkchi2 = 1; - } else { - calib.trkchi2 = 2; - }; - } else { - calib.trkchi2 = 1; - }; + if ( imtrack != 0 ) { + if ( trk.chi2[0]>trk.chi2[1] ) { + calib.trkchi2 = 1; } else { - calib.good2 = false; - for (Int_t e = 0; e<5 ; e++){ - calib.al_p[e][0] = 0.; - calib.al_p[e][1] = 0.; - }; + calib.trkchi2 = 2; }; } else { - calib.good2 = false; - for (Int_t e = 0; e<5 ; e++){ - calib.al_p[e][0] = 0.; - calib.al_p[e][1] = 0.; - }; + calib.trkchi2 = 1; }; - // - // calibrate calorimeter data and retrieve level2 informations - // - // ====> NOTICE: in the case of no tracks or not good tracker events the program will fill the ntuple with zero values; <==== - // ====> in case of multiple tracks the program will calculate variables using the state vector and rigidity of the first track stored <==== - // - printf("Event %i \n",i); - caloerr = calolevel2core(i,b,otr,tree,clevel2,evento,calib,calcalibfile); - // - // - // - if ( caloerr ){ - if ( i > 0 ){ - printf("\nTrying to close the file anyway...\n"); - if ( paw ){ - name.str(""); - name << tmpname; - cclosecalol2((char *)name.str().c_str()); - } else { - hfile->Write(); - hfile->Close(); - }; - printf("...done!\n"); - } else { - printf("\nERROR: output file is empty! \n"); - stringstream remfile; - remfile.str(""); - remfile << "rm -f " << file2.str().c_str(); - gSystem->Exec(remfile.str().c_str()); - }; - printf("\nERROR: processed %i out of %i events\n",i,(int)nevents); - printf("\nERROR: an error occurred during processing!\n\n Exiting...\n\n"); - goto theend; - }; - Int_t ciccio; - ciccio = 1; - }; - if ( (i+1) < trnevents ){ - printf(" WARNING: no more calorimeter data. \n"); - if ( nevents > trnevents ) calolost += (snevents-trnevents); - if ( nevents < trnevents ) trklost += (-snevents+trnevents); - printf(" Sync info: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n",trklost,calolost); - retval = 9; + } else { + calib.good2 = false; + for (Int_t e = 0; e<5 ; e++){ + calib.al_p[e][0] = 0.; + calib.al_p[e][1] = 0.; + }; + }; + } else { + calib.good2 = false; + for (Int_t e = 0; e<5 ; e++){ + calib.al_p[e][0] = 0.; + calib.al_p[e][1] = 0.; + }; }; // - // close rootple/ntuple files + // calibrate calorimeter data and retrieve level2 informations // - closeandexit: - if ( paw ){ - name.str(""); - name << tmpname; - cclosecalol2((char *)name.str().c_str()); - } else { - hfile->Write(); - hfile->Close(); - }; - if ( isrootple ){ - trfile->Close(); - } else { - ccltrklev2(trk); - }; - printf("\n"); - printf(" Finished, exiting... \n"); - printf("\n"); - if ( trklost || calolost ){ - printf("\n Sync info over the whole file: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n",trklost,calolost); - }; - printf("\n"); - theend: + // ====> NOTICE: in the case of no tracks or not good tracker events the program will fill the ntuple with zero values; <==== + // ====> in case of multiple tracks the program will calculate variables using the state vector and rigidity of the first track stored <==== // - // go back to the starting path and unload fortran libraries + caloerr = calolevel2core(i,b,otr,tree,clevel2,evento,calib,calcalibfile); // - gSystem->ChangeDirectory(path); + // + // + if ( caloerr ){ + if ( i > 0 ){ + printf("\nTrying to close the file anyway...\n"); + if ( paw ){ + name.str(""); + name << tmpname; + cclosecalol2((char *)name.str().c_str()); + } else { + hfile->Write(); + hfile->Close(); + }; + printf("...done!\n"); + } else { + printf("\nERROR: output file is empty! \n"); + stringstream remfile; + remfile.str(""); + remfile << "rm -f " << file2.str().c_str(); + gSystem->Exec(remfile.str().c_str()); + }; + printf("\nERROR: processed %i out of %i events\n",i,(int)nevents); + printf("\nERROR: an error occurred during processing!\n\n Exiting...\n\n"); + goto theend; + }; + Int_t ciccio; + ciccio = 1; + }; + if ( (i+1) < trnevents ){ + printf(" WARNING: no more calorimeter data. \n"); + if ( nevents > trnevents ) calolost += (snevents-trnevents); + if ( nevents < trnevents ) trklost += (-snevents+trnevents); + printf(" Sync info: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n",trklost,calolost); + retval = 9; + }; + // + // close rootple/ntuple files + // + closeandexit: + if ( paw ){ + name.str(""); + name << tmpname; + cclosecalol2((char *)name.str().c_str()); + } else { + hfile->Write(); + hfile->Close(); + }; + if ( isrootple ){ + trfile->Close(); + } else { + ccltrklev2(trk); + }; + printf("\n"); + printf(" Finished, exiting... \n"); + printf("\n"); + if ( trklost || calolost ){ + printf("\n Sync info over the whole file: \n - tracker packets without calorimeter %i\n - calorimeter packets without tracker %i\n",trklost,calolost); + }; + printf("\n"); + theend: + // + // go back to the starting path and unload fortran libraries + // + gSystem->ChangeDirectory(path); #if defined (__CINT__) - if ( !isrootple ){ - const char *sdir=gSystem->Getenv("PAM_LIB"); - bdir.str(""); - bdir << sdir << "/libcltrklev2_C.so"; - gSystem->Unload(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libcltrklev2.so"; - gSystem->Unload(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libretrklev2_C.so"; - gSystem->Unload(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/libretrklev2.so"; - gSystem->Unload(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/liboptrklev2_C.so"; - gSystem->Unload(bdir.str().c_str()); - bdir.str(""); - bdir << sdir << "/liboptrklev2.so"; - gSystem->Unload(bdir.str().c_str()); - }; - // unloadf77lib(Framework,MySQL); + if ( !isrootple ){ + const char *sdir=gSystem->Getenv("PAM_LIB"); + bdir.str(""); + bdir << sdir << "/libcltrklev2_C.so"; + gSystem->Unload(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libcltrklev2.so"; + gSystem->Unload(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libretrklev2_C.so"; + gSystem->Unload(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/libretrklev2.so"; + gSystem->Unload(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/liboptrklev2_C.so"; + gSystem->Unload(bdir.str().c_str()); + bdir.str(""); + bdir << sdir << "/liboptrklev2.so"; + gSystem->Unload(bdir.str().c_str()); + }; + // unloadf77lib(Framework,MySQL); #endif - return(retval); + return(retval); }