LEVEL2/macros/CaloLEVEL2.cra

//
//   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)
//
//   The only input needed is the path to the directory created by YODA for the data file you want to analyze.
//
//   Changelog:
//
//   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.
//
//   3.17 - 4.00 (2005-11-29): preparing for the final release, small changes in the makefile.
//
//   3.16 - 3.17 (2005-10-13): do not exclude all strips marked "bad" but only the ones with RMS>6.
//
//   3.15 - 3.16 (2005-10-04): now it reads Tracker Data version 2.00. Process self-trigger events assuming R=1000 GV and using the calorimeter fit of the track; 
//                             self-trigger events are still marked as "bad" and there is no energy info, nor any angular correction on the track.
//
//   3.14 - 3.15 (2005-09-20): QX22 and NX22 are not related to the last plane but the 11th! fixed!  
//
//   3.13 - 3.14 (2005-09-09): Still some bugs in 64 bit arch, fixed.
//
//   3.12 - 3.13 (2005-09-08): Bug: variable cbar is filled incorrectly. Fixed.
//
//   3.11 - 3.12 (2005-08-31): Added two variables in the ntuple (thex, they sin of x/y angles as measured by the calorimeter)
//
//   3.10 - 3.11 (2005-08-30): Don't run in MySQL standalone mode even if in FORCE mode.
//
//   3.09 - 3.10 (2005-08-05): 64 bit arch bugs fixed.
//
//   3.08 - 3.09 (2005-07-27): x and y view were inverted when fitting the track using calorimeter information, fixed. (affect variables cbar, cibar, npcfit and varcfit).
//
//   3.07 - 3.08 (2005-07-25): planetot and qmean always empty in rootple, fixed. Added variables on goodness of the calorimeter fit (npcfit, varcfit).
//
//   3.06 - 3.07 (2005-07-21): If not in FORCE mode abort processing if no calibration file is found.
//
//   3.05 - 3.06 (2005-07-19): Infinite loop in case of OBT jumps, fixed. NaN in the ntuple and rootple in some cases, fixed (lowered optimization level in 
//                             the compilation of the ntuple and check for NaN from the tracking routine). Added variables in the ntuple/rootple: planetot, qmean.
//
//   3.04 - 3.05 (2005-07-18): No more class CalorimeterLevel2 since it makes very slow the filling of the rootple (unknown reason). Now the tree is filled
//                             with variables from structure CaLevel2.
//
//   3.03 - 3.04 (2005-07-15): More tuning, now it can compile and it works! Not loaded anymore yodaUtility.c, use clone routines instead.
//
//   3.02 - 3.03 (2005-07-13): Start the tuning for compilation and speed optimization.
//
//   3.01 - 3.02 (2005-07-11): Correction for high energy deposit. Correction in the calorimeter planes position.
//
//   3.00 - 3.01 (2005-07-11): Infinite loop in some cases of lost of sync, fixed. Cannot write ntuples bigger than 3 MB, fixed.
//
//   2.10 - 3.00 (2005-07-05): Include track information.
//
//   2.09 - 2.10 (2005-07-05): Started the tuning for the compilation of the macro. Added include files of ROOT subroutines.
//
//   2.08 - 2.09 (2005-06-29): Simplification in the check of goodness for the events. Integer convention for software code.
//
//   2.07 - 2.08 (2005-06-24): Changed the structure of the release and added the enviromental variables.
//
//   2.06 - 2.07 (2005-06-23): Added boolean variables to ntuple: good, crc, perr, swerr.
//
//   2.05 - 2.06 (2005-06-22): Trigger flag set as default to -1. It will remain -1 if the program is not able to determine who gave the trigger.
//                             Added: new tree (Software)/ntuple (number 2) with one entry only which contains the information about the software used to 
//                             produce the LEVEL2 data. Changed names in commons/struct/ntuples/rootples: evfile -> pro_num, headc -> pkt_num .
//                             Disabled MySQL input variable (still there for testing purpouse).
//
//   2.04 - 2.05 (2005-06-17): Trigger flag is not always correct for S4 triggers, fixed with a workaround.
//
//   2.03 - 2.04 (2005-06-16): Bug, when casting from integer to float the resulting float is WRONG! changed in the rootple OBT, HEADC, EVFILE from real to integer.
//                             Changes also in the fortran libraries and passing structures.
//
//   2.02 - 2.03 (2005-06-15): Bug, in PAW and FORCE mode the existence of file is not checked, finally fixed.
//
//   2.01 - 2.02 (2005-06-07): Bug, the calibration used could be incorrect for some events. Fixed: check that the event time is in the calibration window, 
//                             if not search for another calibration.
//
//   2.00 - 2.01 (2005-06-07): Check if file is older than 050515_007.
//
//   1.04 - 2.00 (2005-06-06): Added support to MySQL Trieste database.
//
//   1.03 - 1.04 (2005-05-31): Now it can unload all the libraries and you can run this script as many times you want without quitting ROOT.
//
//   1.02 - 1.03 (2005-05-27): cannot find libraries and calibration file when running twice in the same ROOT section, fixed. KNOWN BUG: cannot unload cfillcalol2_C!!
//
//   1.01 - 1.02 (2005-05-18): filename corruption in some cases. Workaround for the moment...
//
//   1.00 - 1.01 (2005-05-17): problems when starting with a possible backup calibration file. Fixed.
//
//   0.01 - 1.00 (2005-05-06): first time it works. evfile must be i+1 not i.
//
//   0.01        (2005-05-04): created.
//
#include <fstream>
#include <TTree.h>
#include <TClassEdit.h>
#include <TObject.h>
#include <TList.h>
#include <TSystem.h>
#include <TSystemDirectory.h>
#include <TString.h>
#include <TFile.h>
#include <TClass.h>
#include <TCanvas.h>
#include <TH1.h>
#include <TH1F.h>
#include <TH2D.h>
#include <TLatex.h>
#include <TPad.h>
//
#include <calomysqlstruct.h>
#include <ccal2struct.h>
#include <ctrkstruct.h>
#include <ctrkinclude.h>
//
#if !defined (__CINT__)                         
extern  short int creadmysql(char *, const char *, struct Mystruct &);
extern  void cfillcalol2(struct CaLevel1 &, struct CaLevel2 &);
extern  void creadB(const char []);
extern  void coptrklev2(const char [], struct Tracklev2 &, int &);
extern  void copencalol2(char [], Int_t &, Int_t &, int &);
extern  void cretrklev2(int &, Tracklev2 &);
extern  void cclosecalol2(char []);
extern  void ccltrklev2(struct Tracklev2 &);
extern  const char *pathtocalibration();
#include <event/PamelaRun.h>
#include <event/physics/calorimeter/CalorimeterEvent.h>
#include <event/physics/trigger/TriggerEvent.h>
#include <event/CalibCalPedEvent.h>
#endif
//
#include <caloclasses.h>
#include <CaloFunctions.h>
//
#if !defined (__CINT__)                         
void unloadf77lib(TString Framework, Int_t MySQL){}
#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());
    };    
    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 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
    //
    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 {    
        //
        // for each event check that the calibration we are using are still within calibration limits, if not call the next calibration
        //
        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] ; 
                  };
                };
                //
                // 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);
                        };
                    };
                };
                //
                // 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; nn<upnn; nn++ ){
                                        evento.dexyc[l][m][nn] = de->dexyc[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;
                        };
                    };
                }; 
                //
                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);
}
//
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;
    //
