/[PAMELA software]/tof/flight/ToFReprocessing/ToFreproc.cpp
ViewVC logotype

Diff of /tof/flight/ToFReprocessing/ToFreproc.cpp

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

revision 1.1 by mocchiut, Fri Nov 28 14:58:50 2008 UTC revision 1.7 by mocchiut, Fri Jan 17 15:16:24 2014 UTC
# Line 44  extern "C" int tofl2com(); Line 44  extern "C" int tofl2com();
44  extern "C" int toftrk();  extern "C" int toftrk();
45  ///extern "C" int toftrk(float);  ///extern "C" int toftrk(float);
46  #define rdtofcal rdtofcal_  #define rdtofcal rdtofcal_
47  extern "C" int rdtofcal(char [], int *);  //extern "C" int rdtofcal(char [], int *);
48    extern "C" int rdtofcal(const char *, int *);
49    
50  using namespace std;  using namespace std;
51    
# Line 57  TString   DIR; Line 58  TString   DIR;
58  TString   OUTF;  TString   OUTF;
59  TString   LIST;  TString   LIST;
60  TString   OPTIONS;  TString   OPTIONS;
61    TString   CALIBF;
62    
63  PamLevel2 *pam_event = NULL;  PamLevel2 *pam_event = NULL;
64    
# Line 72  Int_t Loop(TString ddir,TString list, TS Line 74  Int_t Loop(TString ddir,TString list, TS
74    //  read input file/list    //  read input file/list
75    //  --------------------    //  --------------------
76    pam_event = new PamLevel2(ddir,list,options);    pam_event = new PamLevel2(ddir,list,options);
77      pam_event->SetSELLI(2);
78    TTree::SetMaxTreeSize(1000*Long64_t(2000000000));    TTree::SetMaxTreeSize(1000*Long64_t(2000000000));
79   //    //
80    TString outfile_t =0;    TString outfile_t = "";
81    outfile_t = OUTF;    outfile_t = OUTF;
82    outfile_t.Append("b.root");    outfile_t.Append(".root");
83    TFile *outft = (TFile*)gROOT->FindObject(outfile_t);    TFile *outft = (TFile*)gROOT->FindObject(outfile_t);
84    if (outft) outft->Close();    if (outft) outft->Close();
85    outft = new TFile(outfile_t,"RECREATE");    outft = new TFile(outfile_t,"RECREATE");
# Line 88  Int_t Loop(TString ddir,TString list, TS Line 91  Int_t Loop(TString ddir,TString list, TS
91    //    //
92    outft->cd();    outft->cd();
93    ToFLevel2 *tof = new ToFLevel2();    ToFLevel2 *tof = new ToFLevel2();
94      ToFdEdx *tofdedx = new ToFdEdx();
95    TTree *toft = new TTree("ToF","PAMELA Level2 ToF data");    TTree *toft = new TTree("ToF","PAMELA Level2 ToF data");
96    toft->SetAutoSave(900000000000000LL);              toft->SetAutoSave(900000000000000LL);          
97    tof->Set();//ELENA **TEMPORANEO?**    tof->Set();//ELENA **TEMPORANEO?**
# Line 97  Int_t Loop(TString ddir,TString list, TS Line 101  Int_t Loop(TString ddir,TString list, TS
101    //    //
102    Int_t ntrkentry = 0;    Int_t ntrkentry = 0;
103    Int_t npmtentry = 0;    Int_t npmtentry = 0;
104      Float_t xleft=0;
105      Float_t xright=0;
106      Float_t yleft=0;
107      Float_t yright=0;
108    //    //
109    ULong64_t nevents = pam_event->GetEntries();    ULong64_t nevents = pam_event->GetEntries();
110    printf("\n\n Running on %llu events \n\n",nevents);    printf("\n\n Running on %llu events \n\n",nevents);
# Line 145  Int_t Loop(TString ddir,TString list, TS Line 153  Int_t Loop(TString ddir,TString list, TS
153    for(ULong64_t iev=0; iev<nevents; iev++){          for(ULong64_t iev=0; iev<nevents; iev++){      
154    
155    
156  //==================================================================      //==================================================================
157    
158    
159      //  for(ULong64_t iev=0; iev<50; iev++){            //  for(ULong64_t iev=0; iev<50; iev++){      
# Line 166  Int_t Loop(TString ddir,TString list, TS Line 174  Int_t Loop(TString ddir,TString list, TS
174          //      printf(" rhtime %u myrun %i mysrun %i \n",pam_event->GetRunInfo()->RUNHEADER_TIME,(int)myrun,(int)mysrun);          //      printf(" rhtime %u myrun %i mysrun %i \n",pam_event->GetRunInfo()->RUNHEADER_TIME,(int)myrun,(int)mysrun);
175          mysrun = myrun;          mysrun = myrun;
176          //          //
177          Int_t error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,201,dbc); // parameters stored in DB in GL_PRAM table  
178    
179    
180    
181            Int_t error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,204,dbc); // parameters stored in DB in GL_PRAM table
182            if ( error<0 ) {
183              return(1);
184            };
185            //
186            tofdedx->ReadParAtt((glparam->PATH+glparam->NAME).Data());
187            printf(" Reading Attenuation file: %s \n",(glparam->PATH+glparam->NAME).Data());        
188    
189            //
190            error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,205,dbc); // parameters stored in DB in GL_PRAM table
191            if ( error<0 ) {
192              return(1);
193            };
194            //
195            tofdedx->ReadParPos((glparam->PATH+glparam->NAME).Data());
196            printf(" Reading desaturation1 file: %s \n",(glparam->PATH+glparam->NAME).Data());
197    
198            //
199            error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,206,dbc); // parameters stored in DB in GL_PRAM table
200            if ( error<0 ) {
201              return(1);
202            };
203            //
204            tofdedx->ReadParBBneg((glparam->PATH+glparam->NAME).Data());
205            printf(" Reading BBneg file: %s \n",(glparam->PATH+glparam->NAME).Data());
206    
207            //
208            error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,207,dbc); // parameters stored in DB in GL_PRAM table
209            if ( error<0 ) {
210              return(1);
211            };
212            //
213            tofdedx->ReadParBBpos((glparam->PATH+glparam->NAME).Data());
214            printf(" Reading BBpos file: %s \n",(glparam->PATH+glparam->NAME).Data());
215    
216            //
217            error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,208,dbc); // parameters stored in DB in GL_PRAM table
218            if ( error<0 ) {
219              return(1);
220            };
221            //
222            tofdedx->ReadParDesatBB((glparam->PATH+glparam->NAME).Data());
223            printf(" Reading desaturation2 file: %s \n",(glparam->PATH+glparam->NAME).Data());
224            
225            tofdedx->CheckConnectors(pam_event->GetRunInfo()->RUNHEADER_TIME,glparam,dbc);
226            
227            
228            error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,201,dbc); // parameters stored in DB in GL_PRAM table
229          if ( error<0 ) {          if ( error<0 ) {
230            return(1);            return(1);
231          };          };
# Line 175  Int_t Loop(TString ddir,TString list, TS Line 234  Int_t Loop(TString ddir,TString list, TS
234          //          //
235          if ( (UInt_t)glparam->TO_TIME != (UInt_t)4294967295UL ) defcal = false;          if ( (UInt_t)glparam->TO_TIME != (UInt_t)4294967295UL ) defcal = false;
236          //          //
237          Int_t nlen = (Int_t)(glparam->PATH+glparam->NAME).Length();          TString pippo=(glparam->PATH+glparam->NAME).Data();
238          rdtofcal((char *)(glparam->PATH+glparam->NAME).Data(),&nlen);          //Int_t nlen = (Int_t)(glparam->PATH+glparam->NAME).Length();
239            //rdtofcal((char *)(glparam->PATH+glparam->NAME).Data(),&nlen);
240            Int_t nlen = (Int_t)pippo.Length();
241            rdtofcal(pippo.Data(),&nlen);
242          //          //
       };        
       //  
   
