/[PAMELA software]/DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp
ViewVC logotype

Diff of /DarthVader/CalorimeterLevel2/src/CaloLevel0.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.1 by mocchiut, Sun Mar 18 20:22:54 2007 UTC revision 1.25 by mocchiut, Fri Dec 12 16:18:17 2008 UTC
# Line 56  CaloLevel0::CaloLevel0(){ Line 56  CaloLevel0::CaloLevel0(){
56    clevel1 = &clevel1_;    clevel1 = &clevel1_;
57    clevel2 = &clevel2_;    clevel2 = &clevel2_;
58    //    //
59    //  extern struct FlEventi eventi_;
60    //   extern struct FlGruppo gruppo_;
61    //   extern struct FlGruppo2 gruppo2_;
62    //   extern struct FlGruppo4 gruppo4_;
63    //   extern struct FlTaglioen taglioen_;
64    //   extern struct FlAngolo angolo_;
65    //   extern struct FlWhere where_;
66    //   extern struct FlGeneral general_;
67    //   extern struct FlCh ch_;
68    //   extern struct FlCalofit calofit_;
69    //   extern struct FlPawcd pawcd_;
70    //   extern struct FlQuestd questd_;
71    //   eventi = &eventi_;
72    //   gruppo = &gruppo_;
73    //   gruppo2 = &gruppo2_;
74    //   gruppo4 = &gruppo4_;
75    //   taglioen = &taglioen_;
76    //   angolo = &angolo_;
77    //   where = &where_;
78    //   general = &general_;
79    //   ch = &ch_;
80    //   calofit = &calofit_;
81    //   pawcd = &pawcd_;
82    //   questd = &questd_;
83      //
84    trkseqno = 0;    trkseqno = 0;
85    ClearStructs();    ClearStructs();
86    //    //
# Line 64  CaloLevel0::CaloLevel0(){ Line 89  CaloLevel0::CaloLevel0(){
89    memset(mip, 0, 2*22*96*sizeof(Float_t));    memset(mip, 0, 2*22*96*sizeof(Float_t));
90    memset(base, 0, 2*22*6*sizeof(Float_t));    memset(base, 0, 2*22*6*sizeof(Float_t));
91    memset(sbase, 0, 2*22*6*sizeof(Float_t));    memset(sbase, 0, 2*22*6*sizeof(Float_t));
92      memset(obadmask, 0, 2*22*96*sizeof(Int_t));
93      memset(obadpulsemask, 0, 2*22*6*sizeof(Int_t));
94      memset(ctprecor, 0, 2*22*6*sizeof(Float_t));
95      memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
96      memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
97    calopar1 = true;    calopar1 = true;
98    calopar2 = true;    calopar2 = true;
99    calopar3 = true;    calopar3 = true;
100      calopar4 = true;
101      calopar5 = true;
102    crosst = true;    crosst = true;
103    ftcalopar1 = 0;    ftcalopar1 = 0;
104    ttcalopar1 = 0;    ttcalopar1 = 0;
# Line 74  CaloLevel0::CaloLevel0(){ Line 106  CaloLevel0::CaloLevel0(){
106    ttcalopar2 = 0;    ttcalopar2 = 0;
107    ftcalopar3 = 0;    ftcalopar3 = 0;
108    ttcalopar3 = 0;    ttcalopar3 = 0;
109      ftcalopar4 = 0;
110      ttcalopar4 = 0;
111      ftcalopar5 = 0;
112      ttcalopar5 = 0;
113  }  }
114    
115  void CaloLevel0::SetCrossTalk(Bool_t ct){  void CaloLevel0::SetCrossTalk(Bool_t ct){
116    crosst = ct;    crosst = ct;
117  };  }
118    
119    void CaloLevel0::SetCrossTalkType(Bool_t ct){
120      ctground = ct;
121    }
122    
123    void CaloLevel0::SetCrossTalkType(Int_t ct){
124      if ( ct == 0 ) ctground = true;
125      if ( ct == 1 ){
126        ctground = false;
127        noselfct = false;
128      };
129      if ( ct == 2 ){
130        ctground = false;
131        noselfct = true;
132      };
133    }
134    
135  void CaloLevel0::SetVerbose(Bool_t ct){  void CaloLevel0::SetVerbose(Bool_t ct){
136    verbose = ct;    verbose = ct;
137  };  }
138    
139  /**  /**
140   * Initialize CaloLevel0 object   * Initialize CaloLevel0 object
141  **/  **/
142  void CaloLevel0::ProcessingInit(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){  void CaloLevel0::ProcessingInit(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
143      if ( !dbc->IsConnected() ) throw -116;  
144      this->InitDo(dbc,hs,sgnl,l0tree,isdeb,isverb);
145    }
146    
147    /**
148     * Initialize CaloLevel0 object
149    **/
150    void CaloLevel0::ProcessingInit(GL_TABLES *glt, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
151      //
152      const TString host = glt->CGetHost();
153      const TString user = glt->CGetUser();
154      const TString psw = glt->CGetPsw();
155      TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
156      if ( !dbc->IsConnected() ) throw -116;  
157      this->InitDo(dbc,hs,sgnl,l0tree,isdeb,isverb);
158      dbc->Close();
159      delete dbc;
160    }
161    
162    
163    void CaloLevel0::InitDo(TSQLServer *dbc, UInt_t hs, Int_t &sgnl, TTree *l0tree, Bool_t isdeb, Bool_t isverb){
164      stringstream myquery;
165      myquery.str("");
166      myquery << "SET time_zone='+0:00'";
167      dbc->Query(myquery.str().c_str());
168    //    //
169    debug = isdeb;    debug = isdeb;
170    verbose = isverb;    verbose = isverb;
# Line 169  void CaloLevel0::ProcessingInit(TSQLServ Line 246  void CaloLevel0::ProcessingInit(TSQLServ
246    //    //
247  }  }
248    
249  Int_t CaloLevel0::ChkCalib(TSQLServer *dbc, UInt_t atime){  Int_t CaloLevel0::ChkCalib(GL_TABLES *glt, UInt_t atime){
250    Int_t sgnl = 0;    Int_t sgnl = 0;
251    for ( Int_t s = 0; s < 4; s++){    for ( Int_t s = 0; s < 4; s++){
252      if ( atime > totime[s] ){      if ( atime > totime[s] ){
253        if ( !dbc->IsConnected() ) throw -116;        sgnl = Update(glt,atime,s);
       sgnl = Update(dbc,atime,s);  
254        if ( sgnl < 0 ) return(sgnl);        if ( sgnl < 0 ) return(sgnl);
255      };      };
256    };    };
257    return(sgnl);    return(sgnl);
258  }  }
259    
260  Int_t CaloLevel0::ChkParam(TSQLServer *dbc, UInt_t runheader){  Int_t CaloLevel0::ChkParam(TSQLServer *dbc, UInt_t runheader, Bool_t mechal){
261      Int_t sig = this->ChkParamDo(dbc,runheader,mechal);
262      return(sig);
263    }
264    
265    Int_t CaloLevel0::ChkParam(GL_TABLES *glt, UInt_t runheader, Bool_t mechal){
266      const TString host = glt->CGetHost();
267      const TString user = glt->CGetUser();
268      const TString psw = glt->CGetPsw();
269      TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
270      if ( !dbc->IsConnected() ) throw -116;
271      stringstream myquery;
272      myquery.str("");
273      myquery << "SET time_zone='+0:00'";
274      dbc->Query(myquery.str().c_str());
275      //
276      Int_t sig = this->ChkParamDo(dbc,runheader,mechal);
277      dbc->Close();
278      delete dbc;
279      return(sig);
280    }
281    
282    Int_t CaloLevel0::ChkParamDo(TSQLServer *dbc, UInt_t runheader, Bool_t mechal){
283      //
284    stringstream calfile;    stringstream calfile;
285    stringstream bmfile;    stringstream bmfile;
286    stringstream aligfile;    stringstream aligfile;
# Line 191  Int_t CaloLevel0::ChkParam(TSQLServer *d Line 290  Int_t CaloLevel0::ChkParam(TSQLServer *d
290    GL_PARAM *glparam = new GL_PARAM();    GL_PARAM *glparam = new GL_PARAM();
291    //    //
292    if ( calopar1 || ( ttcalopar1 != 0 && ttcalopar1 < runheader ) ){    if ( calopar1 || ( ttcalopar1 != 0 && ttcalopar1 < runheader ) ){
293        //
294        if ( debug ) printf(" calopar1 %i ftcalopar1 %u ttcalopar1 %u runheader %u \n",calopar1,ftcalopar1,ttcalopar1,runheader);
295        //
296        if ( calopar1 ){
297          //
298          // determine where I can find calorimeter ADC to MIP conversion file  
299          //
300          if ( verbose ) printf(" Querying DB for calorimeter parameters files...\n");
301          //
302          error = 0;
303          error = glparam->Query_GL_PARAM(runheader,101,dbc);
304          if ( error < 0 ) return(error);
305          //
306          calfile.str("");
307          calfile << glparam->PATH.Data() << "/";
308          calfile << glparam->NAME.Data();
309          //
310          if ( verbose ) printf("\n Using ADC to MIP conversion file: \n %s \n",calfile.str().c_str());
311          f = fopen(calfile.str().c_str(),"rb");
312          if ( !f ){
313            if ( verbose ) printf(" CALORIMETER - ERROR: no ADC to MIP file!\n");
314            return(-105);
315          };
316          //
317          for (Int_t m = 0; m < 2 ; m++ ){
318            for (Int_t k = 0; k < 22; k++ ){
319              for (Int_t l = 0; l < 96; l++ ){
320                fread(&mip[m][k][l],sizeof(mip[m][k][l]),1,f);
321                if ( debug ) printf(" %f \n",mip[m][k][l]);
322              };
323            };
324          };
325          fclose(f);      
326        };
327        //
328      calopar1 = false;      calopar1 = false;
329      //      //
330      // determine where I can find calorimeter ADC to MIP conversion file        // flight extra corrections:
331      //      //
332      if ( verbose ) printf(" Querying DB for calorimeter parameters files...\n");      if ( verbose ) printf(" Querying DB for calorimeter flight ADC to MIP files...\n");
333      //      //
334      error = 0;      error = 0;
335      error = glparam->Query_GL_PARAM(runheader,101,dbc);      error = glparam->Query_GL_PARAM(runheader,110,dbc);
336      if ( error < 0 ) return(error);      if ( error < 0 ) return(error);
337      //      //
338      calfile.str("");      calfile.str("");
# Line 207  Int_t CaloLevel0::ChkParam(TSQLServer *d Line 341  Int_t CaloLevel0::ChkParam(TSQLServer *d
341      ftcalopar1 = glparam->FROM_TIME;      ftcalopar1 = glparam->FROM_TIME;
342      ttcalopar1 = glparam->TO_TIME;      ttcalopar1 = glparam->TO_TIME;
343      //      //
344      if ( verbose ) printf("\n Using ADC to MIP conversion file: \n %s \n",calfile.str().c_str());      if ( verbose ) printf("\n Using ADC to MIP special conversion file: \n %s \n",calfile.str().c_str());
345      f = fopen(calfile.str().c_str(),"rb");      ifstream spfile;
346      if ( !f ){      spfile.open(calfile.str().c_str());
347        if ( verbose ) printf(" CALORIMETER - ERROR: no ADC to MIP file!\n");      if ( !spfile ){
348        return(-105);        if ( verbose ) printf(" CALORIMETER - ERROR: no special calibration file!\n");
349          return(-123);
350      };      };
351        //  
352        Int_t vview = 0;
353        Int_t vplane = 0;
354        Int_t vstrip = 0;
355        Float_t vval = 0.;
356        while ( spfile >> vview && spfile >> vplane && spfile >> vstrip && spfile >> vval){
357          if ( debug ) printf(" Setting ADC to MIP conversion factor: view %i plane %i strip %i mip %f  \n",vview,vplane,vstrip,vval);
358          mip[vview][vplane][vstrip] = vval;
359        };
360      //      //
     for (Int_t m = 0; m < 2 ; m++ ){  
       for (Int_t k = 0; k < 22; k++ ){  
         for (Int_t l = 0; l < 96; l++ ){  
           fread(&mip[m][k][l],sizeof(mip[m][k][l]),1,f);  
         };  
       };  
     };  
     fclose(f);  
361    };    };
362    //    //
363      //
364    if ( calopar2 || ( ttcalopar2 != 0 && ttcalopar2 < runheader ) ){    if ( calopar2 || ( ttcalopar2 != 0 && ttcalopar2 < runheader ) ){
365        //
366        if ( debug ) printf(" calopar2 %i ftcalopar2 %u ttcalopar2 %u runheader %u \n",calopar2,ftcalopar2,ttcalopar2,runheader);
367      calopar2 = false;      calopar2 = false;
368      //      //
369      // determine where I can find calorimeter alignment file        // determine where I can find calorimeter alignment file  
# Line 240  Int_t CaloLevel0::ChkParam(TSQLServer *d Line 379  Int_t CaloLevel0::ChkParam(TSQLServer *d
379      ftcalopar2 = glparam->FROM_TIME;      ftcalopar2 = glparam->FROM_TIME;
380      ttcalopar2 = glparam->TO_TIME;      ttcalopar2 = glparam->TO_TIME;
381      //      //
382      if ( verbose ) printf("\n Using alignment file: \n %s \n",aligfile.str().c_str());      if ( verbose ) printf("\n Using parameter file: \n %s \n",aligfile.str().c_str());
383      f = fopen(aligfile.str().c_str(),"rb");      f = fopen(aligfile.str().c_str(),"rb");
384      if ( !f ){      if ( !f ){
385        if ( verbose ) printf(" CALORIMETER - ERROR: no alignement file!\n");        if ( verbose ) printf(" CALORIMETER - ERROR: no parameter file!\n");
386        return(-106);        return(-106);
387      };      };
388      //      //
389      fread(&clevel1->xalig,sizeof(clevel1->xalig),1,f);      if ( !mechal ){
390      if ( debug ) printf(" xalig = %f \n",clevel1->xalig);        //
391      fread(&clevel1->yalig,sizeof(clevel1->yalig),1,f);        fread(&clevel1->xalig,sizeof(clevel1->xalig),1,f);
392      if ( debug ) printf(" yalig = %f \n",clevel1->yalig);        if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
393      fread(&clevel1->zalig,sizeof(clevel1->zalig),1,f);        fread(&clevel1->yalig,sizeof(clevel1->yalig),1,f);
394      if ( debug ) printf(" zalig = %f \n",clevel1->zalig);        if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
395          fread(&clevel1->zalig,sizeof(clevel1->zalig),1,f);
396          if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
397        } else {
398          if ( verbose ) printf("\n Using MECHANICAL alignement parameters \n");
399          //
400          CaloStrip cs = CaloStrip();
401          cs.UseMechanicalAlig();
402          clevel1->xalig = cs.GetXalig();
403          if ( debug ) printf(" xalig = %f \n",clevel1->xalig);
404          clevel1->yalig = cs.GetYalig();
405          if ( debug ) printf(" yalig = %f \n",clevel1->yalig);
406          clevel1->zalig = cs.GetZalig();
407          if ( debug ) printf(" zalig = %f \n",clevel1->zalig);
408          //
409          Float_t tmp = 0;
410          fread(&tmp,sizeof(clevel1->xalig),1,f);
411          fread(&tmp,sizeof(clevel1->yalig),1,f);
412          fread(&tmp,sizeof(clevel1->zalig),1,f);
413          //      clevel1->zalig = -265.82;
414          //    
415        };
416      fread(&clevel1->emin,sizeof(clevel1->emin),1,f);        fread(&clevel1->emin,sizeof(clevel1->emin),1,f);  
417      if ( debug ) printf(" signal threshold = %f \n",clevel1->emin);      if ( debug ) printf(" signal threshold = %f \n",clevel1->emin);
418      //      //
# Line 262  Int_t CaloLevel0::ChkParam(TSQLServer *d Line 422  Int_t CaloLevel0::ChkParam(TSQLServer *d
422    // Load offline bad strip mask    // Load offline bad strip mask
423    //    //
424    if ( calopar3 || ( ttcalopar3 != 0 && ttcalopar3 < runheader ) ){    if ( calopar3 || ( ttcalopar3 != 0 && ttcalopar3 < runheader ) ){
425        if ( debug ) printf(" calopar3 %i ftcalopar3 %u ttcalopar3 %u runheader %u \n",calopar3,ftcalopar3,ttcalopar3,runheader);
426      calopar3 = false;      calopar3 = false;
427      //      //
428      // determine where I can find calorimeter alignment file        // determine where I can find calorimeter alignment file  
# Line 308  Int_t CaloLevel0::ChkParam(TSQLServer *d Line 469  Int_t CaloLevel0::ChkParam(TSQLServer *d
469      badfile.close();      badfile.close();
470    };    };
471    //    //
472      // calopar4
473      //
474      if ( calopar4 || ( ttcalopar4 != 0 && ttcalopar4 < runheader ) ){
475        //
476        if ( debug ) printf(" calopar4 %i ftcalopar4 %u ttcalopar4 %u runheader %u \n",calopar4,ftcalopar4,ttcalopar4,runheader);
477        //
478        calopar4 = false;
479        //
480        // flight extra corrections:
481        //
482        if ( verbose ) printf(" Querying DB for calorimeter max rms file...\n");
483        //
484        error = 0;
485        error = glparam->Query_GL_PARAM(runheader,109,dbc);
486        if ( error < 0 ) return(error);
487        //
488        calfile.str("");
489        calfile << glparam->PATH.Data() << "/";
490        calfile << glparam->NAME.Data();
491        ftcalopar4 = glparam->FROM_TIME;
492        ttcalopar4 = glparam->TO_TIME;
493        //
494        if ( verbose ) printf("\n Using calorimeter max rms file: \n %s \n",calfile.str().c_str());
495        ifstream spfile;
496        spfile.open(calfile.str().c_str());
497        if ( !spfile ){
498          if ( verbose ) printf(" CALORIMETER - ERROR: no max rms file!\n");
499          return(-124);
500        };
501        //  
502        Int_t vview = 0;
503        Int_t vplane = 0;
504        Int_t vval = 0;
505        for (Int_t l=0; l<2; l++){
506          for (Int_t m=0; m<22; m++){
507            maxrms[l][m] = 26;
508          };
509        };
510        while ( spfile >> vview && spfile >> vplane && spfile >> vval){
511          if ( debug ) printf(" Setting view %i plane %i max rms %i  \n",vview,vplane,vval);
512          maxrms[vview][vplane] = vval;
513        };
514        spfile.close();
515        //
516      };
517      //
518      // calopar5
519      //
520      if ( calopar5 || ( ttcalopar5 != 0 && ttcalopar5 < runheader ) ){
521        //
522        if ( debug ) printf(" calopar5 %i ftcalopar5 %u ttcalopar5 %u runheader %u \n",calopar5,ftcalopar5,ttcalopar5,runheader);
523        //
524        calopar5 = false;
525        //
526        // flight extra corrections:
527        //
528        if ( verbose ) printf(" Querying DB for calorimeter noise to signal threshold file...\n");
529        //
530        error = 0;
531        error = glparam->Query_GL_PARAM(runheader,111,dbc);
532        if ( error < 0 ) return(error);
533        //
534        calfile.str("");
535        calfile << glparam->PATH.Data() << "/";
536        calfile << glparam->NAME.Data();
537        ftcalopar5 = glparam->FROM_TIME;
538        ttcalopar5 = glparam->TO_TIME;
539        //
540        if ( verbose ) printf("\n Using calorimeter noise to signal threshold file: \n %s \n",calfile.str().c_str());
541        ifstream spfile;
542        spfile.open(calfile.str().c_str());
543        if ( !spfile ){
544          if ( verbose ) printf(" CALORIMETER - ERROR: no noise to signal threshold file!\n");
545          return(-125);
546        };
547        //  
548        Int_t vview = 0;
549        Int_t vplane = 0;
550        Int_t vstrip = 0;
551        Float_t vval = 0.;
552        for (Int_t l=0; l<2; l++){
553          for (Int_t m=0; m<22; m++){
554            for (Int_t n=0; n<96; n++){
555              memin[l][m][n] = 0.7;
556            };
557          };
558        };
559        while ( spfile >> vview && spfile >> vplane && spfile >> vstrip && spfile >> vval){
560          if ( vstrip == -1 ){
561            for (Int_t ll=0; ll<96; ll++){
562              if ( debug ) printf(" Setting view %i plane %i strip %i noise to signal ratio %f  \n",vview,vplane,ll,vval);
563              memin[vview][vplane][ll] = vval;
564            };
565          } else {
566            if ( debug ) printf(" Setting view %i plane %i strip %i noise to signal ratio %f  \n",vview,vplane,vstrip,vval);
567            memin[vview][vplane][vstrip] = vval;
568          };
569        };
570        spfile.close();
571        //
572      };
573      //
574      //
575    delete glparam;    delete glparam;
576    //    //
577    return(0);    return(0);
578  }  }
579    
580    Int_t CaloLevel0::CalcCrossTalkCorr(TSQLServer *dbc, UInt_t runheader, Bool_t ctusetable){
581      Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,ctusetable);
582      return(sig);  
583    };
584    
585    Int_t CaloLevel0::CalcCrossTalkCorr(TSQLServer *dbc, UInt_t runheader){
586      Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,true);
587      return(sig);
588    }
589    
590  void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){  Int_t CaloLevel0::CalcCrossTalkCorr(GL_TABLES *glt, UInt_t runheader, Bool_t usetable){
591      Float_t minstrip = 100000.;    const TString host = glt->CGetHost();
592      Float_t rms = 0.;    const TString user = glt->CGetUser();
593      base[l][m][pre] = 0.;    const TString psw = glt->CGetPsw();
594      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){    TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
595          if ( calgood[l][m][e] == 0. && obadmask[l][m][e] == 0 &&  dexy[l][m][e]-calped[l][m][e] < minstrip &&  dexy[l][m][e] > 0.) {    if ( !dbc->IsConnected() ) throw -116;
596              minstrip = dexy[l][m][e]-calped[l][m][e];    stringstream myquery;
597              rms = calthr[l][m][pre];    myquery.str("");
598          };    myquery << "SET time_zone='+0:00'";
599      dbc->Query(myquery.str().c_str());
600      //
601      Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,usetable);
602      dbc->Close();
603      delete dbc;
604      //
605      return(sig);
606      //  
607    };
608    
609    Int_t CaloLevel0::CalcCrossTalkCorr(GL_TABLES *glt, UInt_t runheader){
610      const TString host = glt->CGetHost();
611      const TString user = glt->CGetUser();
612      const TString psw = glt->CGetPsw();
613      TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
614      if ( !dbc->IsConnected() ) throw -116;
615      stringstream myquery;
616      myquery.str("");
617      myquery << "SET time_zone='+0:00'";
618      dbc->Query(myquery.str().c_str());
619      //
620      Int_t sig = CalcCrossTalkCorrDo(dbc,runheader,true);
621      dbc->Close();
622      delete dbc;
623      //
624      return(sig);
625      //
626    }
627    
628    Int_t CaloLevel0::CalcCrossTalkCorrDo(TSQLServer *dbc, UInt_t runheader, Bool_t usetable){
629      //
630      if ( ctground ) return(0);
631      //
632      Int_t error = 0;
633      GL_PARAM *glparam = new GL_PARAM();  
634      //
635      // determine where I can find file with offline bad pulser mask
636      //
637      stringstream bmfile;
638      error = 0;
639      error = glparam->Query_GL_PARAM(runheader,105,dbc);
640      if ( error < 0 ) return(error);
641      //
642      bmfile.str("");
643      bmfile << glparam->PATH.Data() << "/";
644      bmfile << glparam->NAME.Data();
645      //
646      ifstream badfile;
647      if ( verbose ) printf("\n Using bad pulser offline mask file: \n %s \n\n",bmfile.str().c_str());
648      badfile.open(bmfile.str().c_str());
649      if ( !badfile ){
650        if ( verbose ) printf(" CALORIMETER - ERROR: no bad pulser offline mask file!\n");
651        return(-115);
652      };
653      //
654      Bool_t isdone = false;
655      Int_t bad = 0;
656      Int_t view = 1;
657      Int_t pre = 0;
658      Int_t plane = 21;
659      while ( !isdone ) {
660        badfile >> bad;
661        obadpulsemask[view][plane][pre] = bad;
662        if ( debug && bad ) printf(" SETTING view %i plane %i pre %i  BAD = %i \n",view,plane,pre,bad);
663        pre++;
664        if ( pre > 5 ){
665          pre = 0;
666          plane--;
667          if ( plane < 0 ){
668            plane = 21;
669            view--;
670          };
671          if ( view < 0 ) isdone = true;
672      };      };
673      if ( minstrip != 100000. ) {    };
674          Float_t strip6s = 0.;    //
675          for (Int_t e = pre*16; e < (pre+1)*16 ; e++){    badfile.close();
676              if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {    if ( !usetable ){
677                  strip6s += 1.;      //
678                  base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);      // Let's start with cross-talk correction calculation
679        //
680        GL_CALOPULSE_CALIB *glp = new GL_CALOPULSE_CALIB();
681        Float_t adcp[2][22][96];
682        Float_t adcpcal[2][22][96];
683        memset(adcp , 0, 2*22*96*sizeof(Float_t));
684        memset(adcpcal , 0, 2*22*96*sizeof(Float_t));
685        //
686        UInt_t pampli = 0;  
687        for (Int_t s=0; s<4; s++){
688          //
689          // Save into matrix adcp the values of the highest pulse calibration (pulse amplitude = 2)
690          //
691          pampli = 2;
692          error = 0;
693          error = glp->Query_GL_CALOPULSE_CALIB(runheader,s,pampli,dbc);
694          if ( error < 0 ){
695            if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
696            return(error);
697          };
698          //
699          UInt_t idcalib = glp->ID_ROOT_L0;
700          UInt_t fromtime = glp->FROM_TIME;
701          UInt_t calibno = glp->EV_ROOT;
702          //
703          // determine path and name and entry of the calibration file
704          //
705          GL_ROOT *glroot = new GL_ROOT();  
706          if ( verbose ) printf("\n");
707          if ( verbose ) printf(" ** SECTION %i **\n",s);
708          //
709          error = 0;
710          error = glroot->Query_GL_ROOT(idcalib,dbc);
711          if ( error < 0 ){
712            if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
713            return(error);
714          };
715          //        
716          stringstream name;
717          name.str("");
718          name << glroot->PATH.Data() << "/";
719          name << glroot->NAME.Data();
720          //
721          TString fcalname = (TString)name.str().c_str();
722          ifstream myfile;
723          myfile.open(fcalname.Data());
724          if ( !myfile ){    
725            return(-107);
726          };
727          myfile.close();
728          //
729          TFile *File = new TFile(fcalname.Data());
730          if ( !File ) return(-108);
731          TTree *tr = (TTree*)File->Get("CalibCalPulse2");
732          if ( !tr ) return(-119);
733          //
734          TBranch *calo = tr->GetBranch("CalibCalPulse2");
735          //
736          pamela::CalibCalPulse2Event *ce = 0;
737          tr->SetBranchAddress("CalibCalPulse2", &ce);
738          //
739          Long64_t ncalibs = calo->GetEntries();
740          //
741          if ( !ncalibs ) return(-110);
742          //
743          calo->GetEntry(calibno);
744          if ( verbose ) printf(" PULSE2 using entry %u from file %s",calibno,fcalname.Data());
745          //
746          // retrieve calibration table
747          //
748          if ( ce->pstwerr[s] && ce->pperror[s] == 0 && ce->unpackError == 0 ){
749            for ( Int_t d=0 ; d<11 ;d++  ){
750              for ( Int_t j=0; j<96 ;j++){
751                if ( s == 2 ){
752                  adcp[0][2*d+1][j] = ce->calpuls[3][d][j];
753              };              };
754              //              if ( s == 3 ){
755              //  compression                adcp[0][2*d][j] = ce->calpuls[1][d][j];
756              //              };
757              if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {              if ( s == 0 ){
758                  dexyc[l][m][e] = 0.;                adcp[1][2*d][j] = ce->calpuls[0][d][j];
759              } else {              };
760                  dexyc[l][m][e] = dexy[l][m][e];              if ( s == 1 ){
761                  adcp[1][2*d+1][j] = ce->calpuls[2][d][j];
762                };
763              };
764            };
765          } else {
766            if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
767            return(-111);
768          };
769          //
770          File->Close();
771          delete glroot;
772          //
773          // Save into matrix adcpcal the calibrated values of the pulse calibration (subtraction of pulse amplitude = 0 relative to the pulse2 calibration used)
774          //
775          pampli = 0;
776          error = 0;
777          error = glp->Query_GL_CALOPULSE_CALIB(fromtime,s,pampli,dbc);
778          if ( error < 0 ){
779            if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
780            return(error);
781          };
782          //
783          idcalib = glp->ID_ROOT_L0;
784          calibno = glp->EV_ROOT;
785          //
786          // determine path and name and entry of the calibration file
787          //
788          glroot = new GL_ROOT();  
789          if ( verbose ) printf("\n");
790          if ( verbose ) printf(" ** SECTION %i **\n",s);
791          //
792          error = 0;
793          error = glroot->Query_GL_ROOT(idcalib,dbc);
794          if ( error < 0 ){
795            if ( verbose ) printf(" CALORIMETER - ERROR: error from GLTables\n");
796            return(error);
797          };
798          //        
799          name.str("");
800          name << glroot->PATH.Data() << "/";
801          name << glroot->NAME.Data();
802          //
803          fcalname = (TString)name.str().c_str();
804          myfile.open(fcalname.Data());
805          if ( !myfile ){    
806            return(-107);
807          };
808          myfile.close();
809          //
810          TFile *File1 = new TFile(fcalname.Data());
811          if ( !File1 ) return(-108);
812          TTree *tr1 = (TTree*)File1->Get("CalibCalPulse1");
813          if ( !tr1 ) return(-120);
814          //
815          TBranch *calo1 = tr1->GetBranch("CalibCalPulse1");
816          //
817          pamela::CalibCalPulse1Event *ce1 = 0;
818          tr1->SetBranchAddress("CalibCalPulse1", &ce1);
819          //
820          ncalibs = calo1->GetEntries();
821          //
822          if ( !ncalibs ) return(-110);
823          //
824          calo1->GetEntry(calibno);
825          if ( verbose ) printf(" PULSE1 using entry %u from file %s",calibno,fcalname.Data());
826          //
827          // retrieve calibration table
828          //
829          if ( ce1->pstwerr[s] && ce1->pperror[s] == 0 && ce1->unpackError == 0 ){
830            for ( Int_t d=0 ; d<11 ;d++  ){
831              for ( Int_t j=0; j<96 ;j++){
832                if ( s == 2 ){
833                  adcpcal[0][2*d+1][j] = adcp[0][2*d+1][j] - ce1->calpuls[3][d][j];
834              };              };
835                if ( s == 3 ){
836                  adcpcal[0][2*d][j] = adcp[0][2*d][j] - ce1->calpuls[1][d][j];
837                };
838                if ( s == 0 ){
839                  adcpcal[1][2*d][j] = adcp[1][2*d][j] - ce1->calpuls[0][d][j];
840                };
841                if ( s == 1 ){
842                  adcpcal[1][2*d+1][j] = adcp[1][2*d+1][j] - ce1->calpuls[2][d][j];
843                };
844              };
845          };          };
846          if ( strip6s >= 9. ){              } else {
847              Double_t arro = base[l][m][pre]/strip6s;          if ( verbose ) printf(" CALORIMETER - ERROR: problems finding a good calibration in this file! \n\n ");
848              Float_t deci = 1000.*((float)arro - float(int(arro)));                              return(-111);
849              if ( deci < 500. ) {        };
850                  arro = double(int(arro));        //
851          File1->Close();
852          //
853          delete glroot;
854          //
855        };// loop on the four sections
856        //
857        //
858        delete glp;
859        //
860        // Ok, now we can try to calculate the cross-talk correction for each pre-amplifier
861        //
862        for ( Int_t v=0; v<2; v++){
863          if ( debug ) printf(" \n\n NEW VIEW \n");
864          for ( Int_t p=0; p<22; p++){
865            for ( Int_t npre=0; npre<6; npre++){
866              ctprecor[v][p][npre] = 1000.;
867              ctneigcor[v][p][npre] = 1000.;
868              Int_t str0=npre*16;
869              Int_t str16= -1 + (1+npre)*16;
870              //
871              UInt_t neigc = 0;
872              UInt_t farc = 0;
873              UInt_t pulsc = 0;
874              Float_t sigpulsed = 0.;
875              Float_t neigbase = 0.;
876              Float_t farbase = 0.;
877              //
878              // Loop over the strip of the pre and sum all signal far away from pulsed strip, signal in the neighbour(s) strip(s) and save the pulsed signal
879              // moreover count the number of strips in each case
880              //
881              for (Int_t s=str0; s<=str16; s++){
882                if ( adcpcal[v][p][s] > 10000.){
883                  sigpulsed = adcpcal[v][p][s];
884                  pulsc++;
885                  if ( s > str0 ){
886                    neigbase += adcpcal[v][p][s-1];
887                    neigc++;
888                    farbase -= adcpcal[v][p][s-1];
889                    farc--;
890                  };
891                  if ( s < str16 ){
892                    neigbase += adcpcal[v][p][s+1];
893                    neigc++;
894                    farbase -= adcpcal[v][p][s+1];
895                    farc--;
896                  };
897              } else {              } else {
898                  arro = 1. + double(int(arro));                farc++;
899                  farbase += adcpcal[v][p][s];
900              };              };
901              base[l][m][pre] = arro;            };
902          } else {            //
903              base[l][m][pre] = 31000.;            // Now calculate the corrections
904              for (Int_t e = pre*16; e < (pre+1)*16 ; e++){            //
905                  dexyc[l][m][e] = dexy[l][m][e];            Float_t avefarbase = 0.;
906              if ( farc ) avefarbase = farbase/(Float_t)farc;
907              Float_t aveneigbase = 0.;
908              if ( neigc ) aveneigbase = neigbase/(Float_t)neigc;
909              //
910              if ( pulsc == 1 && farc && neigc ){
911                ctprecor[v][p][npre] = -avefarbase/(sigpulsed+fabs(avefarbase));
912                ctneigcor[v][p][npre] = fabs(aveneigbase-avefarbase)/(sigpulsed+fabs(avefarbase));    
913                if ( debug ) printf(" Cross-talk correction View %i Plane %i Pre %i : pre-correction: %f neighbour strips correction %f \n",v,p,npre,ctprecor[v][p][npre],ctneigcor[v][p][npre]);
914              } else {
915                //
916                // did not find the pulsed strip or more than one pulsed strip found!
917                //
918                if ( debug ) printf(" Problems finding the cross-talk corrections: \n View %i Plane %i Pre %i number of pulsed strip %i \n Average faraway baseline %f number of strips %i Average neighbour baseline %f number of neighbour strips %i \n",v,p,npre,pulsc,avefarbase,farc,aveneigbase,neigc);
919                //
920              };
921            };
922            if ( debug ) printf(" \n ==================== \n");
923          };
924        };
925      } else {
926        //
927        // use pre-amply table
928        //    
929        //
930        // determine where I can find file with offline neighbour correction table
931        //
932        stringstream bmfile2;
933        error = 0;
934        error = glparam->Query_GL_PARAM(runheader,106,dbc);
935        if ( error < 0 ) return(error);
936        //
937        bmfile2.str("");
938        bmfile2 << glparam->PATH.Data() << "/";
939        bmfile2 << glparam->NAME.Data();
940        //
941        ifstream badfile2;
942        if ( verbose ) printf("\n Using pre-amply neighbour crosstalk table file: \n %s \n\n",bmfile2.str().c_str());
943        badfile2.open(bmfile2.str().c_str());
944        if ( !badfile2 ){
945          if ( verbose ) printf(" CALORIMETER - ERROR: no pre-amply neighbour crosstalk table file!\n");
946          return(-121);
947        };
948        //  
949        Int_t vview = 0;
950        Int_t vplane = 0;
951        Int_t vpre = 0;
952        Float_t vcorr = 0.;
953        while ( badfile2 >> vview && badfile2 >> vplane && badfile2 >> vpre && badfile2 >> vcorr){
954          if ( debug ) printf(" Pre-amply neighbour correction: view %i plane %i pre %i correction %f  \n",vview,vplane,vpre,vcorr);
955          ctneigcor[vview][vplane][vpre] = vcorr;
956        };
957        //
958        // determine where I can find file with offline SECOND neighbour correction table
959        //
960        stringstream bmfile3;
961        error = 0;
962        error = glparam->Query_GL_PARAM(runheader,107,dbc);
963        if ( error < 0 ) return(error);
964        //
965        bmfile3.str("");
966        bmfile3 << glparam->PATH.Data() << "/";
967        bmfile3 << glparam->NAME.Data();
968        //
969        ifstream badfile3;
970        if ( verbose ) printf("\n Using pre-amply second neighbour crosstalk table file: \n %s \n\n",bmfile3.str().c_str());
971        badfile3.open(bmfile3.str().c_str());
972        if ( !badfile3 ){
973          if ( verbose ) printf(" CALORIMETER - ERROR: no pre-amply second neighbour crosstalk table file!\n");
974          return(-122);
975        };
976        //  
977        Int_t pview = 0;
978        Int_t pplane = 0;
979        Int_t ppre = 0;
980        Float_t pcorr = 0.;
981        while ( badfile3 >> pview && badfile3 >> pplane && badfile3 >> ppre && badfile3 >> pcorr){
982          if ( debug ) printf(" Pre-amply second neighbour correction: view %i plane %i pre %i correction %f  \n",pview,pplane,ppre,-pcorr);
983          ctprecor[pview][pplane][ppre] = -pcorr; // data are saved as negatives in the file
984        };
985        //
986        // determine where to find the file containing the Silicon crosstalk correction table
987        //
988        stringstream bmfile4;
989        error = 0;
990        error = glparam->Query_GL_PARAM(runheader,108,dbc);
991        if ( error < 0 ) return(error);
992        //
993        bmfile4.str("");
994        bmfile4 << glparam->PATH.Data() << "/";
995        bmfile4 << glparam->NAME.Data();
996        //
997        ifstream badfile4;
998        if ( verbose ) printf("\n Using Silicon crosstalk table file: \n %s \n\n",bmfile4.str().c_str());
999        badfile4.open(bmfile4.str().c_str());
1000        if ( !badfile4 ){
1001          if ( verbose ) printf(" CALORIMETER - ERROR: no Silicon crosstalk table file!\n");
1002          return(-125);
1003        };
1004        //  
1005        Int_t spview = 0;
1006        Int_t spplane = 0;
1007        Int_t psil = 0;
1008        Float_t spcorr = 0.;
1009        memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
1010        while ( badfile4 >> spview && badfile4 >> spplane && badfile4 >> psil && badfile4 >> spcorr){
1011          if ( debug ) printf(" Silicon correction: view %i plane %i silicon %i correction %f  \n",spview,spplane,psil,-spcorr);
1012          ctsicor[spview][spplane][psil] = -spcorr; // data are saved as negatives in the file
1013        };
1014        //
1015      };
1016      //
1017      delete glparam;
1018      //
1019      // Check the calculated corrections
1020      //
1021      Int_t opre=0;
1022      Int_t ppre=0;
1023      Bool_t found = false;
1024      for ( Int_t v=0; v<2; v++){
1025        for ( Int_t p=0; p<22; p++){
1026          for ( Int_t npre=0; npre<6; npre++){
1027            if ( ctprecor[v][p][npre] == 1000. || ctneigcor[v][p][npre] == 1000. || obadpulsemask[v][p][npre] != 0 ){
1028              if ( debug ) printf(" Cross-talk correction CHANGED for view %i Plane %i Pre %i\n    BEFORE: pre-correction: %f neighbour strips correction %f \n",v,p,npre,ctprecor[v][p][npre],ctneigcor[v][p][npre]);
1029              if ( npre%2 ){
1030                opre = npre-1;
1031              } else {
1032                opre = npre+1;          
1033              };
1034              if ( ctprecor[v][p][opre] == 1000. || ctneigcor[v][p][opre] == 1000. || obadpulsemask[v][p][opre] != 0 ){
1035                ppre=0;
1036                found = false;
1037                while ( ppre < 6 ){
1038                  if ( ctprecor[v][p][ppre] != 1000. && ctneigcor[v][p][ppre] != 1000. && !obadpulsemask[v][p][ppre] ){
1039                    found = true;
1040                    ctprecor[v][p][npre] = ctprecor[v][p][ppre];
1041                    ctneigcor[v][p][npre] = ctneigcor[v][p][ppre];
1042                    break;
1043                  };
1044                  ppre++;
1045                };
1046                if ( !found ){
1047                  if ( verbose ) printf(" WARNING: cannot find a good cross-talk correction for view %i plane %i pre %i \n Setting to default values 0.002 0.002\n",v,p,npre);
1048                  ctprecor[v][p][npre] = 0.002;
1049                  ctneigcor[v][p][npre] = 0.002;
1050              };              };
1051              } else {
1052                ctprecor[v][p][npre] = ctprecor[v][p][opre];
1053                ctneigcor[v][p][npre] = ctneigcor[v][p][opre];
1054              };
1055              if ( debug ) printf("    AFTER: pre-correction: %f neighbour strips correction %f \n",ctprecor[v][p][npre],ctneigcor[v][p][npre]);
1056          };          };
1057      } else {        };
1058        };
1059      };  
1060      //
1061      return(0);
1062    }
1063    
1064    void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre){
1065      Int_t n = 0;
1066      Float_t q = 0;
1067      this->FindBaseCompress(l,m,pre,n,q);
1068    }
1069    
1070    void CaloLevel0::FindBaseCompress(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
1071      for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1072        dexy[l][m][e] = dexyc[l][m][e];
1073      };  
1074      this->FindBaseRaw(l,m,pre,nst,qp);
1075    }
1076    
1077    void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre){
1078      Int_t n = 0;
1079      Float_t q = 0;
1080      this->FindBaseRaw(l,m,pre,n,q);
1081    }
1082    
1083    void CaloLevel0::FindBaseRaw(Int_t l, Int_t m, Int_t pre, Int_t &nst, Float_t &qp){
1084      //
1085      Float_t minstrip = 100000.;
1086      Float_t rms = 0.;
1087      Int_t process = 0;
1088      Int_t onlmask[16];
1089      memset(onlmask, 0, 16*sizeof(Int_t));
1090      //
1091      while ( process < 2 ){
1092        //
1093        minstrip = 100000.;
1094        rms = 0.;
1095        base[l][m][pre] = 0.;
1096        qp = 0.;
1097        //
1098        Int_t spos = -1;
1099        Int_t ee = 0;
1100        for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1101          if ( calgood[l][m][e] == 0. && obadmask[l][m][e] == 0 &&  dexy[l][m][e]-calped[l][m][e] < minstrip &&  dexy[l][m][e] > 0. && onlmask[ee] == 0 ) {
1102            minstrip = dexy[l][m][e]-calped[l][m][e];
1103            rms = calthr[l][m][pre];
1104            spos = ee;
1105          };
1106          ee++;
1107          qp += (dexy[l][m][e]-calped[l][m][e]-sbase[l][m][e]);
1108        };
1109        //
1110        if ( debug && l==0 ){
1111          printf("\n BASELINE CALCULATION for view %i pl %i pre %i: \n => minstrip %f rms %f \n => qp = %f \n",l,m,pre,minstrip,rms,qp);
1112        };
1113        if ( minstrip != 100000. ) {
1114          Float_t strip6s = 0.;
1115          for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1116            if ( (dexy[l][m][e]-calped[l][m][e]) >= minstrip && (dexy[l][m][e]-calped[l][m][e]) <= (minstrip+rms) ) {
1117              strip6s += 1.;
1118              base[l][m][pre] += (dexy[l][m][e] - calped[l][m][e]);
1119            };
1120            //
1121            //  compression
1122            //
1123            //      if ( abs((int)(dexy[l][m][e]-calped[l][m][e])) <= (minstrip+rms) ) {
1124            //        dexyc[l][m][e] = 0.;
1125            //      } else {
1126            dexyc[l][m][e] = dexy[l][m][e];
1127              //    };
1128          };
1129          //
1130          if ( strip6s == 1. && process < 1 ){
1131            onlmask[spos] = 1;
1132            process++;
1133            if ( debug ) printf(" Warning, only one strip to calculate baseline: minstrip %f rms %f spos %i l %i m %i pre %i \n",minstrip,rms,spos,l,m,pre);
1134            continue;
1135          };
1136          process += 2;
1137          nst = (Int_t)strip6s;
1138          //
1139          if ( debug ){
1140            printf(" strip6s %f \n",strip6s);
1141          };
1142          //        if ( strip6s >= 9. ){      
1143          if ( (strip6s >= 2. && process == 2) || (strip6s >= 9. && process > 2) ){    
1144            //if ( (strip6s >= 4. && process == 2) || (strip6s >= 9. && process > 2) ){        
1145            Double_t arro = base[l][m][pre]/strip6s;
1146            Float_t deci = 1000.*((float)arro - float(int(arro)));                
1147            if ( deci < 500. ) {
1148              arro = double(int(arro));
1149            } else {
1150              arro = 1. + double(int(arro));
1151            };
1152            base[l][m][pre] = arro;
1153            //
1154            // if too few strips were used to determine the baseline check if it is comparable with the previous event, if not mark it as bad
1155            //
1156            if ( debug && process > 2 ) printf(" AGH low strip value was discarded process %i strip6s %f minstrip %f rms %f spos %i\n",process,strip6s,minstrip,rms,spos);
1157            if ( debug ) printf(" Calculated baseline: base %f sbase-0.02*qp %f \n",base[l][m][pre],(-qp*0.02+sbase[l][m][pre]));
1158            //
1159            if ( strip6s < 4 && base[l][m][pre] > (-0.015*qp+sbase[l][m][pre]) && sbase[l][m][pre] > 0. ){
1160              if ( debug ) printf(" Suspicious calculated baseline: base %f sbase-0.02*qp %f strip6s %i \n",base[l][m][pre],(-qp*0.02+sbase[l][m][pre]),(Int_t)strip6s);
1161              base[l][m][pre] = 31000.;    
1162              for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1163                dexyc[l][m][e] = dexy[l][m][e];
1164              };      
1165            };
1166          } else {
1167          base[l][m][pre] = 31000.;          base[l][m][pre] = 31000.;
1168            for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1169              dexyc[l][m][e] = dexy[l][m][e];
1170            };
1171          };
1172        } else {
1173          process += 2;
1174          base[l][m][pre] = 31000.;
1175          for (Int_t e = pre*16; e < (pre+1)*16 ; e++){
1176            dexyc[l][m][e] = dexy[l][m][e];
1177          };
1178      };      };
1179      };
1180  }  }
1181    
1182  Int_t CaloLevel0::Calibrate(Int_t ei){  Int_t CaloLevel0::Calibrate(Int_t ei){
# Line 369  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1187  Int_t CaloLevel0::Calibrate(Int_t ei){
1187    //    //
1188    // if it was not a selftrigger event, could it ever been a selftrigger event? if so trigty = 3.    // if it was not a selftrigger event, could it ever been a selftrigger event? if so trigty = 3.
1189    //    //
1190      clevel2->nsatstrip = 0.;
1191    Int_t val = 0;    Int_t val = 0;
1192    Int_t del = 1100;    Int_t del = 1100;
1193      for (Int_t sec = 0; sec < 4; sec++){
1194        for (Int_t dsec = 0; dsec < 7; dsec++){
1195          val = (Int_t)de->calselftrig[sec][dsec];
1196          del = delay(val);
1197          clevel2->selfdelay[sec][dsec] = del;
1198        };
1199      };
1200      val = 0;
1201      del = 1100;
1202    if ( clevel2->trigty != 2. ){    if ( clevel2->trigty != 2. ){
1203      Bool_t ck = false;      Bool_t ck = false;
1204      for (Int_t sec = 0; sec < 4; sec++){      for (Int_t sec = 0; sec < 4; sec++){
# Line 406  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1234  Int_t CaloLevel0::Calibrate(Int_t ei){
1234    Float_t ener;    Float_t ener;
1235    Int_t doneb = 0;    Int_t doneb = 0;
1236    Int_t donec = 0;    Int_t donec = 0;
1237    Int_t ck = 0;    Int_t ck[2][22][6];
1238      memset(ck, 0, 2*22*6*sizeof(Int_t));
1239    Int_t ipre = 0;    Int_t ipre = 0;
1240    Int_t ip[3] = {0};    //  Int_t ip[3] = {0};
1241      Int_t ip[3] = {0,0,0};
1242      Int_t ipp = 0;
1243    Float_t base0, base1, base2;    Float_t base0, base1, base2;
1244    base0 = 0.;    base0 = 0.;
1245    base1 = 0.;    base1 = 0.;
1246    base2 = 0.;    base2 = 0.;
1247    Float_t qpre[6] = {0.,0.,0.,0.,0.,0.};    Float_t qpre[2][22][6];
1248      memset(qpre, 0, 2*22*6*sizeof(Float_t));
1249    Float_t ene[96];    Float_t ene[96];
1250    Int_t chdone[4] = {0,0,0,0};    Int_t chdone[4] = {0,0,0,0};
1251    Int_t pe = 0;    Int_t pe = 0;
# Line 421  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1253  Int_t CaloLevel0::Calibrate(Int_t ei){
1253    Float_t ener0 = 0.;    Float_t ener0 = 0.;
1254    Float_t cbase0 = 0.;    Float_t cbase0 = 0.;
1255    Bool_t pproblem = false;    Bool_t pproblem = false;
1256      Bool_t negbase = false;
1257    //    //
1258    Float_t tim = 0.;    Float_t tim = 0.;
1259    Int_t plo = 0;    Int_t plo = 0;
# Line 432  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1265  Int_t CaloLevel0::Calibrate(Int_t ei){
1265    for (Int_t l = 0; l < 2; l++){    for (Int_t l = 0; l < 2; l++){
1266      for (Int_t m = 0; m < 22; m++){      for (Int_t m = 0; m < 22; m++){
1267        //        //
1268        // determine the section number        // determine the section number  
1269        //        //
1270          negbase = false;
1271        se = 5;        se = 5;
1272        if (l == 0 && m%2 == 0) se = 3;        if (l == 0 && m%2 == 0) se = 3;
1273        if (l == 0 && m%2 != 0) se = 2;        if (l == 0 && m%2 != 0) se = 2;
1274        if (l == 1 && m%2 == 0) se = 1;        if (l == 1 && m%2 != 0) se = 1;
1275        if (l == 1 && m%2 != 0) se = 0;                  if (l == 1 && m%2 == 0) se = 0;          
1276        //        //
1277        // determine what kind of event we are going to analyze        // determine what kind of event we are going to analyze
1278        //        //
# Line 459  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1293  Int_t CaloLevel0::Calibrate(Int_t ei){
1293          };          };
1294          clevel2->perr[se] = 0;          clevel2->perr[se] = 0;
1295          if ( de->perror[se] != 0 ){          if ( de->perror[se] != 0 ){
1296            clevel2->perr[se] = 1;            clevel2->perr[se] = (Int_t)de->perror[se];
1297            pe++;            pe++;
1298          };          };
1299          clevel2->swerr[se] = 0;          clevel2->swerr[se] = 0;
# Line 475  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1309  Int_t CaloLevel0::Calibrate(Int_t ei){
1309          pre = -1;          pre = -1;
1310          //          //
1311          for (Int_t nn = 0; nn < 96; nn++){                            for (Int_t nn = 0; nn < 96; nn++){                  
1312            ene[nn] = 0.;            //      ene[nn] = 0.;
1313            dexy[l][m][nn] = de->dexy[l][m][nn] ;            dexy[l][m][nn] = de->dexy[l][m][nn] ;
1314            dexyc[l][m][nn] = de->dexyc[l][m][nn] ;            dexyc[l][m][nn] = de->dexyc[l][m][nn] ;
1315          };          };
# Line 484  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1318  Int_t CaloLevel0::Calibrate(Int_t ei){
1318          //          //
1319          pre = -1;          pre = -1;
1320          cbase0 = 0.;          cbase0 = 0.;
1321            Int_t nstt[2];
1322            Float_t rqp[2];
1323          for (Int_t i = 0; i < 3; i++){          for (Int_t i = 0; i < 3; i++){
1324              nstt[0] = 1000;
1325              nstt[1] = 1000;
1326              rqp[0] = 0.;
1327              rqp[1] = 0.;
1328            for (Int_t j = 0; j < 2; j++){            for (Int_t j = 0; j < 2; j++){
1329              pre = j + i*2;              pre = j + i*2;
1330              //              //
1331              // baseline check and calculation              // baseline check and calculation
1332              //              //
1333              if ( !isRAW ) {              if ( !isRAW ){
1334                base[l][m][pre] = de->base[l][m][pre] ;                  //
1335                  // if it is a compress event with fully transmitted pre try to calculate the baseline
1336                  //
1337                  if ( de->base[l][m][pre] != 0. && de->base[l][m][pre]<31000. ) {
1338                    base[l][m][pre] = de->base[l][m][pre] ;  
1339                  } else {
1340                    FindBaseCompress(l,m,pre,nstt[j],rqp[j]);
1341                  };
1342                cbase0 += base[l][m][pre];                cbase0 += base[l][m][pre];
1343              } else {              } else {
1344                //                //
1345                // if it is a raw event and we haven't checked                // if it is a raw event calculate the baseline.
               // yet, calculate the baseline.  
1346                //                //
1347                FindBaseRaw(l,m,pre);                FindBaseRaw(l,m,pre,nstt[j],rqp[j]);
1348                cbase0 += base[l][m][pre];                cbase0 += base[l][m][pre];
1349              };              };
1350            };            };
1351              //
1352              // if we are able to calculate the baseline with more than 3 strips on one pre and not in the other one choose the pre with more calculated strips
1353              //
1354              if ( nstt[0] < 4 && nstt[1] >= 4 && nstt[0] != 1000 && nstt[1] != 1000 ) base[l][m][pre-1] = 31000.;
1355              if ( nstt[0] >= 4 && nstt[1] < 4 && nstt[0] != 1000 && nstt[1] != 1000 ) base[l][m][pre] = 31000.;
1356              //      //
1357              //      // if we are NOT able to calculate the baseline with more than 3 strips on both pres take the baseline (if any) of the one which has less energy
1358              //      //
1359              //      if ( nstt[0] < 4 && nstt[1] < 4 ){
1360              //        if ( rqp[0] >= rqp[1] ) base[l][m][pre-1] = 31000.;
1361              //        if ( rqp[0] < rqp[1] ) base[l][m][pre] = 31000.;
1362              //      };
1363          };          };
1364          //          //
1365          // run over strips          // run over strips
# Line 512  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1370  Int_t CaloLevel0::Calibrate(Int_t ei){
1370            ip[i] = 0;            ip[i] = 0;
1371            for (Int_t n = i*32 ; n < (i+1)*32 ; n++){                            for (Int_t n = i*32 ; n < (i+1)*32 ; n++){                
1372              if (n%16 == 0) {              if (n%16 == 0) {
               ck = 0;  
1373                done = 0;                done = 0;
1374                doneb = 0;                doneb = 0;
1375                donec = 0;                donec = 0;
1376                pre++;                pre++;
1377                qpre[pre] = 0.;                ck[l][m][pre] = 0;
1378                  qpre[l][m][pre] = 0.;
1379              };              };
1380              //              //
1381              // baseline check and calculation              // baseline check and calculation
# Line 526  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1384  Int_t CaloLevel0::Calibrate(Int_t ei){
1384              //              //
1385              if ( !done ){              if ( !done ){
1386                if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){                if ( (base[l][m][pre] == 31000. || base[l][m][pre] == 0.) ){
1387                  ck = 1;                  ck[l][m][pre] = 1;
1388                  if (pre%2 == 0) {                  if (pre%2 == 0) {
1389                    ip[i] = pre + 1;                    ip[i] = pre + 1;
1390                  } else {                  } else {
1391                    ip[i] = pre - 1;                    ip[i] = pre - 1;
1392                  };                  };
1393                  if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0.) ){                  if ( (base[l][m][ip[i]] == 31000. || base[l][m][ip[i]] == 0. || !crosst ) ){
1394                    //                    //
1395                    ck = 2;                    ck[l][m][pre] = 2;
1396                    if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {                    if ( sbase[l][m][pre] == 31000. || sbase[l][m][pre] == 0. ) {
1397                      ck = 3;                      ck[l][m][pre] = 3;
1398                    };                    };
1399                  };                  };
                 done = 1;  
1400                };                };
1401                  done = 1;
1402              };                                };                  
1403              //              //
1404              // CALIBRATION ALGORITHM              // CALIBRATION ALGORITHM
1405              //              //
1406              if ( !doneb ){              if ( !doneb ){
1407                switch (ck) {                if ( debug ) printf(" ck[l][m][pre] is %i \n",ck[l][m][pre]);
1408                  switch (ck[l][m][pre]) {
1409                case 0:                case 0:
1410                  base0 = base[l][m][pre];                  base0 = base[l][m][pre];
1411                  base2 = calbase[l][m][pre];                  base2 = calbase[l][m][pre];
1412                    if ( debug ) printf(" base0 = base l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,base[l][m][pre],calbase[l][m][pre]);
1413                  break;                  break;
1414                case 1:                case 1:
1415                  base0 = base[l][m][ip[i]];                  base0 = base[l][m][ip[i]];
1416                  base2 = calbase[l][m][ip[i]];                  base2 = calbase[l][m][ip[i]];
1417                    if ( debug ) printf(" base0 = base l%i m%i ip(i)%i = %f base2 = calbase l m ip(i) = %f \n",l,m,ip[i],base[l][m][ip[i]],calbase[l][m][ip[i]]);
1418                  break;                  break;
1419                case 2:                case 2:
1420                  base0 = sbase[l][m][pre];                  base0 = sbase[l][m][pre];
1421                  base2 = calbase[l][m][pre];                  base2 = calbase[l][m][pre];    
1422                    if ( debug ) printf(" base0 = sbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,sbase[l][m][pre],calbase[l][m][pre]);
1423                  break;                  break;
1424                case 3:                case 3:
1425                  base0 = calbase[l][m][pre];                  base0 = calbase[l][m][pre];
1426                  base2 = calbase[l][m][pre];                  base2 = calbase[l][m][pre];
1427                    if ( debug ) printf(" base0 = calbase l%i m%i pre%i = %f base2 = calbase l m pre = %f \n",l,m,pre,calbase[l][m][pre],calbase[l][m][pre]);
1428                  break;                  break;
1429                };                };
1430                base1 = calbase[l][m][pre];                base1 = calbase[l][m][pre];
# Line 570  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1433  Int_t CaloLevel0::Calibrate(Int_t ei){
1433              ener = dexyc[l][m][n];              ener = dexyc[l][m][n];
1434              ener0 += ener;              ener0 += ener;
1435              clevel1->estrip[n][m][l] = 0.;              clevel1->estrip[n][m][l] = 0.;
1436                if ( de->base[l][m][pre] < 0 ) negbase = true;
1437              if ( base0>0 && base0 < 30000. ){              if ( base0>0 && base0 < 30000. ){
1438                if ( !donec && (base0 - base1 + base2) != 0. ){                //
1439                  sbase[l][m][pre] = base0 - base1 + base2;                // save the baseline only if the energy release is "small"
1440                  //
1441                  if ( !donec && (base0 + base1 - base2) != 0. && (n+1)%16==0 ){
1442                    if ( qpre[l][m][pre] < 200. ) sbase[l][m][pre] = base0 + base1 - base2;
1443                  donec = 1;                  donec = 1;
1444                };                };
1445                if ( ener > 0. ){                if ( ener > 0. ){
# Line 580  Int_t CaloLevel0::Calibrate(Int_t ei){ Line 1447  Int_t CaloLevel0::Calibrate(Int_t ei){
1447                  //                  //
1448                  // 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)                  // 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)
1449                  //                  //
1450                  qpre[pre] += clevel1->estrip[n][m][l];                  if ( clevel1->estrip[n][m][l] > 0. ) qpre[l][m][pre] += clevel1->estrip[n][m][l];
1451                };                  //
1452              };                  //
           };  
           if ( crosst ){  
             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] += (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] +=  qpre[ipre] * 0.00478;  
1453                };                };
1454              };              };
1455            };            };
1456          };          };
1457          //          //
1458          if ( ener0 == 0. && cbase0 == 0. && !pproblem && clevel2->perr[se] == 0){          // check if there were problems with 5.7 or glitches in the power supply
1459            //
1460            if ( ((ener0 == 0. && cbase0 == 0.) || negbase ) && !pproblem && clevel2->perr[se] == 0){
1461            if ( verbose ) printf(" L0 entry %i : calorimeter power problems! event marked as bad perr %f swerr %X view %i plane %i \n",ei,de->perror[se],de->stwerr[se],l,m);            if ( verbose ) printf(" L0 entry %i : calorimeter power problems! event marked as bad perr %f swerr %X view %i plane %i \n",ei,de->perror[se],de->stwerr[se],l,m);
1462            pproblem = true;            pproblem = true;
1463            pe++;            pe++;
1464          };          };
1465          //          //
1466          Int_t j4 = -4;        } else {
1467          Int_t jjj = -3;          for (Int_t nn = 0; nn < 96; nn++){                  
1468          Int_t jj = -2;            clevel1->estrip[nn][m][l] = 0.;
1469          for (Int_t j = 0 ; j < 100 ; j++){          };
1470            jj++;        };
1471            jjj++;      };
1472            j4++;    };
1473            if ( j < 96 ) ene[j] = clevel1->estrip[j][m][l];    //
1474            if ( crosst ){    // run over views and planes to apply crosstalk corrections
1475              if ( jj >= 0 && jj < 96 ){    //
1476                if ( jj != 0 && jj != 32 && jj != 64 ) ene[jj-1] += -clevel1->estrip[jj][m][l] * 0.01581;    for (Int_t l = 0; l < 2; l++){
1477                if ( jj != 31 && jj != 63 && jj != 95 ) ene[jj+1] += -clevel1->estrip[jj][m][l] * 0.01581;                              for (Int_t m = 0; m < 22; m++){
1478              };        //
1479              if ( jjj >= 0 && jjj < 96 ){        // determine the section number  
1480                if ( jjj != 0 && jjj != 32 && jjj != 64 ) clevel1->estrip[jjj-1][m][l] +=  -ene[jjj] * 0.01581;        //
1481                if ( jjj != 31 && jjj != 63 && jjj != 95 ) clevel1->estrip[jjj+1][m][l] +=  -ene[jjj] * 0.01581;        se = 5;
1482              };        if (l == 0 && m%2 == 0) se = 3;
1483            };        if (l == 0 && m%2 != 0) se = 2;
1484            if ( j4 >= 0 && j4 < 96 ){        if (l == 1 && m%2 != 0) se = 1;
1485          if (l == 1 && m%2 == 0) se = 0;          
1486          //
1487          // check for any error in the event
1488          //
1489          if ( clevel2->crc[se] == 0 && (clevel1->good2 == 1 || clevel2->trigty >= 2) ){
1490            //
1491            // Cross-talk corrections  
1492            //
1493            if ( crosst ){
1494              //
1495              // energy on silicon ladders
1496              //
1497              Float_t qsi[3];
1498              qsi[0] = qpre[l][m][0]+qpre[l][m][1];
1499              qsi[1] = qpre[l][m][2]+qpre[l][m][3];
1500              qsi[2] = qpre[l][m][4]+qpre[l][m][5];      
1501              //
1502              for ( pre = 1; pre < 6; pre += 2 ){  
1503                Int_t ladder = (pre - 1)/2;
1504              //              //
1505              // NOTICE: THE FOLLOWING LINE EXCLUDE ALL STRIPS FOR WHICH THE RMS*4 IS GREATER THAN 26 !!! <=============== IMPORTANT! =================>              // If the noselfct flag is set the strip doesn't suffer the self crosstalk due to electronics so we must subtract some energy
1506              //              //
1507              if ( obadmask[l][m][j4] == 1 || clevel1->estrip[j4][m][l] <= clevel1->emin || calrms[l][m][j4] > 26 ){              if ( noselfct ){
1508                clevel1->estrip[j4][m][l] = 0.;                for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1509                    ipre = j/16 ;  
1510                    if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] -=  clevel1->estrip[j][m][l] * ctprecor[l][m][ipre];
1511                  };
1512              };              };
1513              //              //
1514              // code and save the energy for each strip in svstrip              // Using the neighbour pre baseline
1515              //              //
1516              if ( clevel1->estrip[j4][m][l] > clevel1->emin ){              if (ck[l][m][pre] == 1 || ck[l][m][pre-1] == 1){
1517                //                //
1518                tim = 100000.;                // pre-amplifier effect on baseline when using the neighbour pre (ck=1)
1519                plo = m;                //
1520                fbi = 0;                if (ck[l][m][pre] == 1){
1521                if ( clevel1->estrip[j4][m][l] > 0.99995 ){                  ipre = pre;
1522                  tim = 10000.;                  ipp = pre - 1;
1523                  plo = m;                } else {
1524                  fbi = 1;                  ipre = pre - 1;
1525                };                  ipp = pre;
               if ( clevel1->estrip[j4][m][l] > 9.9995 ){  
                 tim = 1000.;  
                 plo = 22 + m;  
                 fbi = 1;  
               };  
               if ( clevel1->estrip[j4][m][l] > 99.995 ){  
                 tim = 100.;  
                 plo = 22 + m;  
                 fbi = 0;  
1526                };                };
1527                if ( clevel1->estrip[j4][m][l] > 999.95 ){                Int_t it = 0;
1528                  tim = 10.;                Float_t nqpre = 0.;
1529                  plo = 44 + m;                //
1530                  fbi = 0;                if ( debug ) printf(" CK1 Limit for while: 0.07 \n");
1531                  for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1532                    if ( !ctground ){
1533                      if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += - qpre[l][m][ipp] * ctprecor[l][m][ipp];
1534                    } else {
1535                      if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += - qpre[l][m][ipp] * 0.00478;
1536                    };
1537                    if ( clevel1->estrip[j][m][l] > 0. ) nqpre += clevel1->estrip[j][m][l] ;
1538                };                };
1539                if ( clevel1->estrip[j4][m][l] > 9999.5 ){                qpre[l][m][ipre] = nqpre;
1540                  tim = 1.;                nqpre = 0.;
1541                  plo = 66 + m;                Float_t deltaqpre = qpre[l][m][ipre];          
1542                  fbi = 0;                //
1543                  // these values are empirically determined, usually the routine converge due to deltaqsi and the latest applied correction is based on less than 1 mip
1544                  //
1545                  while ( it < 10 && deltaqpre > 0.07 ){
1546                    nqpre = 0.;
1547                    for (Int_t j = ipre*16 ; j < (ipre+1)*16 ; j++){
1548                      if ( !ctground ){
1549                        if ( debug ) printf(" CK1 pre correction: iteration %i deltaqpre %f  ctprecor %f  TOTAL CORRECTION %f \n",it,deltaqpre,ctprecor[l][m][ipre],deltaqpre * ctprecor[l][m][ipre]);
1550                        if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += deltaqpre * ctprecor[l][m][ipre];
1551                      } else {
1552                        if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] += deltaqpre * 0.00478;
1553                      };
1554                      if ( clevel1->estrip[j][m][l] > 0. ) nqpre += clevel1->estrip[j][m][l] ;
1555                    };
1556                    if ( ctground ) it = 100;
1557                    it++;
1558                    deltaqpre = nqpre - qpre[l][m][ipre];
1559                    if ( debug ) printf(" CK1 BEFORE: qpre %f \n",qpre[l][m][ipre]);  
1560                    qpre[l][m][ipre] = nqpre;
1561                    if ( debug ) printf(" CK1 AFTER: qpre %f \n",qpre[l][m][ipre]);  
1562                };                };
1563                //                //
1564                cle = (Int_t)lroundf(tim*clevel1->estrip[j4][m][l]);              };
1565                //
1566                // No baseline calculation due to high energy release
1567                //
1568                if (ck[l][m][pre] == 2 && ck[l][m][pre-1] == 2){
1569                  //
1570                  // y^
1571                  //  |
1572                  //  |    6 7 8
1573                  //  |    3 4 5  
1574                  //  |    0 1 2
1575                  //  | --------------------------------------> x
1576                //                //
1577                  Int_t si1 = 0;
1578                  Int_t si2 = 0;
1579                  Int_t si3 = 0;
1580                if ( l == 0 ){                if ( l == 0 ){
1581                  //                  if ( ladder == 0 ){
1582                  // +-PPSSmmmm.mmmm                    si1 = 0;
1583                  //                    si2 = 3;
1584                  svstrip[istrip] = fbi*1000000000 + plo*10000000 + j4*100000 + cle;                    si3 = 6;
1585                    };
1586                    if ( ladder == 1 ){
1587                      si1 = 1;
1588                      si2 = 4;
1589                      si3 = 7;
1590                    };
1591                    if ( ladder == 2 ){
1592                      si1 = 2;
1593                      si2 = 5;
1594                      si3 = 8;
1595                    };
1596                } else {                } else {
1597                  svstrip[istrip] = -(fbi*1000000000 + plo*10000000 + j4*100000 + cle);                  if ( ladder == 0 ){
1598                      si1 = 0;
1599                      si2 = 1;
1600                      si3 = 2;
1601                    };
1602                    if ( ladder == 1 ){
1603                      si1 = 3;
1604                      si2 = 4;
1605                      si3 = 5;
1606                    };
1607                    if ( ladder == 2 ){
1608                      si1 = 6;
1609                      si2 = 7;
1610                      si3 = 8;
1611                    };
1612                };                };
1613                //                //
1614                //              if ( ei >= -770 ) printf(" j %i l %i m %i estrip %f \n",j4,l,m,clevel1->estrip[j4][m][l]);                // Find the energy distribution along the considered plane looking at the two sandwiching plane of the other view.
               //              if ( ei >= -770 ) printf(" num lim %i fbi %i tim %f plo %i cle %i \n",numeric_limits<Int_t>::max(),fbi,tim,plo,cle);  
               //              if ( ei >= -770 ) printf(" svstrip %i \n",svstrip[istrip]);  
1615                //                //
1616                istrip++;                Float_t sied[3] = {0.,0.,0.};
1617                  Int_t othv = !l;
1618                  Int_t othpl1 = m - 1;
1619                  Int_t othpl2 = m + 1;
1620                  Float_t oprof[3] = {0.,0.,0.};
1621                  for(Int_t s=0; s<3; s++){
1622                    for(Int_t t=(s*32); t<32*(s + 1); t++){      
1623                      if ( othpl1 > -1 ) {
1624                        oprof[s] += clevel1->estrip[othv][othpl1][t];
1625                      };
1626                      if ( othpl2 < 22 ) {
1627                        oprof[s] += clevel1->estrip[othv][othpl2][t];
1628                      };          
1629                    };
1630                  };
1631                  Float_t otote = fabs(oprof[0]) + fabs(oprof[1]) + fabs(oprof[2]);
1632                  for(Int_t g=0; g<3; g++){
1633                    if ( otote > 0. ){
1634                      sied[g] = fabs(oprof[g])/otote;
1635                    } else {
1636                      sied[g] = 1./3.;
1637                    };
1638                  };
1639                  //
1640                  //
1641                  //
1642                  Int_t it = 0;
1643                  Int_t jpre = 0;
1644                  Float_t nqsi = 0.;
1645                  Float_t snqsi = qsi[ladder];
1646                  Float_t nqpre[2] = {0.,0.};
1647                  Float_t deltaqsi = qsi[ladder];
1648                  Float_t deltaqpre[2];
1649                  deltaqpre[0] = qpre[l][m][pre-1];
1650                  deltaqpre[1] = qpre[l][m][pre];
1651                  //
1652                  if ( debug ) printf(" Limit for while: 0.07 it < 10 \n");
1653                  //
1654                  // these values are empirically determined, usually the routine converge due to deltaqsi and the latest applied correction is based on less than 1 mip
1655                  //
1656                  while ( it < 10 && (deltaqsi > 0.07 || deltaqpre[0] > 0.07 || deltaqpre[1] > 0.07) ){
1657                    nqsi = 0.;
1658                    nqpre[0] = 0.;
1659                    nqpre[1] = 0.;
1660                    for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1661                      ipre = 0;
1662                      if ( j > (ladder*32)+15 ) ipre = 1;
1663                      jpre = j/16 ;
1664                      //
1665                      // Silicon effect on the baseline when using the same pre previous baseline (ck = 2) + pre-amply effect
1666                      //          
1667                      if ( !ctground ){
1668                        if ( debug ) printf(" silicon correction: iteration %i deltaqsi[%i] %f ctsicor %f %f %f sied %f %f %f si %i %i %i TOTAL CORRECTION %f \n",it,ladder,deltaqsi,ctsicor[l][m][si1],ctsicor[l][m][si2],ctsicor[l][m][si3],sied[0],sied[1],sied[2],si1,si2,si3,deltaqsi * (ctsicor[l][m][si1] * sied[0] + ctsicor[l][m][si2] * sied[1] + ctsicor[l][m][si3] * sied[2]));
1669                        if ( debug ) printf(" pre correction: iteration %i deltaqpre[0] %f deltaqpre[1] %f ctprecor %f  TOTAL CORRECTION %f \n",it,deltaqpre[0],deltaqpre[1],ctprecor[l][m][jpre],deltaqpre[ipre] * ctprecor[l][m][jpre]);
1670                        if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] +=  (deltaqsi * (ctsicor[l][m][si1] * sied[0] + ctsicor[l][m][si2] * sied[1] + ctsicor[l][m][si3] * sied[2])/mip[l][m][j]) + deltaqpre[ipre] * ctprecor[l][m][jpre];
1671                      } else {
1672                        if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] +=  0. + qpre[l][m][jpre] * 0.00478; // no correction
1673                      };
1674                      if ( clevel1->estrip[j][m][l] > 0. ) nqsi += clevel1->estrip[j][m][l] ;
1675                      if ( clevel1->estrip[j][m][l] > 0. ) nqpre[ipre] += clevel1->estrip[j][m][l] ;
1676                    };      
1677                    if ( ctground ) it = 100;
1678                    deltaqsi = nqsi-snqsi;
1679                    snqsi = nqsi;
1680                    it++;
1681                    deltaqpre[0] = nqpre[0] - qpre[l][m][pre-1];
1682                    deltaqpre[1] = nqpre[1] - qpre[l][m][pre];
1683                    if ( debug ) printf(" BEFORE: qpre 0 %f qpre 1 %f \n",qpre[l][m][pre-1],qpre[l][m][pre]);  
1684                    qpre[l][m][pre-1] = nqpre[0];
1685                    qpre[l][m][pre] = nqpre[1];    
1686                    if ( debug ) printf(" AFTER: qpre 0 %f qpre 1 %f \n",qpre[l][m][pre-1],qpre[l][m][pre]);  
1687                  };
1688                  //
1689                  //
1690                  //
1691    //            for (Int_t j = ladder*32 ; j < (ladder+1)*32 ; j++){
1692    //              ipre = j/16 ;          
1693    //              //
1694    //              // pre-amplifier effect on baseline when using the same pre previous event baseline (ck=2)
1695    //              //
1696    //              if ( !ctground ){
1697    //                if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] +=  qpre[l][m][ipre] * ctprecor[l][m][ipre];
1698    //              } else {
1699    //                if ( clevel1->estrip[j][m][l] != 0. ) clevel1->estrip[j][m][l] +=  qpre[l][m][ipre] * 0.00478;
1700    //              };
1701    //            };
1702              };              };
1703            };            };
1704          };          };
1705          //        };
1706        } else {        //
1707          for (Int_t nn = 0; nn < 96; nn++){                          Int_t j4 = -4;
1708            clevel1->estrip[nn][m][l] = 0.;        Int_t jjj = -3;
1709          Int_t jj = -2;
1710          Int_t jjpre = -1;
1711          Int_t jjjpre = -1;
1712          memset(ene, 0, 96*sizeof(Float_t));
1713          for (Int_t j = 0 ; j < 100 ; j++){          
1714            jj++;
1715            jjj++;
1716            j4++;
1717            if ( j < 96 ) ene[j] = clevel1->estrip[j][m][l];
1718            if ( crosst ){
1719              //
1720              // "Real" crosstalk effect on the neighbour strips respect to the one which have seen the energy deposit
1721              //
1722              if ( jj >= 0 && jj < 96 ){
1723                if ( !ctground ){
1724                  if ( jj%16 == 0 ) jjpre++;              
1725                  if ( jj != 0 && jj != 32 && jj != 64 && ene[jj-1] != 0. ) ene[jj-1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];
1726                  if ( jj != 31 && jj != 63 && jj != 95 && ene[jj+1] != 0. ) ene[jj+1] += -clevel1->estrip[jj][m][l] * ctneigcor[l][m][jjpre];                      
1727                } else {
1728                  if ( jj != 0 && jj != 32 && jj != 64 && ene[jj-1] != 0. ) ene[jj-1] += -clevel1->estrip[jj][m][l] * 0.01581;
1729                  if ( jj != 31 && jj != 63 && jj != 95 && ene[jj+1] != 0. ) ene[jj+1] += -clevel1->estrip[jj][m][l] * 0.01581;                    
1730                };
1731              };
1732              if ( jjj >= 0 && jjj < 96 ){
1733                if ( !ctground ){
1734                  if ( jjj%16 == 0 ) jjjpre++;            
1735                  if ( jjj != 0 && jjj != 32 && jjj != 64 && clevel1->estrip[jjj-1][m][l] != 0. ) clevel1->estrip[jjj-1][m][l] +=  -ene[jjj] * ctneigcor[l][m][jjjpre];
1736                  if ( jjj != 31 && jjj != 63 && jjj != 95 && clevel1->estrip[jjj+1][m][l] !=0. ) clevel1->estrip[jjj+1][m][l] +=  -ene[jjj] * ctneigcor[l][m][jjjpre];
1737                } else {
1738                  if ( jjj != 0 && jjj != 32 && jjj != 64 && clevel1->estrip[jjj-1][m][l] != 0. ) clevel1->estrip[jjj-1][m][l] +=  -ene[jjj] * 0.01581;
1739                  if ( jjj != 31 && jjj != 63 && jjj != 95 && clevel1->estrip[jjj+1][m][l] != 0. ) clevel1->estrip[jjj+1][m][l] +=  -ene[jjj] * 0.01581;
1740                };
1741              };
1742            };
1743            if ( j4 >= 0 && j4 < 96 ){
1744              //
1745              // CALOLEVEL1 CODING AND FILLING
1746              //
1747              //
1748              // NOTICE: THE FOLLOWING LINE EXCLUDE ALL STRIPS FOR WHICH THE RMS*4 IS GREATER THAN 26 !!! <=============== IMPORTANT! =================>
1749              //
1750              if ( obadmask[l][m][j4] == 1 || clevel1->estrip[j4][m][l] <= clevel1->emin || clevel1->estrip[j4][m][l] <= memin[l][m][j4] || calrms[l][m][j4] > maxrms[l][m] ){
1751                clevel1->estrip[j4][m][l] = 0.;
1752              };
1753              //
1754              // code and save the energy for each strip in svstrip    
1755              //
1756              if ( clevel1->estrip[j4][m][l] > clevel1->emin ){
1757                //
1758                Float_t savel1 = clevel1->estrip[j4][m][l];
1759                //      if ( dexyc[l][m][j4] == 32767. ){
1760                if ( dexyc[l][m][j4] > 32000. ){
1761                  savel1 += 5000.;
1762                  clevel2->nsatstrip += 1.;
1763                };
1764                //
1765                tim = 100000.;
1766                plo = m;
1767                fbi = 0;
1768                if ( savel1 > 0.99995 ){
1769                  tim = 10000.;
1770                  plo = m;
1771                  fbi = 1;
1772                };
1773                if ( savel1 > 9.9995 ){
1774                  tim = 1000.;
1775                  plo = 22 + m;
1776                  fbi = 1;
1777                };
1778                if ( savel1 > 99.995 ){
1779                  tim = 100.;
1780                  plo = 22 + m;
1781                  fbi = 0;
1782                };
1783                if ( savel1 > 999.95 ){
1784                  tim = 10.;
1785                  plo = 44 + m;
1786                  fbi = 0;
1787                };
1788                if ( savel1 > 9999.5 ){
1789                  tim = 1.;
1790                  plo = 66 + m;
1791                  fbi = 0;
1792                };
1793                //
1794                cle = (Int_t)lroundf(tim*savel1);
1795                //
1796                if ( l == 0 ){
1797                  //
1798                  // +-PPSSmmmm.mmmm
1799                  //
1800                  svstrip[istrip] = fbi*1000000000 + plo*10000000 + j4*100000 + cle;
1801                } else {
1802                  svstrip[istrip] = -(fbi*1000000000 + plo*10000000 + j4*100000 + cle);
1803                };
1804                //
1805                istrip++;
1806              };
1807          };          };
1808        };        };
1809          //
1810      };      };
1811    };    };
1812      //
1813      // store goodness flag
1814      //
1815    if ( !pe  ){    if ( !pe  ){
1816      clevel2->good = 1;      clevel2->good = 1;
1817    } else {    } else {
1818      clevel2->good = 0;      clevel2->good = 0;
1819    };    };
1820      //
1821      // done
1822      //
1823    return(0);    return(0);
1824  }  }
1825    
# Line 749  void CaloLevel0::FillTrkVar(CaloLevel2 * Line 1868  void CaloLevel0::FillTrkVar(CaloLevel2 *
1868      memcpy(t_ca->tibar,clevel2->cibar,sizeof(clevel2->cibar));      memcpy(t_ca->tibar,clevel2->cibar,sizeof(clevel2->cibar));
1869      memcpy(t_ca->tbar,clevel2->cbar,sizeof(clevel2->cbar));      memcpy(t_ca->tbar,clevel2->cbar,sizeof(clevel2->cbar));
1870      memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));      memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1871        memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1872      ca->varcfit[2] = clevel2->varcfit[0];      ca->varcfit[2] = clevel2->varcfit[0];
1873      ca->varcfit[3] = clevel2->varcfit[1];      ca->varcfit[3] = clevel2->varcfit[1];
1874      ca->npcfit[2] = clevel2->npcfit[0];      ca->npcfit[2] = clevel2->npcfit[0];
# Line 789  void CaloLevel0::FillCommonVar(CaloLevel Line 1909  void CaloLevel0::FillCommonVar(CaloLevel
1909    memcpy(ca->swerr,clevel2->swerr,sizeof(clevel2->swerr));    memcpy(ca->swerr,clevel2->swerr,sizeof(clevel2->swerr));
1910    memcpy(ca->crc,clevel2->crc,sizeof(clevel2->crc));    memcpy(ca->crc,clevel2->crc,sizeof(clevel2->crc));
1911    ca->nstrip = (Int_t)clevel2->nstrip;    ca->nstrip = (Int_t)clevel2->nstrip;
1912      ca->nsatstrip = (Int_t)clevel2->nsatstrip;
1913    ca->qtot = clevel2->qtot;    ca->qtot = clevel2->qtot;
1914    //  ca->impx = clevel2->impx;    //  ca->impx = clevel2->impx;
1915    //  ca->impy = clevel2->impy;    //  ca->impy = clevel2->impy;
# Line 801  void CaloLevel0::FillCommonVar(CaloLevel Line 1922  void CaloLevel0::FillCommonVar(CaloLevel
1922    ca->selen = clevel2->selen;    ca->selen = clevel2->selen;
1923    memcpy(ca->qq,clevel2->qq,sizeof(clevel2->qq));    memcpy(ca->qq,clevel2->qq,sizeof(clevel2->qq));
1924    memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));    memcpy(ca->planemax,clevel2->planemax,sizeof(clevel2->planemax));
1925      memcpy(ca->selfdelay,clevel2->selfdelay,sizeof(clevel2->selfdelay));
1926    ca->varcfit[0] = clevel2->varcfit[0];    ca->varcfit[0] = clevel2->varcfit[0];
1927    ca->varcfit[1] = clevel2->varcfit[1];    ca->varcfit[1] = clevel2->varcfit[1];
1928    ca->npcfit[0] = clevel2->npcfit[0];    ca->npcfit[0] = clevel2->npcfit[0];
# Line 824  void CaloLevel0::ClearStructs(){ Line 1946  void CaloLevel0::ClearStructs(){
1946    ClearCommonVar();    ClearCommonVar();
1947  }  }
1948    
1949    void CaloLevel0::Delete(Option_t *t){
1950      if ( de ) delete de;
1951      delete this;
1952    }
1953    
1954    
1955  void CaloLevel0::RunClose(){  void CaloLevel0::RunClose(){
1956    l0tr->Delete();    l0tr->Delete();
1957    ClearStructs();    ClearStructs();
# Line 832  void CaloLevel0::RunClose(){ Line 1960  void CaloLevel0::RunClose(){
1960    memset(dexyc, 0, 2*22*96*sizeof(Float_t));    memset(dexyc, 0, 2*22*96*sizeof(Float_t));
1961    memset(base, 0, 2*22*6*sizeof(Float_t));    memset(base, 0, 2*22*6*sizeof(Float_t));
1962    memset(sbase, 0, 2*22*6*sizeof(Float_t));    memset(sbase, 0, 2*22*6*sizeof(Float_t));
1963      memset(ctprecor, 0, 2*22*6*sizeof(Float_t));
1964      memset(ctsicor, 0, 2*22*9*sizeof(Float_t));
1965      memset(ctneigcor, 0, 2*22*6*sizeof(Float_t));
1966    //    //
1967  }  }
1968    
# Line 876  void CaloLevel0::ClearCommonVar(){ Line 2007  void CaloLevel0::ClearCommonVar(){
2007    clevel2->wartrig = 0.;    clevel2->wartrig = 0.;
2008    clevel2->good = 0;    clevel2->good = 0;
2009    clevel2->nstrip = 0.;    clevel2->nstrip = 0.;
2010      clevel2->nsatstrip = 0.;
2011    clevel2->qtot = 0.;    clevel2->qtot = 0.;
2012    //  clevel2->impx = 0.;    //  clevel2->impx = 0.;
2013    //  clevel2->impy = 0.;    //  clevel2->impy = 0.;
# Line 891  void CaloLevel0::ClearCommonVar(){ Line 2023  void CaloLevel0::ClearCommonVar(){
2023    memset(clevel2->varcfit, 0, 4*sizeof(Float_t));    memset(clevel2->varcfit, 0, 4*sizeof(Float_t));
2024    memset(clevel2->npcfit, 0, 4*sizeof(Int_t));    memset(clevel2->npcfit, 0, 4*sizeof(Int_t));
2025    memset(clevel2->planemax, 0, 2*sizeof(Int_t));    memset(clevel2->planemax, 0, 2*sizeof(Int_t));
2026      memset(clevel2->selfdelay, 0, 4*7*sizeof(Int_t));
2027    memset(clevel2->fmode, 0, 2*sizeof(Int_t));    memset(clevel2->fmode, 0, 2*sizeof(Int_t));
2028    memset(clevel2->cibar, 0, 2*22*sizeof(Int_t));    memset(clevel2->cibar, 0, 2*22*sizeof(Int_t));
2029    memset(clevel2->cbar, 0, 2*22*sizeof(Float_t));    memset(clevel2->cbar, 0, 2*22*sizeof(Float_t));
# Line 947  void CaloLevel0::ClearCalibVals(Int_t s) Line 2080  void CaloLevel0::ClearCalibVals(Int_t s)
2080    return;    return;
2081  }  }
2082    
2083  Int_t CaloLevel0::Update(TSQLServer *dbc, UInt_t atime, Int_t s){  Int_t CaloLevel0::Update(GL_TABLES *glt, UInt_t atime, Int_t s){
2084    //    //
2085      const TString host = glt->CGetHost();
2086      const TString user = glt->CGetUser();
2087      const TString psw = glt->CGetPsw();
2088      TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
2089      if ( !dbc->IsConnected() ) throw -116;
2090      stringstream myquery;
2091      myquery.str("");
2092      myquery << "SET time_zone='+0:00'";
2093      dbc->Query(myquery.str().c_str());
2094    Int_t sgnl = 0;    Int_t sgnl = 0;
2095    //    //
2096    GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();    GL_CALO_CALIB *glcalo = new GL_CALO_CALIB();
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.25
  ViewVC Help
Powered by ViewVC 1.1.23