#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();
#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.
    //
#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());
    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 << "/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 defined (__CINT__)  
    Int_t MySQLload = 0;
#endif
    //
    // 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);
#endif
        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);
#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");
#if defined (__CINT__)                          
        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");
#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(3);          
            } else {
                printf("\n OK, I will create it!\n\n");
            };
        };
    };
    //    
    // 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 {
        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);
#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());
#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");
#if defined (__CINT__)                          
            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();
    } 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");       
        };
        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());
#if defined (__CINT__)                          
            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;
    };

    //
    // 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";
        };
        //
        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());
#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);
        };
    } 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); 
    } 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
        //
        // 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);
        };
        //
        // 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++;
                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; e<trk.ntrk ; e++){
            if ( trk.image[e] != 0 ) imtrack++;
          };
          if ( (trk.ntrk-(imtrack/2)) == 1 ){ 
            calib.good2 = trk.good2;
            for (Int_t e = 0; e<5 ; e++){
              calib.al_p[e][0] = trk.al[0][e];
              if ( imtrack != 0 ) calib.al_p[e][1] = trk.al[1][e];
            };
            if ( imtrack != 0 ) {
              if ( trk.chi2[0]>trk.chi2[1] ) {
                calib.trkchi2 = 1;
              } else {
                calib.trkchi2 = 2;
              };
            } else {
              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.;
            };              
          };
        } else {
          calib.good2 = false;
          for (Int_t e = 0; e<5 ; e++){
            calib.al_p[e][0] = 0.;
            calib.al_p[e][1] = 0.;
          };        
        };
        //
        // 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;
    };
    //
    // 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);
#endif
    return(retval);
}