243    
244  //================================================================          tofdedx->Clear();
 //==================================================================  
 //---  Define absolute time  
        Float_t tabs=pam_event->GetOrbitalInfo()->absTime;  
   
 //==================================================================  
 Float_t dedx_corr_m[2000][48],dedx_corr[48];  
 Double_t mtime[2000],t1,t2,tm;  
 Float_t yhelp1,yhelp2,slope,inter,thelp1,thelp2;  
 Float_t xmean1,xwidth1;  
   
 Int_t ical,ii,j,jj;  
   
 if (iev==0) {  
   
 ical=0;  // counter set to zero if first-time reading  
   
 //-----------------------------------------------------------  
 // Here I  read the dEdx_korr parameters  
 //-----------------------------------------------------------  
   
 jj=0;  
   
 ifstream fin("adcmonitor.7th.100k.cut.dat");  
   
 while (! fin.eof()) {  
 fin>>t1>>tm>>t2;  
 //cout<<jj<<" "<<tm<<endl;  
 cout << setiosflags(ios::fixed)  << setw(10) << setprecision(3) << tm << endl;  
 mtime[jj]=tm;  
 for (ii=0; ii<48;ii++) {  
 fin>>j>>xmean1>>xwidth1;  
 dedx_corr_m[jj][ii]=xmean1;  
                        }  
 jj=jj+1;  
 }  
245    
246  fin.close();        };      
247          //
248    
249    
250  while (tabs > mtime[ical]) {        //================================================================
251  ical = ical+1;        //==================================================================
252            }        //---  Define absolute time
253  ical = ical-1;        UInt_t tabs=pam_event->GetOrbitalInfo()->absTime;
254  cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;        
255  } // if iev==0...        if ( !(iev%100000) ) printf(" ATIME %u re %u \n",(Int_t)tabs,(UInt_t)iev);
256    
257  //==================================================================        //==================================================================
258  //==  End of first time reading & filling the array        Float_t dedx_corr_m[2000][48];//,dedx_corr[48];
259  //==================================================================        Double_t mtime[2000],t1,t2,tm;
260          //      Float_t yhelp1,yhelp2,slope,inter,thelp1,thelp2;
261  //==================================================================        Float_t yhelp1,yhelp2;
262  //==  if time is outside time limits:        Float_t xmean1,xwidth1;
263  //==================================================================  
264          Int_t ical=0;
265  if (tabs<mtime[ical] || tabs>mtime[ical+1]) {        Int_t ii,j,jj;
266  cout<<"Checking Time Limits!"<<endl;  
267  ical=0;        if (iev==0) {
268  while (tabs > mtime[ical]) {  
269  ical = ical+1;          ical=0;  // counter set to zero if first-time reading
270            }  
271  ical = ical-1;          //-----------------------------------------------------------
272  cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;          // Here I  read the dEdx_korr parameters
273                                             };          //-----------------------------------------------------------
274    
275  //==================================================================          jj=0;
276  //== interpolate betwen time limits          printf(" READING NEW CALIBRATION FILE: %s \n",CALIBF.Data());
277  //==================================================================  
278            ifstream fin(CALIBF.Data());
279         thelp1 = mtime[ical];          //cout << "topolino" << endl;
280         thelp2 = mtime[ical+1];          while (! fin.eof()) {
281              fin>>t1>>tm>>t2;
282         for (ii=0; ii<48;ii++) {            //cout<<jj<<" "<<tm<<endl;
283         yhelp1 = dedx_corr_m[ical][ii];            //cout << setiosflags(ios::fixed)  << setw(10) << setprecision(3) << tm << endl;
284         yhelp2 = dedx_corr_m[ical+1][ii];            mtime[jj]=tm;
285         slope  = (yhelp2-yhelp1)/(thelp2-thelp1);            for (ii=0; ii<48;ii++) {
286         inter  = yhelp1 - slope*thelp1;              fin>>j>>xmean1>>xwidth1;
287         dedx_corr[ii] = slope*tabs + inter;              dedx_corr_m[jj][ii]=xmean1;
288  //       if (ii==0) cout<<thelp1<<" "<<thelp2<<" "<<tabs<<" "<<yhelp1<<" "<<yhelp2<<" "<<dedx_corr[0]<<endl;            }
289                         }            jj=jj+1;
290              //      printf(" kk %i \n",jj);
291            }
292            //      printf(" 1ical %i \n",ical);
293            //cout << "pippo" << endl;
294            fin.close();
295    
296    
297            while (tabs<mtime[ical] || tabs > mtime[ical+1]) {
298              //      printf(" ical %i \n",ical);
299              ical = ical+1;
300            }
301            //      ical = ical-1;
302            //cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;
303          } // if iev==0...
304    
305          //==================================================================
306          //==  End of first time reading & filling the array
307          //==================================================================
308    
309          //==================================================================
310          //==  if time is outside time limits:
311          //==================================================================
312    
313          if (tabs<mtime[ical] || tabs>mtime[ical+1]) {
314            cout<<"Checking Time Limits!"<<endl;
315            ical=0;
316            while (tabs > mtime[ical+1] || tabs<mtime[ical]) {
317              ical = ical+1;
318            }
319            //      ical = ical-1;
320            //cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;
321          };
322    
323  //================================================================        //      printf(" 2ical %i \n",ical);
324  //================================================================        //==================================================================
325          //== interpolate betwen time limits
326          //==================================================================
327    
328          //      thelp1 = mtime[ical];
329          //      thelp2 = mtime[ical+1];
330    
331          for (ii=0; ii<48;ii++) {
332            yhelp1 = fabs(dedx_corr_m[ical][ii]);
333            //      yhelp1 = 6.;
334            if ( yhelp1 < 0.1 ) yhelp1 = 4.;
335            yhelp2 = fabs(dedx_corr_m[ical+1][ii]);
336            //      yhelp2 = 6.;
337            if ( yhelp2 < 0.1 ) yhelp2 = 4.;
338            //      slope  = (yhelp2-yhelp1)/(thelp2-thelp1);
339            //      inter  = yhelp1 - slope*thelp1;
340            //      dedx_corr[ii] = slope*tabs + inter;
341            //       if (ii==0) cout<<thelp1<<" "<<thelp2<<" "<<tabs<<" "<<yhelp1<<" "<<yhelp2<<" "<<dedx_corr[0]<<endl;
342          }
343    
344          //================================================================
345          //================================================================
346          //      printf("cpippo \n");
347    
348        // process tof data        // process tof data
349        //        //
350        for (Int_t hh=0; hh<12;hh++){        for (Int_t hh=0; hh<12;hh++){
351          for (Int_t kk=0; kk<4;kk++){          for (Int_t kk=0; kk<4;kk++){
352                  adc[kk][hh] = 4095;            adc[kk][hh] = 4095;
353                  tdc[kk][hh] = 4095;            tdc[kk][hh] = 4095;
354                  tdcc[kk][hh] = 4095.;            tdcc[kk][hh] = 4095.;
355                  tofinput_.adc[hh][kk] = 4095;            tofinput_.adc[hh][kk] = 4095;
356                  tofinput_.tdc[hh][kk] = 4095;            tofinput_.tdc[hh][kk] = 4095;
357          };          };
358        };        };
359  //      memset(adc, 0, 12*4*sizeof(Int_t));        //      memset(adc, 0, 12*4*sizeof(Int_t));
360  //      memset(tdc, 0, 12*4*sizeof(Int_t));        //      memset(tdc, 0, 12*4*sizeof(Int_t));
361        Int_t gg = 0;        Int_t gg = 0;
362        Int_t hh = 0;        Int_t hh = 0;
363        Int_t adcf[48];        Int_t adcf[48];
# Line 289  cout<<"abs time "<<tabs<<" limit low "<< Line 367  cout<<"abs time "<<tabs<<" limit low "<<
367        for (Int_t pm=0; pm < tofl2->ntrk() ; pm++){        for (Int_t pm=0; pm < tofl2->ntrk() ; pm++){
368          ToFTrkVar *ttf = tofl2->GetToFTrkVar(pm);          ToFTrkVar *ttf = tofl2->GetToFTrkVar(pm);
369          for ( Int_t nc=0; nc < ttf->npmttdc; nc++){          for ( Int_t nc=0; nc < ttf->npmttdc; nc++){
370                  if ( (ttf->tdcflag).At(nc) != 0 ) tdcf[(ttf->pmttdc).At(nc)] = 1;            if ( (ttf->tdcflag).At(nc) != 0 ) tdcf[(ttf->pmttdc).At(nc)] = 1;
371          };          };
372          for ( Int_t nc=0; nc < ttf->npmtadc; nc++){          for ( Int_t nc=0; nc < ttf->npmtadc; nc++){
373                  if ( (ttf->adcflag).At(nc) != 0 ) adcf[(ttf->pmtadc).At(nc)] = 1;            if ( (ttf->adcflag).At(nc) != 0 ) adcf[(ttf->pmtadc).At(nc)] = 1;
374          };          };
375        };        };
376        //        //
# Line 300  cout<<"abs time "<<tabs<<" limit low "<< Line 378  cout<<"abs time "<<tabs<<" limit low "<<
378          ToFPMT *pmt = tofl2->GetToFPMT(pm);          ToFPMT *pmt = tofl2->GetToFPMT(pm);
379          tofl2->GetPMTIndex(pmt->pmt_id, gg, hh);          tofl2->GetPMTIndex(pmt->pmt_id, gg, hh);
380          if ( adcf[pmt->pmt_id] == 0 ){          if ( adcf[pmt->pmt_id] == 0 ){
381                  tofinput_.adc[gg][hh] = (int)pmt->adc;            tofinput_.adc[gg][hh] = (int)pmt->adc;
382                  adc[hh][gg] = (int)pmt->adc;            adc[hh][gg] = (int)pmt->adc;
383          };          };
384          if ( tdcf[pmt->pmt_id] == 0 ){          if ( tdcf[pmt->pmt_id] == 0 ){
385                  tofinput_.tdc[gg][hh] = (int)pmt->tdc;            tofinput_.tdc[gg][hh] = (int)pmt->tdc;
386                  tdc[hh][gg] = (int)pmt->tdc;            tdc[hh][gg] = (int)pmt->tdc;
387          };          };
388          tdcc[hh][gg] = (float)pmt->tdc_tw;          tdcc[hh][gg] = (float)pmt->tdc_tw;
389          //        Int_t pppid = tof->GetPMTid(hh,gg);          //        Int_t pppid = tof->GetPMTid(hh,gg);
390          //      printf(" pm %i pmt_id %i pppid %i hh %i gg %i tdcc %f tdc %f adc %f \n",pm,pmt->pmt_id,pppid,hh,gg,pmt->tdc_tw,pmt->tdc,pmt->adc);          //      printf(" pm %i pmt_id %i pppid %i hh %i gg %i tdcc %f tdc %f adc %f \n",pm,pmt->pmt_id,pppid,hh,gg,pmt->tdc_tw,pmt->tdc,pmt->adc);
391        };        };
392          for (Int_t hh=0; hh<12;hh++){
393            for (Int_t kk=0; kk<4;kk++){      
394  //      for (Int_t pm=0; pm <48 ; pm++){            tofdedx->Init(kk,hh,adc[kk][hh]);
395  //      tof->GetPMTIndex(pm, gg, hh);          };
396  //      tofinput_.tdc[hh][gg] = (int)500.;        };
397  //      tofinput_.adc[hh][gg] = (int)500.;        //      for (Int_t pm=0; pm <48 ; pm++){
398  //      tdc[hh][gg] = (int)500.;        //        tof->GetPMTIndex(pm, gg, hh);
399  //      adc[hh][gg] = (int)500.;        //        tofinput_.tdc[hh][gg] = (int)500.;
400  //      //      printf(" hh %i gg %i tdc %f adc %f \n",hh,gg,pmt->tdc,pmt->adc);        //        tofinput_.adc[hh][gg] = (int)500.;
401  //    };        //        tdc[hh][gg] = (int)500.;
402          //        adc[hh][gg] = (int)500.;
403          //        //      printf(" hh %i gg %i tdc %f adc %f \n",hh,gg,pmt->tdc,pmt->adc);
404          //    };
405        //        //
406        for (Int_t hh=0; hh<5;hh++){        for (Int_t hh=0; hh<5;hh++){
407          tofinput_.patterntrig[hh]=pam_event->GetTrigLevel2()->patterntrig[hh];          tofinput_.patterntrig[hh]=pam_event->GetTrigLevel2()->patterntrig[hh];
# Line 368  cout<<"abs time "<<tabs<<" limit low "<< Line 449  cout<<"abs time "<<tabs<<" limit low "<<
449          t_tof->beta[kk] = tofoutput_.betatof_a[kk];              t_tof->beta[kk] = tofoutput_.betatof_a[kk];    
450        }        }
451        //        //
       t_tof->npmtadc = 0;  
       for (Int_t hh=0; hh<12;hh++){  
         for (Int_t kk=0; kk<4;kk++){  
           if ( tofoutput_.adctof_c[hh][kk] < 1000 ){  
             //      t_tof->dedx.AddAt(tofoutput_.adctof_c[hh][kk],t_tof->npmtadc); // EMILIANO  
             pmt_id = tof->GetPMTid(kk,hh);  
             t_tof->dedx.AddAt((tofoutput_.adctof_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO  
             t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);  
             t_tof->adcflag.AddAt(tofoutput_.adcflagtof[hh][kk],t_tof->npmtadc); // gf: Jan 09/07  
             t_tof->npmtadc++;  
           };  
         };  
       };  
       //  
452        memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));        memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
453        memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));        memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
454        memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));        memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
455        memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));        memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
456        //        //      {
457        new(t[ntrkentry]) ToFTrkVar(*t_tof);  //      Float_t xtof_temp[6]={0.,t_tof->xtofpos[0],t_tof->xtofpos[1],0.,0.,t_tof->xtofpos[2]};
458        ntrkentry++;        //      Float_t ytof_temp[6]={t_tof->ytofpos[0],0.,0.,t_tof->ytofpos[1],t_tof->ytofpos[2],0.};
459        t_tof->Clear();  //      tofdedx->Process(pam_event->GetOrbitalInfo()->absTime,t_tof->beta[12], (Float_t *)xtof_temp,(Float_t *)ytof_temp);
460        //  //       }
461        //  //       t_tof->npmtadc = 0;
462        //  
463        t_pmt->Clear();  
464        //  //       for (Int_t hh=0; hh<12;hh++){
465    //      for (Int_t kk=0; kk<4;kk++){
466    //        pmt_id = tof->GetPMTid(kk,hh);
467    //        if ( tofdedx->GetdEdx_pmt(pmt_id)>-1. ){
468    //          t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*36./pow(dedx_corr[pmt_id],2)),t_tof->npmtadc);
469    //          t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
470    //          t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
471    //          t_tof->npmtadc++;
472    //        };
473    //      };
474    //       };
475         {
476    
477            Float_t xtof_temp[6]={100.,100.,100.,100.,100.,100.};
478            Float_t ytof_temp[6]={100.,100.,100.,100.,100.,100.};
479    
480            if(t_tof->xtofpos[0]<100. && t_tof->ytofpos[0]<100.){
481              xtof_temp[1]=t_tof->xtofpos[0];
482              ytof_temp[0]=t_tof->ytofpos[0];
483            }else if(t_tof->xtofpos[0]>=100. && t_tof->ytofpos[0]<100.){
484              ytof_temp[0]=t_tof->ytofpos[0];
485              tof->GetPaddleGeometry(0,(Int_t)log2(tof->tof_j_flag[0]),xleft, xright, yleft, yright);
486              xtof_temp[1]=xleft+2.55;
487            }else if(t_tof->ytofpos[0]>=100. && t_tof->xtofpos[0]<100.){
488              xtof_temp[1]=t_tof->xtofpos[0];
489              tof->GetPaddleGeometry(1,(Int_t)log2(tof->tof_j_flag[1]),xleft, xright, yleft, yright);
490              ytof_temp[0]=yleft+2.75;
491            }
492    
493            if(t_tof->xtofpos[1]<100. && t_tof->ytofpos[1]<100.){
494              xtof_temp[2]=t_tof->xtofpos[1];
495              ytof_temp[3]=t_tof->ytofpos[1];
496            }else if(t_tof->xtofpos[1]>=100. && t_tof->ytofpos[1]<100.){
497              ytof_temp[3]=t_tof->ytofpos[1];
498              tof->GetPaddleGeometry(3,(Int_t)log2(tof->tof_j_flag[3]),xleft, xright, yleft, yright);
499              xtof_temp[2]=xleft+4.5;
500            }else if(t_tof->ytofpos[1]>=100. && t_tof->xtofpos[1]<100.){
501              xtof_temp[2]=t_tof->xtofpos[1];
502              tof->GetPaddleGeometry(2,(Int_t)log2(tof->tof_j_flag[2]),xleft, xright, yleft, yright);
503              ytof_temp[3]=yleft+3.75;
504            }
505    
506            if(t_tof->xtofpos[2]<100. && t_tof->ytofpos[2]<100.){
507              xtof_temp[5]=t_tof->xtofpos[2];
508              ytof_temp[4]=t_tof->ytofpos[2];
509            }else if(t_tof->xtofpos[2]>=100. && t_tof->ytofpos[2]<100.){
510              ytof_temp[4]=t_tof->ytofpos[2];
511              tof->GetPaddleGeometry(4,(Int_t)log2(tof->tof_j_flag[4]),xleft, xright, yleft, yright);
512              xtof_temp[5]=xleft+3;
513            }else if(t_tof->ytofpos[2]>=100. && t_tof->xtofpos[2]<100.){
514              xtof_temp[5]=t_tof->xtofpos[2];
515              tof->GetPaddleGeometry(5,(Int_t)log2(tof->tof_j_flag[5]),xleft, xright, yleft, yright);
516              ytof_temp[4]=yleft+2.5;
517            }
518    //       Float_t xtof_temp[6]={0.,t_tof->xtofpos[0],t_tof->xtofpos[1],0.,0.,t_tof->xtofpos[2]};
519    //       Float_t ytof_temp[6]={t_tof->ytofpos[0],0.,0.,t_tof->ytofpos[1],t_tof->ytofpos[2],0.};
520    //      tofdedx->Process(atime,t_tof->beta[12], (Float_t *)xtof_temp,(Float_t *)ytof_temp);
521            tofdedx->Process(pam_event->GetOrbitalInfo()->absTime,t_tof->beta[12], (Float_t *)xtof_temp,(Float_t *)ytof_temp);
522            t_tof->npmtadc = 0;
523            for (Int_t hh=0; hh<12;hh++){
524              for (Int_t kk=0; kk<4;kk++){
525                pmt_id = tof->GetPMTid(kk,hh);
526                Int_t Iplane=-1;
527                Int_t Ipaddle=-1;
528                //      Int_t IpaddleT=-1;
529                tof->GetPMTPaddle(pmt_id, Iplane, Ipaddle);
530                tof->GetPaddleGeometry(Iplane,Ipaddle,xleft,xright,yleft,yright);
531                if ( tofdedx->GetdEdx_pmt(pmt_id)>-1. &&((xtof_temp[Iplane]>=xleft&&xtof_temp[Iplane]<=xright) || (ytof_temp[Iplane]>=yleft&&ytof_temp[Iplane]<=yright)) ){  //attenzione:qui va inserito un controllo sulla traccia tof o sulle variabili di posizione !!!!
532                  //t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*4./dedx_corr[pmt_id]),t_tof->npmtadc);
533                  //t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*dedx_corr[pmt_id]/4.),t_tof->npmtadc);//annullo wolfrizzazione
534                  t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)),t_tof->npmtadc);//annullo wolfrizzazione
535                  t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
536                  t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
537                  t_tof->npmtadc++;
538                };
539              };
540            };
541          };
542    
543    
544    
545    
546    
547    
548    
549    
550    
551        
552         //       for (Int_t hh=0; hh<12;hh++){
553         //         for (Int_t kk=0; kk<4;kk++){
554         //           if ( tofoutput_.adctof_c[hh][kk] < 1000 ){
555         //             //      t_tof->dedx.AddAt(tofoutput_.adctof_c[hh][kk],t_tof->npmtadc); // EMILIANO
556         //             pmt_id = tof->GetPMTid(kk,hh);
557         //             t_tof->dedx.AddAt((tofoutput_.adctof_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO
558         //             t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
559         //             t_tof->adcflag.AddAt(tofoutput_.adcflagtof[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
560         //             t_tof->npmtadc++;
561         //           };
562         //         };
563         //       };
564         //
565        
566         //
567         new(t[ntrkentry]) ToFTrkVar(*t_tof);
568         ntrkentry++;      
569         t_tof->Clear();
570         //
571         //
572         //
573         t_pmt->Clear();
574         //
575        for (Int_t hh=0; hh<12;hh++){        for (Int_t hh=0; hh<12;hh++){
576          for (Int_t kk=0; kk<4;kk++){          for (Int_t kk=0; kk<4;kk++){
577           // new WM            // new WM
578  //        if ( tofoutput_.tdc_c[hh][kk] < 4095 || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){                //      if ( tofoutput_.tdc_c[hh][kk] < 4095 || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){    
579            if ( tdcc[kk][hh] < 4095. || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){                    if ( tdcc[kk][hh] < 4095. || adc[kk][hh] < 4095  || tdc[kk][hh] < 4095 ){        
580              //              //
581              t_pmt->pmt_id = tof->GetPMTid(kk,hh);              t_pmt->pmt_id = tof->GetPMTid(kk,hh);
# Line 431  cout<<"abs time "<<tabs<<" limit low "<< Line 608  cout<<"abs time "<<tabs<<" limit low "<<
608            for (Int_t e = 0; e < 5 ; e++){            for (Int_t e = 0; e < 5 ; e++){
609              tofinput_.al_pp[e] = ptt->al[e];              tofinput_.al_pp[e] = ptt->al[e];
610            };                  };      
611    
612              // new input for 9th reduction: tracker dEdx
613              tofinput_.trkmip = ptt->GetDEDX();
614    
615            //            //
616            // Get tracker related variables for this track            // Get tracker related variables for this track
617            //            //
618            toftrk();            toftrk();
619  //        toftrk(thelp);            //      toftrk(thelp);
620            //            //
621            // Copy values in the class from the structure (we need to use a temporary class to store variables).            // Copy values in the class from the structure (we need to use a temporary class to store variables).
622            //            //
# Line 454  cout<<"abs time "<<tabs<<" limit low "<< Line 635  cout<<"abs time "<<tabs<<" limit low "<<
635              t_tof->beta[kk] = tofoutput_.beta_a[kk];                  t_tof->beta[kk] = tofoutput_.beta_a[kk];    
636            };            };
637            //            //
638              memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
639              memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
640              memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
641              memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
642              //
643              tofdedx->Process(pam_event->GetOrbitalInfo()->absTime,t_tof->beta[12], (Float_t *)t_tof->xtr_tof,(Float_t *)t_tof->ytr_tof);
644            t_tof->npmtadc = 0;            t_tof->npmtadc = 0;
645    
646            for (Int_t hh=0; hh<12;hh++){            for (Int_t hh=0; hh<12;hh++){
647              for (Int_t kk=0; kk<4;kk++){              for (Int_t kk=0; kk<4;kk++){
648                if ( tofoutput_.adc_c[hh][kk] < 1000 ){                pmt_id = tof->GetPMTid(kk,hh);
649               //     t_tof->dedx.AddAt(tofoutput_.adc_c[hh][kk],t_tof->npmtadc); // EMILIANO                Int_t Iplane=-1;
650                  pmt_id = tof->GetPMTid(kk,hh);                Int_t Ipaddle=-1;
651                  t_tof->dedx.AddAt((tofoutput_.adc_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO                Int_t IpaddleT=-1;
652                  tof->GetPMTPaddle(pmt_id, Iplane, Ipaddle);
653                  IpaddleT=tof->GetPaddleIdOfTrack(t_tof->xtr_tof[Iplane],t_tof->ytr_tof[Iplane], Iplane,0.0);
654                  if ( tofdedx->GetdEdx_pmt(pmt_id) > -1. && Ipaddle==IpaddleT ){
655                    //t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*4./dedx_corr[pmt_id]),t_tof->npmtadc);
656                    //t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*dedx_corr[pmt_id]/4.),t_tof->npmtadc);//annullo wolfrizzazione
657                    t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)),t_tof->npmtadc);//annullo wolfrizzazione
658                  t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);                  t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
659                  t_tof->adcflag.AddAt(tofoutput_.adcflag[hh][kk],t_tof->npmtadc); // gf: Jan 09/07                  t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
660                  t_tof->npmtadc++;                  t_tof->npmtadc++;
661                };                };
662                  
663              };              };
664            };            };
665    
666    
667    //        for (Int_t hh=0; hh<12;hh++){
668    //          for (Int_t kk=0; kk<4;kk++){
669    //            pmt_id = tof->GetPMTid(kk,hh);
670    //            if ( tofdedx->GetdEdx_pmt(pmt_id) > -1. ){
671    //              t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*36./pow(dedx_corr[pmt_id],2)),t_tof->npmtadc);
672    //              t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
673    //              t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
674    //              printf(" nt %i npmtadc %i dedx %f dedx corr %f\n",nt,t_tof->npmtadc,(tofdedx->GetdEdx_pmt(pmt_id)*36./pow(dedx_corr[pmt_id],2)),dedx_corr[pmt_id]);
675    //              t_tof->npmtadc++;
676    //            };
677                  
678    //          };
679    //        };
680    //        t_tof->npmtadc = 0;
681    //        for (Int_t hh=0; hh<12;hh++){
682    //          for (Int_t kk=0; kk<4;kk++){
683    //            if ( tofoutput_.adc_c[hh][kk] < 1000 ){
684    //              //          t_tof->dedx.AddAt(tofoutput_.adc_c[hh][kk],t_tof->npmtadc); // EMILIANO
685    //              pmt_id = tof->GetPMTid(kk,hh);
686    //              t_tof->dedx.AddAt((tofoutput_.adc_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO
687    //              t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
688    //              t_tof->adcflag.AddAt(tofoutput_.adcflag[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
689    //              t_tof->npmtadc++;
690    //            };
691    //          };
692    //        };
693            //            //
694            memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));  
           memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));  
           memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));  
           memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));  
