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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (show annotations) (download)
Mon Nov 23 09:38:05 2009 UTC (15 years ago) by mocchiut
Branch: MAIN
Changes since 1.2: +129 -25 lines
New dE/dx calibration implemented NB: it requires now -calibFile option to run

1 #if !defined(__CINT__) || defined(__MAKECINT__)
2
3
4 #include <fstream>
5 #include <sstream>
6 #include <string>
7 #include <math.h>
8 #include <iostream>
9 #include <TTree.h>
10 #include <TClassEdit.h>
11 #include <TObject.h>
12 #include <TList.h>
13 #include <TSystem.h>
14 #include <TSystemDirectory.h>
15 #include <TString.h>
16 #include <TFile.h>
17 #include <TClass.h>
18 #include <TCanvas.h>
19 #include <TH1.h>
20 #include <TH1F.h>
21 #include <TH2D.h>
22 #include <TLatex.h>
23 #include <TPad.h>
24 #include <TPaveLabel.h>
25 #include <TChain.h>
26 #include <TGraph.h>
27 #include <TMultiGraph.h>
28 #include <TLine.h>
29 #include <TStyle.h>
30 #include <TF1.h>
31 #include <TH1F.h>
32 #include <TTree.h>
33 #include <TCanvas.h>
34 #include <PamLevel2.h>
35 #include <iomanip>
36 //
37 #include <ToFCore.h>
38 //
39 // Declaration of the core fortran routines
40 //
41 #define tofl2com tofl2com_
42 extern "C" int tofl2com();
43 #define toftrk toftrk_
44 extern "C" int toftrk();
45 ///extern "C" int toftrk(float);
46 #define rdtofcal rdtofcal_
47 extern "C" int rdtofcal(char [], int *);
48
49 using namespace std;
50
51 #endif
52
53 //===================================
54 // global variables
55 //===================================
56 TString DIR;
57 TString OUTF;
58 TString LIST;
59 TString OPTIONS;
60 TString CALIBF;
61
62 PamLevel2 *pam_event = NULL;
63
64 //==========================================
65 //000000000000000000000000000000000000000000
66 //==========================================
67 Int_t Loop(TString ddir,TString list, TString options){
68 //
69 //
70 extern struct ToFInput tofinput_;
71 extern struct ToFOutput tofoutput_;
72 // --------------------
73 // read input file/list
74 // --------------------
75 pam_event = new PamLevel2(ddir,list,options);
76 pam_event->SetSELLI(2);
77 TTree::SetMaxTreeSize(1000*Long64_t(2000000000));
78 //
79 TString outfile_t = "";
80 outfile_t = OUTF;
81 outfile_t.Append(".root");
82 TFile *outft = (TFile*)gROOT->FindObject(outfile_t);
83 if (outft) outft->Close();
84 outft = new TFile(outfile_t,"RECREATE");
85 if(outft->IsZombie()){
86 cout << "Output file could not be created\n";
87 return 1;
88 };
89 cout << "Created output file: "<<outft->GetName()<<endl;
90 //
91 outft->cd();
92 ToFLevel2 *tof = new ToFLevel2();
93 ToFdEdx *tofdedx = new ToFdEdx();
94 TTree *toft = new TTree("ToF","PAMELA Level2 ToF data");
95 toft->SetAutoSave(900000000000000LL);
96 tof->Set();//ELENA **TEMPORANEO?**
97 toft->Branch("ToFLevel2","ToFLevel2",&tof);
98 //
99 pam_event->CreateCloneTrees(outft);
100 //
101 Int_t ntrkentry = 0;
102 Int_t npmtentry = 0;
103 //
104 ULong64_t nevents = pam_event->GetEntries();
105 printf("\n\n Running on %llu events \n\n",nevents);
106 //
107 TFile *tfile = TFile::Open(list.Data(),"READ");
108 TTree *toftr = (TTree*)tfile->Get("ToF");
109 ToFLevel2 *tofl2 = new ToFLevel2();
110 toftr->SetBranchAddress("ToFLevel2", &tofl2);
111 ULong64_t ntofev = toftr->GetEntries();
112 printf(" ToF events are %llu \n",ntofev);
113 //
114 //
115 GL_PARAM *glparam = new GL_PARAM();
116 TrkLevel2 *trk = new TrkLevel2();
117 //
118 TString host;
119 TString user;
120 TString psw;
121 const char *pamdbhost=gSystem->Getenv("PAM_DBHOST");
122 const char *pamdbuser=gSystem->Getenv("PAM_DBUSER");
123 const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW");
124 if ( !pamdbhost ) pamdbhost = "";
125 if ( !pamdbuser ) pamdbuser = "";
126 if ( !pamdbpsw ) pamdbpsw = "";
127 if ( strcmp(pamdbhost,"") ) host = pamdbhost;
128 if ( strcmp(pamdbuser,"") ) user = pamdbuser;
129 if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw;
130 //
131 //
132 TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
133 if ( !dbc->IsConnected() ) return 1;
134 stringstream myquery;
135 myquery.str("");
136 myquery << "SET time_zone='+0:00'";
137 dbc->Query(myquery.str().c_str());
138 glparam->Query_GL_PARAM(1,1,dbc); // parameters stored in DB in GL_PRAM table
139 trk->LoadField(glparam->PATH+glparam->NAME);
140 //
141 Long64_t myrun = 0;
142 Long64_t mysrun = -1;
143 Int_t adc[4][12];
144 Int_t tdc[4][12];
145 Float_t tdcc[4][12];
146 //
147 Bool_t defcal = true;
148 for(ULong64_t iev=0; iev<nevents; iev++){
149
150
151 //==================================================================
152
153
154 // for(ULong64_t iev=0; iev<50; iev++){
155 //
156 if ( !(iev%1000) ) printf(" %iK \n",(int)iev/1000);
157 pam_event->Clear();
158 //
159 if( pam_event->GetEntry(iev) ){
160 //
161 pam_event->FillCloneTrees();
162 toftr->GetEntry(iev);
163 //
164 // variable to save information about the tof calibration used
165 //
166 //
167 myrun=pam_event->GetRunInfo()->ID;
168 if ( myrun != mysrun ){
169 // printf(" rhtime %u myrun %i mysrun %i \n",pam_event->GetRunInfo()->RUNHEADER_TIME,(int)myrun,(int)mysrun);
170 mysrun = myrun;
171 //
172
173
174
175
176 Int_t error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,204,dbc); // parameters stored in DB in GL_PRAM table
177 if ( error<0 ) {
178 return(1);
179 };
180 //
181 tofdedx->ReadParAtt((glparam->PATH+glparam->NAME).Data());
182
183 //
184 error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,205,dbc); // parameters stored in DB in GL_PRAM table
185 if ( error<0 ) {
186 return(1);
187 };
188 //
189 tofdedx->ReadParPos((glparam->PATH+glparam->NAME).Data());
190
191 //
192 error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,206,dbc); // parameters stored in DB in GL_PRAM table
193 if ( error<0 ) {
194 return(1);
195 };
196 //
197 tofdedx->ReadParBBneg((glparam->PATH+glparam->NAME).Data());
198
199 //
200 error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,207,dbc); // parameters stored in DB in GL_PRAM table
201 if ( error<0 ) {
202 return(1);
203 };
204 //
205 tofdedx->ReadParBBpos((glparam->PATH+glparam->NAME).Data());
206
207 //
208 error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,208,dbc); // parameters stored in DB in GL_PRAM table
209 if ( error<0 ) {
210 return(1);
211 };
212 //
213 tofdedx->ReadParDesatBB((glparam->PATH+glparam->NAME).Data());
214
215
216 tofdedx->CheckConnectors(pam_event->GetRunInfo()->RUNHEADER_TIME,glparam,dbc);
217
218
219
220 error=glparam->Query_GL_PARAM(pam_event->GetRunInfo()->RUNHEADER_TIME,201,dbc); // parameters stored in DB in GL_PRAM table
221 if ( error<0 ) {
222 return(1);
223 };
224 //
225 printf(" Reading ToF parameter file: %s \n",(glparam->PATH+glparam->NAME).Data());
226 //
227 if ( (UInt_t)glparam->TO_TIME != (UInt_t)4294967295UL ) defcal = false;
228 //
229 Int_t nlen = (Int_t)(glparam->PATH+glparam->NAME).Length();
230 rdtofcal((char *)(glparam->PATH+glparam->NAME).Data(),&nlen);
231 //
232 };
233 //
234
235
236 //================================================================
237 //==================================================================
238 //--- Define absolute time
239 UInt_t tabs=pam_event->GetOrbitalInfo()->absTime;
240
241 printf(" ATIME %u re %u \n",(Int_t)tabs,(UInt_t)iev);
242
243 //==================================================================
244 Float_t dedx_corr_m[2000][48],dedx_corr[48];
245 Double_t mtime[2000],t1,t2,tm;
246 Float_t yhelp1,yhelp2,slope,inter,thelp1,thelp2;
247 Float_t xmean1,xwidth1;
248
249 Int_t ical,ii,j,jj;
250
251 if (iev==0) {
252
253 ical=0; // counter set to zero if first-time reading
254
255 //-----------------------------------------------------------
256 // Here I read the dEdx_korr parameters
257 //-----------------------------------------------------------
258
259 jj=0;
260
261 ifstream fin(CALIBF.Data());
262
263 while (! fin.eof()) {
264 fin>>t1>>tm>>t2;
265 // cout<<jj<<" "<<tm<<endl;
266 cout << setiosflags(ios::fixed) << setw(10) << setprecision(3) << tm << endl;
267 mtime[jj]=tm;
268 for (ii=0; ii<48;ii++) {
269 fin>>j>>xmean1>>xwidth1;
270 dedx_corr_m[jj][ii]=xmean1;
271 }
272 jj=jj+1;
273 // printf(" kk %i \n",jj);
274 }
275 // printf(" 1ical %i \n",ical);
276
277 fin.close();
278
279
280 while (tabs<mtime[ical] || tabs > mtime[ical+1]) {
281 // printf(" ical %i \n",ical);
282 ical = ical+1;
283 }
284 // ical = ical-1;
285 cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;
286 } // if iev==0...
287
288 //==================================================================
289 //== End of first time reading & filling the array
290 //==================================================================
291
292 //==================================================================
293 //== if time is outside time limits:
294 //==================================================================
295
296 if (tabs<mtime[ical] || tabs>mtime[ical+1]) {
297 cout<<"Checking Time Limits!"<<endl;
298 ical=0;
299 while (tabs > mtime[ical+1] || tabs<mtime[ical]) {
300 ical = ical+1;
301 }
302 // ical = ical-1;
303 cout<<"abs time "<<tabs<<" limit low "<<mtime[ical]<<" limit up "<<mtime[ical+1]<<" ical "<<ical<<endl;
304 };
305
306 // printf(" 2ical %i \n",ical);
307 //==================================================================
308 //== interpolate betwen time limits
309 //==================================================================
310
311 thelp1 = mtime[ical];
312 thelp2 = mtime[ical+1];
313
314 for (ii=0; ii<48;ii++) {
315 yhelp1 = dedx_corr_m[ical][ii];
316 // yhelp1 = 6.;
317 if ( yhelp1 < 0.1 ) yhelp1 = 6.;
318 yhelp2 = dedx_corr_m[ical+1][ii];
319 // yhelp2 = 6.;
320 if ( yhelp2 < 0.1 ) yhelp2 = 6.;
321 slope = (yhelp2-yhelp1)/(thelp2-thelp1);
322 inter = yhelp1 - slope*thelp1;
323 dedx_corr[ii] = slope*tabs + inter;
324 // if (ii==0) cout<<thelp1<<" "<<thelp2<<" "<<tabs<<" "<<yhelp1<<" "<<yhelp2<<" "<<dedx_corr[0]<<endl;
325 }
326
327 //================================================================
328 //================================================================
329 // printf("cpippo \n");
330
331 // process tof data
332 //
333 for (Int_t hh=0; hh<12;hh++){
334 for (Int_t kk=0; kk<4;kk++){
335 adc[kk][hh] = 4095;
336 tdc[kk][hh] = 4095;
337 tdcc[kk][hh] = 4095.;
338 tofinput_.adc[hh][kk] = 4095;
339 tofinput_.tdc[hh][kk] = 4095;
340 };
341 };
342 // memset(adc, 0, 12*4*sizeof(Int_t));
343 // memset(tdc, 0, 12*4*sizeof(Int_t));
344 Int_t gg = 0;
345 Int_t hh = 0;
346 Int_t adcf[48];
347 memset(adcf, 0, 48*sizeof(Int_t));
348 Int_t tdcf[48];
349 memset(tdcf, 0, 48*sizeof(Int_t));
350 for (Int_t pm=0; pm < tofl2->ntrk() ; pm++){
351 ToFTrkVar *ttf = tofl2->GetToFTrkVar(pm);
352 for ( Int_t nc=0; nc < ttf->npmttdc; nc++){
353 if ( (ttf->tdcflag).At(nc) != 0 ) tdcf[(ttf->pmttdc).At(nc)] = 1;
354 };
355 for ( Int_t nc=0; nc < ttf->npmtadc; nc++){
356 if ( (ttf->adcflag).At(nc) != 0 ) adcf[(ttf->pmtadc).At(nc)] = 1;
357 };
358 };
359 //
360 for (Int_t pm=0; pm < tofl2->npmt() ; pm++){
361 ToFPMT *pmt = tofl2->GetToFPMT(pm);
362 tofl2->GetPMTIndex(pmt->pmt_id, gg, hh);
363 if ( adcf[pmt->pmt_id] == 0 ){
364 tofinput_.adc[gg][hh] = (int)pmt->adc;
365 adc[hh][gg] = (int)pmt->adc;
366 };
367 if ( tdcf[pmt->pmt_id] == 0 ){
368 tofinput_.tdc[gg][hh] = (int)pmt->tdc;
369 tdc[hh][gg] = (int)pmt->tdc;
370 };
371 tdcc[hh][gg] = (float)pmt->tdc_tw;
372 // Int_t pppid = tof->GetPMTid(hh,gg);
373 // 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);
374 };
375 for (Int_t hh=0; hh<12;hh++){
376 for (Int_t kk=0; kk<4;kk++){
377 tofdedx->Init(kk,hh,adc[kk][hh]);
378 };
379 };
380 // for (Int_t pm=0; pm <48 ; pm++){
381 // tof->GetPMTIndex(pm, gg, hh);
382 // tofinput_.tdc[hh][gg] = (int)500.;
383 // tofinput_.adc[hh][gg] = (int)500.;
384 // tdc[hh][gg] = (int)500.;
385 // adc[hh][gg] = (int)500.;
386 // // printf(" hh %i gg %i tdc %f adc %f \n",hh,gg,pmt->tdc,pmt->adc);
387 // };
388 //
389 for (Int_t hh=0; hh<5;hh++){
390 tofinput_.patterntrig[hh]=pam_event->GetTrigLevel2()->patterntrig[hh];
391 };
392 //
393 tof->Clear();
394 Int_t pmt_id = 0;
395 ToFPMT *t_pmt = new ToFPMT();
396 if(!(tof->PMT))tof->PMT = new TClonesArray("ToFPMT",12); //ELENA
397 TClonesArray &tpmt = *tof->PMT;
398 ToFTrkVar *t_tof = new ToFTrkVar();
399 if(!(tof->ToFTrk))tof->ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA
400 TClonesArray &t = *tof->ToFTrk;
401 //
402 //
403 // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract common and track-related variables.
404 //
405 npmtentry = 0;
406 //
407 ntrkentry = 0;
408 //
409 // Calculate tracks informations from ToF alone
410 //
411 tofl2com();
412 //
413 memcpy(tof->tof_j_flag,tofoutput_.tof_j_flag,6*sizeof(Int_t));
414 //
415 t_tof->trkseqno = -1;
416 //
417 // and now we must copy from the output structure to the level2 class:
418 //
419 t_tof->npmttdc = 0;
420 //
421 for (Int_t hh=0; hh<12;hh++){
422 for (Int_t kk=0; kk<4;kk++){
423 if ( tofoutput_.tofmask[hh][kk] != 0 ){
424 pmt_id = tof->GetPMTid(kk,hh);
425 t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc);
426 t_tof->tdcflag.AddAt(tofoutput_.tdcflagtof[hh][kk],t_tof->npmttdc); // gf: Jan 09/07
427 t_tof->npmttdc++;
428 };
429 };
430 };
431 for (Int_t kk=0; kk<13;kk++){
432 t_tof->beta[kk] = tofoutput_.betatof_a[kk];
433 }
434 //
435 memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
436 memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
437 memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
438 memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
439 {
440 Float_t xtof_temp[6]={0.,t_tof->xtofpos[0],t_tof->xtofpos[1],0.,0.,t_tof->xtofpos[2]};
441 Float_t ytof_temp[6]={t_tof->ytofpos[0],0.,0.,t_tof->ytofpos[1],t_tof->ytofpos[2],0.};
442 tofdedx->Process(pam_event->GetOrbitalInfo()->absTime,t_tof->beta[12], (Float_t *)xtof_temp,(Float_t *)ytof_temp);
443 }
444 t_tof->npmtadc = 0;
445
446
447 for (Int_t hh=0; hh<12;hh++){
448 for (Int_t kk=0; kk<4;kk++){
449 pmt_id = tof->GetPMTid(kk,hh);
450 if ( tofdedx->GetdEdx_pmt(pmt_id)>-1. ){
451 t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*36./pow(dedx_corr[pmt_id],2)),t_tof->npmtadc);
452 t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
453 t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
454 t_tof->npmtadc++;
455 };
456 };
457 };
458 // for (Int_t hh=0; hh<12;hh++){
459 // for (Int_t kk=0; kk<4;kk++){
460 // if ( tofoutput_.adctof_c[hh][kk] < 1000 ){
461 // // t_tof->dedx.AddAt(tofoutput_.adctof_c[hh][kk],t_tof->npmtadc); // EMILIANO
462 // pmt_id = tof->GetPMTid(kk,hh);
463 // t_tof->dedx.AddAt((tofoutput_.adctof_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO
464 // t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
465 // t_tof->adcflag.AddAt(tofoutput_.adcflagtof[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
466 // t_tof->npmtadc++;
467 // };
468 // };
469 // };
470 //
471
472 //
473 new(t[ntrkentry]) ToFTrkVar(*t_tof);
474 ntrkentry++;
475 t_tof->Clear();
476 //
477 //
478 //
479 t_pmt->Clear();
480 //
481 for (Int_t hh=0; hh<12;hh++){
482 for (Int_t kk=0; kk<4;kk++){
483 // new WM
484 // if ( tofoutput_.tdc_c[hh][kk] < 4095 || adc[kk][hh] < 4095 || tdc[kk][hh] < 4095 ){
485 if ( tdcc[kk][hh] < 4095. || adc[kk][hh] < 4095 || tdc[kk][hh] < 4095 ){
486 //
487 t_pmt->pmt_id = tof->GetPMTid(kk,hh);
488 t_pmt->tdc_tw = tofoutput_.tdc_c[hh][kk];
489 t_pmt->adc = (Float_t)adc[kk][hh];
490 t_pmt->tdc = (Float_t)tdc[kk][hh];
491 //
492 new(tpmt[npmtentry]) ToFPMT(*t_pmt);
493 npmtentry++;
494 t_pmt->Clear();
495 };
496 };
497 };
498 //
499 // Calculate track-related variables
500 //
501 if ( pam_event->GetTrkLevel2()->ntrk() > 0 ){
502 //
503 // We have at least one track
504 //
505 //
506 // Run over tracks
507 //
508 for(Int_t nt=0; nt < pam_event->GetTrkLevel2()->ntrk(); nt++){
509 //
510 TrkTrack *ptt = pam_event->GetTrkLevel2()->GetStoredTrack(nt);
511 //
512 // Copy the alpha vector in the input structure
513 //
514 for (Int_t e = 0; e < 5 ; e++){
515 tofinput_.al_pp[e] = ptt->al[e];
516 };
517 //
518 // Get tracker related variables for this track
519 //
520 toftrk();
521 // toftrk(thelp);
522 //
523 // Copy values in the class from the structure (we need to use a temporary class to store variables).
524 //
525 t_tof->npmttdc = 0;
526 for (Int_t hh=0; hh<12;hh++){
527 for (Int_t kk=0; kk<4;kk++){
528 if ( tofoutput_.tofmask[hh][kk] != 0 ){
529 pmt_id = tof->GetPMTid(kk,hh);
530 t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc);
531 t_tof->tdcflag.AddAt(tofoutput_.tdcflag[hh][kk],t_tof->npmttdc); // gf: Jan 09/07
532 t_tof->npmttdc++;
533 };
534 };
535 };
536 for (Int_t kk=0; kk<13;kk++){
537 t_tof->beta[kk] = tofoutput_.beta_a[kk];
538 };
539 //
540 memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos));
541 memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos));
542 memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof));
543 memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof));
544 //
545 tofdedx->Process(pam_event->GetOrbitalInfo()->absTime,t_tof->beta[12], (Float_t *)t_tof->xtr_tof,(Float_t *)t_tof->ytr_tof);
546 t_tof->npmtadc = 0;
547 for (Int_t hh=0; hh<12;hh++){
548 for (Int_t kk=0; kk<4;kk++){
549 pmt_id = tof->GetPMTid(kk,hh);
550 if ( tofdedx->GetdEdx_pmt(pmt_id) > -1. ){
551 t_tof->dedx.AddAt((tofdedx->GetdEdx_pmt(pmt_id)*36./pow(dedx_corr[pmt_id],2)),t_tof->npmtadc);
552 t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
553 t_tof->adcflag.AddAt(0,t_tof->npmtadc); // gf: Jan 09/07
554 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]);
555 t_tof->npmtadc++;
556 };
557
558 };
559 };
560 // t_tof->npmtadc = 0;
561 // for (Int_t hh=0; hh<12;hh++){
562 // for (Int_t kk=0; kk<4;kk++){
563 // if ( tofoutput_.adc_c[hh][kk] < 1000 ){
564 // // t_tof->dedx.AddAt(tofoutput_.adc_c[hh][kk],t_tof->npmtadc); // EMILIANO
565 // pmt_id = tof->GetPMTid(kk,hh);
566 // t_tof->dedx.AddAt((tofoutput_.adc_c[hh][kk]*4./dedx_corr[pmt_id]),t_tof->npmtadc); // EMILIANO
567 // t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc);
568 // t_tof->adcflag.AddAt(tofoutput_.adcflag[hh][kk],t_tof->npmtadc); // gf: Jan 09/07
569 // t_tof->npmtadc++;
570 // };
571 // };
572 // };
573 //
574
575 //
576 // Store the tracker track number in order to be sure to have shyncronized data during analysis
577 //
578 t_tof->trkseqno = nt;
579 //
580 // create a new object for this event with track-related variables
581 //
582 new(t[ntrkentry]) ToFTrkVar(*t_tof);
583 ntrkentry++;
584 t_tof->Clear();
585 //
586 }; // loop on all the tracks
587 };
588 //
589 tof->unpackError = tofl2->unpackError;
590 if ( defcal ){
591 tof->default_calib = 1;
592 } else {
593 tof->default_calib = 0;
594 };
595 //
596 // Fill the rootple
597 //
598 toft->Fill();
599 //
600 // toft->Show();
601 // toftr->Show();
602 //
603 delete t_tof;
604 //
605 //
606 };
607 //
608 };
609 //
610 pam_event->Clear();
611 //
612 outft->cd();
613 pam_event->WriteCloneTrees();
614 toftr->Delete();
615
616 toft->Write();
617 outft->Close();
618 //
619 return 0;
620 //
621 }
622
623 /////////////////////////////////////////////////////////////////
624
625 #if !defined(__CINT__)
626
627 void usage(){
628
629 cout << "------------------------------------------------------------"<<endl;
630 cout << "Loop over events (on one or more Level2-files), applying some selection cuts (defined in My-Selection.cpp), \n";
631 cout << "creating output histograms (defined in My-Histos.cpp) and/or trees with selected events. \n \n ";
632 cout << "USAGE:"<<endl;
633 cout << "-processDir DIR - Level2 data directory \n";
634 cout << "-processList LIST - list of files (.txt) or single file (.root) to be analysed \n";
635 cout << "-outputFile PATH - name of the output file \n";
636 cout << "-calibFile PATH+NAME - name of the calibration file \n";
637 cout << "-NumEvents XXX - number of events to be analysed \n";
638 cout << "--debug, -g - debug mode \n";
639 cout << "--help, -h - print this help \n";
640 cout << "-options [ options ] - options: \n";
641 cout << " fillHistos --> create an output file with histograms \n";
642 cout << " fillTree --> create an output file with trees storing the selected events \n ";
643 cout << " +(-)ALL --> inlcude(exclude) all trees and branches \n " ;
644 cout << " +(-)TRK1 +(-)TRK2 +(-)CAL1 +(-)CAL2 +(-)TOF +(-)TRG +(-)ND +(-)S4 +(-)ORB --> inlcude(exclude) trees and branches \n" ;
645 cout << "------------------------------------------------------------"<<endl;
646 }
647 //
648 int HandleInputPar(int argc, char **argv){
649
650 if(argc>1){
651
652 if(!strcmp(argv[1], "-h") || !strcmp(argv[1], "--help") ){
653 usage();
654 return(1);
655 };
656 // -----------------------
657 // Read input parameters
658 // -----------------------
659 DIR = gSystem->WorkingDirectory();
660 LIST = "";
661 OUTF = "myfile";
662 OPTIONS = "+AUTO -TOF";
663
664 for (int i = 1; i < argc; i++){
665 // -----------------------------------------------------//
666 if (!strcmp(argv[i], "-processDir")){
667 if (++i >= argc) throw -1;
668 DIR = argv[i];
669 cout << "processDir "<<DIR<<endl;
670 continue;
671 }
672 // -----------------------------------------------------//
673 else if (!strcmp(argv[i], "-processList")){
674 if (++i >= argc) throw -1;
675 LIST = argv[i];
676 cout << "processList "<<LIST<<endl;
677 continue;
678 }
679 // -----------------------------------------------------//
680 else if (!strcmp(argv[i], "-outputFile")){
681 if (++i >= argc) throw -1;
682 OUTF = argv[i];
683 cout << "outputFile "<<OUTF<<endl;
684 continue;
685 }
686 // -----------------------------------------------------//
687 // -----------------------------------------------------//
688 else if (!strcmp(argv[i], "-calibFile")){
689 if (++i >= argc) throw -1;
690 CALIBF = argv[i];
691 cout << "calibFile "<<CALIBF<<endl;
692 continue;
693 }
694 // -----------------------------------------------------//
695 else if (!strcmp(argv[i], "-options")){
696 if (++i >= argc) throw -1;
697 OPTIONS = argv[i];
698 if( OPTIONS.Contains("[") ){
699 do{
700 if (++i >= argc) throw -1;
701 OPTIONS.Append(argv[i]);
702 }while(!OPTIONS.Contains("]"));
703 }else cout << "wrong option format --> ignoring " << endl;
704 }
705 else{
706 cout << "Unidentified input parameter. Ignored."<< endl;
707 };
708 };
709 }else{
710 usage();
711 return(1);
712 };
713 // -----------------------
714 // Check input parameters
715 // -----------------------
716
717
718 return(0);
719
720 };
721 //
722
723 int main(int argc, char **argv)
724 {
725
726 if( HandleInputPar(argc,argv) )return(1);
727
728 cout << "OPTIONS "<<OPTIONS<<endl;
729 Loop(DIR,LIST,OPTIONS);
730
731 cout << "Back to main - end"<<endl;
732
733 return 0;
734
735 };
736
737 #endif

  ViewVC Help
Powered by ViewVC 1.1.23