695            //            //
696            // Store the tracker track number in order to be sure to have shyncronized data during analysis            // Store the tracker track number in order to be sure to have shyncronized data during analysis
697            //            //
# Line 526  cout<<"abs time "<<tabs<<" limit low "<< Line 746  cout<<"abs time "<<tabs<<" limit low "<<
746    
747  void usage(){  void usage(){
748    
749      cout << "------------------------------------------------------------"<<endl;    cout << "------------------------------------------------------------"<<endl;
750      cout << "Loop over events (on one or more Level2-files), applying some selection cuts (defined in My-Selection.cpp), \n";    cout << "Loop over events (on one or more Level2-files), applying some selection cuts (defined in My-Selection.cpp), \n";
751      cout << "creating output histograms (defined in My-Histos.cpp) and/or trees with selected events. \n \n ";    cout << "creating output histograms (defined in My-Histos.cpp) and/or trees with selected events. \n \n ";
752      cout << "USAGE:"<<endl;    cout << "USAGE:"<<endl;
753      cout << "-processDir  DIR     -  Level2 data directory \n";    cout << "-processDir  DIR     -  Level2 data directory \n";
754      cout << "-processList LIST    -  list  of files (.txt) or single file (.root) to be analysed \n";    cout << "-processList LIST    -  list  of files (.txt) or single file (.root) to be analysed \n";
755      cout << "-outputFile  PATH    -  name of the output file \n";    cout << "-outputFile  PATH    -  name of the output file \n";
756      cout << "-NumEvents   XXX     -  number of events to be analysed \n";    cout << "-calibFile PATH+NAME -  name of the calibration file \n";
757      cout << "--debug, -g          -  debug mode \n";    cout << "-NumEvents   XXX     -  number of events to be analysed \n";
758      cout << "--help, -h           -  print this help \n";    cout << "--debug, -g          -  debug mode \n";
759      cout << "-options [ options ] -  options: \n";    cout << "--help, -h           -  print this help \n";
760      cout << "                        fillHistos --> create an output file with histograms  \n";    cout << "-options [ options ] -  options: \n";
761      cout << "                        fillTree   --> create an output file with trees storing the selected events \n ";    cout << "                        fillHistos --> create an output file with histograms  \n";
762      cout << "                        +(-)ALL    --> inlcude(exclude) all trees and branches \n "   ;    cout << "                        fillTree   --> create an output file with trees storing the selected events \n ";
763      cout << "                        +(-)TRK1 +(-)TRK2 +(-)CAL1 +(-)CAL2 +(-)TOF +(-)TRG +(-)ND +(-)S4 +(-)ORB --> inlcude(exclude) trees and branches  \n"  ;    cout << "                        +(-)ALL    --> inlcude(exclude) all trees and branches \n "   ;
764      cout << "------------------------------------------------------------"<<endl;    cout << "                        +(-)TRK1 +(-)TRK2 +(-)CAL1 +(-)CAL2 +(-)TOF +(-)TRG +(-)ND +(-)S4 +(-)ORB --> inlcude(exclude) trees and branches  \n"  ;
765      cout << "------------------------------------------------------------"<<endl;
766  }  }
767  //  //
768  int HandleInputPar(int argc, char **argv){  int HandleInputPar(int argc, char **argv){
769    
770      if(argc>1){    if(argc>1){
771                                    
772          if(!strcmp(argv[1], "-h") || !strcmp(argv[1], "--help") ){      if(!strcmp(argv[1], "-h") || !strcmp(argv[1], "--help") ){
773              usage();        usage();
774              return(1);        return(1);
775          };      };
776  // -----------------------      // -----------------------
777  // Read input parameters      // Read input parameters
778  // -----------------------        // -----------------------  
779          DIR     = gSystem->WorkingDirectory();      DIR     = gSystem->WorkingDirectory();
780          LIST    = "";      LIST    = "";
781          OUTF    = "myfile";      OUTF    = "myfile";
782          OPTIONS = "+AUTO -TOF";      OPTIONS = "+AUTO -TOF";
783                                    
784          for (int i = 1; i < argc; i++){              for (int i = 1; i < argc; i++){            
785              // -----------------------------------------------------//        // -----------------------------------------------------//
786              if (!strcmp(argv[i], "-processDir")){        if (!strcmp(argv[i], "-processDir")){
787                  if (++i >= argc) throw -1;          if (++i >= argc) throw -1;
788                  DIR = argv[i];          DIR = argv[i];
789                  cout << "processDir "<<DIR<<endl;          cout << "processDir "<<DIR<<endl;
790                  continue;          continue;
791              }        }
792              // -----------------------------------------------------//        // -----------------------------------------------------//
793              else if (!strcmp(argv[i], "-processList")){        else if (!strcmp(argv[i], "-processList")){
794                  if (++i >= argc) throw -1;          if (++i >= argc) throw -1;
795                  LIST = argv[i];          LIST = argv[i];
796                  cout << "processList "<<LIST<<endl;          cout << "processList "<<LIST<<endl;
797                  continue;          continue;
798              }          }  
799              // -----------------------------------------------------//        // -----------------------------------------------------//
800              else if (!strcmp(argv[i], "-outputFile")){        else if (!strcmp(argv[i], "-outputFile")){
801                  if (++i >= argc) throw -1;          if (++i >= argc) throw -1;
802                  OUTF = argv[i];          OUTF = argv[i];
803                  cout << "outputFile "<<OUTF<<endl;          cout << "outputFile "<<OUTF<<endl;
804                  continue;          continue;
805              }          }  
806              // -----------------------------------------------------//        // -----------------------------------------------------//
807              else if (!strcmp(argv[i], "-options")){        // -----------------------------------------------------//
808                  if (++i >= argc) throw -1;        else if (!strcmp(argv[i], "-calibFile")){
809                  OPTIONS = argv[i];          if (++i >= argc) throw -1;
810                  if( OPTIONS.Contains("[") ){          CALIBF = gSystem->ExpandPathName(argv[i]);
811                      do{          cout << "calibFile "<<CALIBF<<endl;
812                          if (++i >= argc) throw -1;          continue;
813                          OPTIONS.Append(argv[i]);        }  
814                      }while(!OPTIONS.Contains("]"));        // -----------------------------------------------------//
815                  }else cout << "wrong option format --> ignoring " << endl;        else if (!strcmp(argv[i], "-options")){
816              }          if (++i >= argc) throw -1;
817              else{          OPTIONS = argv[i];
818                  cout << "Unidentified input parameter. Ignored."<< endl;          if( OPTIONS.Contains("[") ){
819              };            do{
820          };              if (++i >= argc) throw -1;
821      }else{              OPTIONS.Append(argv[i]);
822          usage();            }while(!OPTIONS.Contains("]"));
823          return(1);          }else cout << "wrong option format --> ignoring " << endl;
824          }
825          else{
826            cout << "Unidentified input parameter. Ignored."<< endl;
827          };
828      };      };
829  // -----------------------    }else{
830  // Check input parameters      usage();
831  // -----------------------        return(1);
832      };
833      // -----------------------
834      // Check input parameters
835      // -----------------------  
836                    
837                                    
838      return(0);    return(0);
839                    
840  };  };
841  //  //
# Line 614  int HandleInputPar(int argc, char **argv Line 843  int HandleInputPar(int argc, char **argv
843  int main(int argc, char **argv)  int main(int argc, char **argv)
844  {  {
845                    
846      if( HandleInputPar(argc,argv) )return(1);    if( HandleInputPar(argc,argv) )return(1);
847                    
848      cout << "OPTIONS "<<OPTIONS<<endl;    cout << "OPTIONS "<<OPTIONS<<endl;
849      Loop(DIR,LIST,OPTIONS);    Loop(DIR,LIST,OPTIONS);
850            
851      cout << "Back to main - end"<<endl;    cout << "Back to main - end"<<endl;
852            
853      return 0;    return 0;
854                    
855  };  };
856    
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.7
  ViewVC Help
Powered by ViewVC 1.1.23