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

Annotation of /DarthVader/CalorimeterLevel2/src/CaloCore.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.60 - (hide annotations) (download)
Tue Sep 16 19:47:32 2014 UTC (10 years, 2 months ago) by mocchiutti
Branch: MAIN
Changes since 1.59: +3 -0 lines
10RED bug fixed: crash when reprocessing without one of the extended tracks

1 mocchiut 1.53 /**
2     * \file CaloCore.cpp
3     * \author Emiliano Mocchiutti
4     */
5 mocchiut 1.1 //
6     // Given a calibration and a data file this program create an ntuple with LEVEL2 calorimeter variables - Emiliano Mocchiutti
7     //
8 mocchiut 1.3 // CaloCore.cxx version 3.05 (2006-05-30)
9 mocchiut 1.1 //
10     // The only input needed is the path to the directory created by YODA for the data file you want to analyze.
11     //
12     // Changelog:
13 mocchiut 1.3 //
14 mocchiut 1.53 // 7.00 (2013-11): Code to call this routine for level1 data only added.
15     //
16 mocchiut 1.11 // 3.08 (2006-11-13): Added high energy nuclei capability and "process all events" capability.
17     //
18 mocchiut 1.4 // 3.04 - 3.05 (2006-05-30): Qlast and nlast are now calculated using 4 (not 8) strips aournd the shower axis. Small bug fixed.
19 mocchiut 1.2 //
20     // 3.03 - 3.04 (2006-05-23): Forgot to put impx and impy in the PAMELA reference system, fixed.
21 mocchiut 1.1 //
22     // 3.02 - 3.03 (2006-05-18): updated to be called in DarthVader. Output dimension are now in cm and in the PAMELA reference system.
23     //
24     // 3.01 - 3.02 (2006-04-21): when copying entries get size of caclone and not of ca... (shouldn't matter). Fixed increasing file dimension bug when reprocessing.
25     // Added variable planemax[2], plane of maximum energy release (x and y) in final output. Use ItoRunInfo instead of RunInfo.
26     //
27     // 3.00 - 3.01 (2006-04-14): fixed small bug in tagging the track used to determine track-related variables, put in caloprocessing the opening of parameters files, fixed
28     // small bug in fortran routines
29     //
30     // 2.01 - 3.00 (2006-04-14): almost everything has changed. Now it can process one, all or some runs, introduced the final CaloLevel2 class+methods and the
31     // working class "CaloProcessing", linked to the preliminary tracker flight software v0r00, reads YODA unique output files,
32     // reduced the number of installed libraries, F77 programs splitted depending on the function contained, introduced the readout and
33     // processing of self-trigger events, if the tracker provides more than one track all calorimeter track-related variables are saved
34     // as many times as the number of tracks (via the TClonesArray object in the level2 rootple) and many other small changes.
35     //
36     // 2.00 - 2.01 (2006-01-26): bug: wrong calculation of baselines in some cases, fixed.
37     //
38     // 1.00 - 2.00 (2006-01-11): use TSQL ROOT classes instead of directly calling MySQL.
39     //
40     // 0.00 - 1.00 (2005-09-14): seems working.
41     //
42     // 0.00 (2005-09-09): clone of CaloLEVEL2.c .
43     //
44     // C/C++ headers
45     //
46     #include <fstream>
47     #include <string.h>
48     //
49     // ROOT headers
50     //
51     #include <TTree.h>
52     #include <TClassEdit.h>
53     #include <TObject.h>
54     #include <TList.h>
55 mocchiut 1.8 #include <TArrayI.h>
56 mocchiut 1.1 #include <TSystem.h>
57     #include <TSystemDirectory.h>
58     #include <TString.h>
59     #include <TFile.h>
60     #include <TClass.h>
61     #include <TCanvas.h>
62     #include <TH1.h>
63     #include <TH1F.h>
64     #include <TH2D.h>
65     #include <TLatex.h>
66     #include <TPad.h>
67     #include <TSQLServer.h>
68     #include <TSQLRow.h>
69     #include <TSQLResult.h>
70     #include <TClonesArray.h>
71     #include <TStreamerInfo.h>
72     //
73 mocchiut 1.20 // This program headers
74     //
75     #include <RunInfo.h>
76     //
77 mocchiut 1.1 // YODA headers
78     //
79     #include <PamelaRun.h>
80     #include <physics/trigger/TriggerEvent.h>
81     //
82     // This program headers
83     //
84     #include <CaloCore.h>
85 mocchiut 1.13 #include <CaloLevel1.h>
86 mocchiut 1.1 #include <CaloLevel2.h>
87 mocchiut 1.22 #include <CaloLevel0.h>
88 mocchiut 1.1 #include <CaloVerl2.h>
89     //
90     // Tracker classes headers and definitions
91     //
92     #include <TrkLevel2.h>
93 mocchiut 1.55 #include <ExtTrack.h> // new tracking code
94 mocchiut 1.1 //
95     using namespace std;
96     //
97 mocchiut 1.47 // IsCaloSelfTrigger method from Nicola (adapted)
98     //
99     //===============================================================================
100     // The IsCaloSelfTrigger functions checks in the trigger pattern register
101     // whether the event has been triggered by the calorimeter or not.
102     // Since the calorimeter self trigger coincidences are slower than both
103     // TOF and S4 triggers the triggering condition is checked for these detectors
104     // before looking for a calorimeter coincidence
105     //===============================================================================
106     //Int_t IsCaloSelfTrigger(TrigLevel2* event)
107     Int_t IsCaloSelfTrigger(UInt_t tc, UInt_t *pt)
108     {
109     // UInt_t tc, *pt;
110     bool S1or, S2or, S3or;
111     bool S1and, S2and, S3and;
112    
113     // tc = event->trigconf;
114     // pt = event->patterntrig;
115    
116     if (tc & CCALIB_ON) return false; // calibration trigger
117    
118     S1or = S11(pt) || S12(pt);
119     S2or = S21(pt) || S22(pt);
120     S3or = S31(pt) || S32(pt);
121     S1and = S11(pt) && S12(pt);
122     S2and = S21(pt) && S22(pt);
123     S3and = S31(pt) && S32(pt);
124    
125     switch (tc & ~(CCALO | CS4 | CCALIB_ON)) {
126     case CTOF1:
127     if (S1or && S2or && S3or) return 0; // TOF1 trigger
128     break;
129     case CTOF2:
130     if (S1and && S2and && S3and) return 0; // TOF2 trigger
131     break;
132     case CTOF3:
133     if (S2or && S3or) return 0; // TOF3 trigger
134     break;
135     case CTOF4:
136     if (S2and && S3and) return 0; // TOF4 trigger
137     break;
138     case CTOF5:
139     if (S12(pt) && S2and) return 0; // TOF5 trigger
140     break;
141     case CTOF6:
142     if (S1or && S3or) return 0; // TOF6 trigger
143     break;
144     case CTOF7:
145     if (S1and && S3and) return 0; // TOF7 trigger
146     break;
147     default: // other trigger configuration
148     break;
149     }
150    
151     if ((tc & CS4) && PS4(pt)) return 2; // S4 trigger
152    
153     if (PCaloST(pt)) return 2; // calorimeter trigger
154    
155     return 2; // nobody gernerated the trigger?
156     }
157    
158     //
159 mocchiut 1.1 // CORE ROUTINE
160     //
161 mocchiut 1.32 int CaloCore(UInt_t run, TFile *file, GL_TABLES *glt, Int_t calargc, char *calargv[]){
162 mocchiut 1.1 //
163     // Set these to true to have a very verbose output.
164     //
165     Bool_t verbose = false;
166     Bool_t debug = false;
167     //
168 mocchiut 1.21 Bool_t crosst = true;
169 mocchiut 1.36 Bool_t ctground = false;
170 mocchiut 1.41 Bool_t usetable = true;
171     Bool_t noselfct = false;
172 mocchiut 1.21 //
173 mocchiut 1.11 Bool_t trackanyway = true;
174     //
175     Float_t rigdefault = 50.;
176 mocchiut 1.35 //
177 mocchiut 1.11 Bool_t hZn = true;
178     //
179     Bool_t withtrk = true;
180     //
181     Bool_t st = true;
182     //
183 mocchiut 1.13 Bool_t getl1 = true;
184     //
185 mocchiut 1.23 Bool_t mechal = false;
186     //
187 mocchiut 1.30 Bool_t checksimu = true;
188     //
189 mocchiut 1.47 Bool_t mask18 = false;
190     //
191 mocchiut 1.53 Bool_t cl1only = false;
192     //
193     Bool_t froml1 = false;
194     //
195     Bool_t delct = false;
196     //
197 mocchiut 1.1 // Output directory is the working directoy.
198     //
199     const char* outdir = gSystem->DirName(gSystem->DirName(file->GetPath()));
200     //
201     Int_t ri = 0;
202 mocchiut 1.26 TString processFolder = Form("calorimeterFolder_%u",run);
203 mocchiut 1.1 if ( calargc > 0 ){
204     ri = 0;
205     while ( ri < calargc ){
206     if ( !strcmp(calargv[ri],"-processFolder") ) {
207     if ( calargc < ri+1 ){
208     throw -3;
209     };
210     processFolder = (TString)calargv[ri+1];
211     ri++;
212     };
213     if ( !strcmp(calargv[ri],"-v") || !strcmp(calargv[ri],"--verbose") ) {
214     verbose = true;
215     };
216     if ( !strcmp(calargv[ri],"-g") || !strcmp(calargv[ri],"--debug") ) {
217 mocchiut 1.13 verbose = true;
218 mocchiut 1.1 debug = true;
219     };
220 mocchiut 1.11 if ( !strcmp(calargv[ri],"--alltracks") ) {
221     trackanyway = true;
222     };
223 mocchiut 1.23 if ( !strcmp(calargv[ri],"--use-default-alig") ) {
224     mechal = true;
225     };
226 mocchiut 1.21 if ( !strcmp(calargv[ri],"--no-crosstalk") ) {
227     crosst = false;
228     };
229 mocchiut 1.47 if ( !strcmp(calargv[ri],"--mask-plane18X") ) {
230     mask18 = true;
231     };
232     if ( !strcmp(calargv[ri],"--unmask-plane18X") ) {
233     mask18 = false;
234     };
235 mocchiut 1.34 if ( !strcmp(calargv[ri],"--flight-crosstalk") ) {
236     ctground = false;
237     };
238 mocchiut 1.41 if ( !strcmp(calargv[ri],"--ct-use-pulse") ) {
239     usetable = false;
240     };
241     if ( !strcmp(calargv[ri],"--ct-use-table") ) {
242     usetable = true;
243     };
244 mocchiut 1.34 if ( !strcmp(calargv[ri],"--ground-crosstalk") ) {
245     ctground = true;
246     };
247 mocchiut 1.41 if ( !strcmp(calargv[ri],"--no-self-crosstalk") ) {
248     noselfct = true;
249     };
250 mocchiut 1.11 if ( !strcmp(calargv[ri],"--no-tracker") ) {
251     withtrk = false;
252     };
253     if ( !strcmp(calargv[ri],"--with-tracker") ) {
254     withtrk = true;
255     };
256     if ( !strcmp(calargv[ri],"--defrig") ) {
257     if ( calargc < ri+1 ){
258     throw -3;
259     };
260     rigdefault = atof(calargv[ri+1]);
261     ri++;
262     };
263     if ( !strcmp(calargv[ri],"--no-alltracks") ) {
264     trackanyway = false;
265     };
266     if ( !strcmp(calargv[ri],"--highZnuclei") ) {
267     hZn = true;
268     };
269     if ( !strcmp(calargv[ri],"--no-highZnuclei") ) {
270     hZn = false;
271     };
272     if ( !strcmp(calargv[ri],"--selftrigger") ) {
273     st = true;
274     };
275     if ( !strcmp(calargv[ri],"--no-selftrigger") ) {
276     st = false;
277     };
278 mocchiut 1.13 if ( !strcmp(calargv[ri],"--no-level1") ) {
279     getl1 = false;
280 mocchiut 1.53 cl1only = false;
281     };
282     if ( !strcmp(calargv[ri],"--delete-calo-tree") ) {
283     delct = true;
284 mocchiut 1.13 };
285 mocchiut 1.41 if ( !strcmp(calargv[ri],"--level1") ) {
286     getl1 = true;
287     };
288 mocchiut 1.30 if ( !strcmp(calargv[ri],"--ignore-h20") ) {
289     checksimu = false;
290     };
291 mocchiut 1.53 if ( !strcmp(calargv[ri],"--level1-only") ) {
292     cl1only = true;
293     withtrk = false;
294     getl1 = true;
295     }
296     if ( !strcmp(calargv[ri],"--from-level1") ) {
297     froml1 = true;
298     }
299 mocchiut 1.11 if ( !strcmp(calargv[ri],"--help") ) {
300     printf("\n\n CALORIMETER HELP CALLED\n\n");
301     printf(" CaloCore options: \n");
302     printf(" -v | --verbose be verbose\n");
303     printf(" -g | --debug be really verbose\n");
304     printf(" --defrig rig rig is the default rigidity in GV to be used to\n");
305     printf(" obtain calorimeter variables in the routines\n");
306     printf(" \"alltracks\" and \"higZnuclei\" [default = 50]\n");
307     printf(" --alltracks fill the track related variables even in the case\n");
308     printf(" of no tracks from tracker and no selftrigger event\n");
309 mocchiut 1.14 printf(" when we have a calorimeter fit for both views [default]\n");
310     printf(" --no-alltracks fill the track related variables only in the case\n");
311     printf(" of a good track from tracker or selftrigger\n");
312 mocchiut 1.11 printf(" --highZnuclei call the routine to analyze high Z nuclei\n");
313     printf(" selftrigger events [default]\n");
314     printf(" --no-highZnuclei do not call the routine to analyze high Z nuclei\n");
315     printf(" selftrigger events\n");
316     printf(" --no-tracker do not use tracker level2\n");
317     printf(" --with-tracker use tracker level2 [default]\n");
318 mocchiut 1.21 printf(" --no-crosstalk do not apply crosstalk corrections\n");
319 mocchiut 1.37 printf(" --ground-crosstalk apply ground crosstalk corrections\n");
320     printf(" --flight-crosstalk apply flight crosstalk corrections [default]\n");
321 mocchiut 1.41 printf(" --no-self-crosstalk do not apply preamplifiers crosstalk corrections to the strip itself\n");
322     printf(" --ct--use-pulse flight crosstalk corrections from pulse calibrations \n");
323     printf(" --ct--use-table flight crosstalk corrections from calibration table [default] \n");
324 mocchiut 1.11 printf(" --selftrigger process selftrigger events [default]\n");
325     printf(" --no-selftrigger skip selftrigger events\n");
326 mocchiut 1.53 printf(" --delete-calo-tree remove the ''Calorimeter'' TTree if found in the file - could cause fatal problems when reprocessing, use it with care!\n");
327     printf(" --level1 save Level1 TBranch [default] \n");
328 mocchiut 1.13 printf(" --no-level1 do not save Level1 TBranch\n");
329 mocchiut 1.53 printf(" --level1-only process level1 only (force --level1 flag)\n");
330     printf(" --from-level1 process level2 data starting from level1 already saved in the file\n");
331 mocchiut 1.30 printf(" --ignore-h20 do not check if it is a file from simulations\n");
332 mocchiut 1.53 printf(" --unmask-plane18X process plane 18X as well [default]\n");
333     printf(" --mask-plane18X do not process plane 18X\n");
334 mocchiut 1.27 printf(" --use-default-alig use default mechanical alignement (defined in CaloLevel1.h)\n");
335 mocchiut 1.11 throw -114;
336     };
337 mocchiut 1.1 ri++;
338     };
339     };
340     //
341 mocchiut 1.13 if ( verbose ){
342     printf("\n");
343     if ( getl1 ) printf(" Saving calorimeter level1 data \n");
344 mocchiut 1.53 if ( cl1only ) printf(" Saving calorimeter level1 data ONLY \n");
345     if ( froml1 ) printf(" Using level1 data saved in the file\n");
346 mocchiut 1.13 if ( st ) printf(" Calling selftrigger subroutine \n");
347     if ( hZn ) printf(" Calling high energy nuclei subroutine \n");
348     if ( trackanyway ) printf(" Filling track related variables for all the possible tracks \n");
349     if ( hZn || trackanyway ) printf(" => default assumed rigidity %f \n",rigdefault);
350     if ( withtrk ) printf(" Using tracker level2 data \n");
351 mocchiut 1.27 if ( mechal ) printf(" Using default mechanical alignement (defined in CaloLevel1.h)\n");
352 mocchiut 1.34 if ( crosst ){
353     printf(" Applying cross-talk corrections \n");
354     if ( ctground ){
355     printf(" => Using ground cross-talk coefficients \n");
356     } else {
357 mocchiut 1.41 printf(" => Using flight cross-talk coefficients... \n");
358     };
359     if ( usetable ){
360     printf(" ... coming from tables \n");
361     } else {
362     printf(" ... coming from pulse calibrations \n");
363 mocchiut 1.34 };
364     } else {
365     printf(" Do not applying cross-talk corrections \n");
366     };
367     if ( !checksimu ) printf(" Check on h20 TTree in level2 file skipped \n");
368 mocchiut 1.13 printf("\n");
369     };
370 mocchiut 1.11 //
371 mocchiut 1.1 // Working filename
372     //
373     TString outputfile;
374     stringstream name;
375     name.str("");
376     name << outdir << "/";
377     //
378     // Variables.
379     //
380     TTree *tracker = 0;
381     TTree *calo = 0;
382     TTree *caloclone = 0;
383     Bool_t reproc = false;
384     Bool_t reprocall = false;
385 mocchiut 1.8 UInt_t nevents = 0;
386 mocchiut 1.1 UInt_t nobefrun = 0;
387     UInt_t noaftrun = 0;
388     UInt_t numbofrun = 0;
389 mocchiut 1.8 UInt_t totnorun = 0;
390 mocchiut 1.1 //
391     Int_t code = 0;
392     Int_t sgnl;
393     //
394     // calorimeter level2 classes
395     //
396 mocchiut 1.13 CaloLevel1 *c1 = 0;
397     CaloLevel1 *c1clone = 0;
398     if ( getl1 ){
399     c1 = new CaloLevel1();
400     c1clone = new CaloLevel1();
401     };
402     //
403     // calorimeter level2 classes
404     //
405 mocchiut 1.1 CaloLevel2 *ca = new CaloLevel2();
406     CaloLevel2 *caclone = new CaloLevel2();
407     //
408     TrkLevel2 *trk = new TrkLevel2();
409 mocchiut 1.8 Int_t nevtrkl2 = 0;
410 mocchiut 1.1 //
411     UInt_t procev = 0;
412     //
413 mocchiut 1.8 // define variables where to store the absolute run header and run trailer times (unsigned long long integers, when set to a number use to store the correct number).
414 mocchiut 1.1 //
415 mocchiut 1.8 UInt_t runheadtime = 0;
416     UInt_t runtrailtime = 0;
417 mocchiut 1.1 UInt_t evfrom = 0;
418     UInt_t evto = 0;
419 mocchiut 1.8 UInt_t totfileentries = 0;
420     UInt_t idRun = 0;
421 mocchiut 1.1 Int_t id_reg_run=-1;
422     stringstream ftmpname;
423     TString fname;
424     //
425     // define variables for opening and reading level0 file
426     //
427     TFile *l0File = 0;
428     TTree *l0tr = 0;
429 mocchiut 1.24 TTree *softinfo = 0;
430 mocchiut 1.1 TBranch *l0head = 0;
431     TBranch *l0calo = 0;
432     TBranch *l0trig = 0;
433 mocchiut 1.8 pamela::EventHeader *eh = 0;
434     pamela::PscuHeader *ph = 0;
435 mocchiut 1.1 pamela::trigger::TriggerEvent *trig = 0;
436 mocchiut 1.24 Int_t yodaver = 0;
437 mocchiut 1.1 //
438     // Define some basic variables
439     //
440 mocchiut 1.22 CaloLevel0 *event = new CaloLevel0(); // NOTICE: very important to call here the constructor!
441 mocchiut 1.21 event->SetCrossTalk(crosst);
442 mocchiut 1.41 // event->SetCrossTalkType(ctground);
443     Int_t ict = -1;
444     if ( ctground ) ict = 0;
445     if ( !ctground ){
446     ict = 1;
447     };
448     if ( !ctground && noselfct ){
449     ict = 2;
450     };
451     event->SetCrossTalkType(ict);
452 mocchiut 1.1 stringstream file2;
453     stringstream file3;
454     stringstream qy;
455     // Bool_t imtrack = false;
456     Bool_t filled = false;
457     //
458 mocchiut 1.8 UInt_t caloevents = 0;
459 mocchiut 1.1 stringstream calfile;
460     stringstream aligfile;
461     //
462     Int_t i = -1;
463     Int_t itr = -1;
464     Int_t badevent = 0;
465     Int_t totevent = 0;
466     //
467 mocchiut 1.8 UInt_t atime = 0;
468 mocchiut 1.1 //
469 mocchiut 1.47 // Int_t S3 = 0;
470     // Int_t S2 = 0;
471     // Int_t S12 = 0;
472     // Int_t S11 = 0;
473 mocchiut 1.1 UInt_t re = 0;
474 mocchiut 1.9 UInt_t jumped = 0;
475 mocchiut 1.1 //
476     TString caloversion;
477     ItoRunInfo *runinfo = 0;
478 mocchiut 1.8 TArrayI *runlist = 0;
479 mocchiut 1.1 //
480 mocchiut 1.47 // Float_t tmptrigty = -1.;
481 mocchiut 1.1 Int_t ntrkentry = 0;
482     GL_PARAM *q4 = new GL_PARAM();
483 mocchiut 1.8 UInt_t tttrkpar1 = 0;
484 mocchiut 1.1 Bool_t trkpar1 = true;
485     GL_ROOT *glroot = new GL_ROOT();
486 mocchiut 1.8 GL_TIMESYNC *dbtime = 0;
487 mocchiut 1.1 //
488 mocchiut 1.15 Long64_t maxsize = 10000000000LL;
489     TTree::SetMaxTreeSize(maxsize);
490 mocchiut 1.1 //
491     //
492     TFile *tempfile = 0;
493     TTree *tempcalo = 0;
494 mocchiut 1.30 Bool_t myfold = false;
495 mocchiut 1.1 stringstream tempname;
496     stringstream calofolder;
497     tempname.str("");
498     tempname << outdir;
499     tempname << "/" << processFolder.Data();
500     calofolder.str("");
501     calofolder << tempname.str().c_str();
502     tempname << "/calotree_run";
503     tempname << run << ".root";
504     //
505     // As a first thing we must check what we have to do: if run = -1 we must process all events in the file has been passed
506     // if run != -1 we must process only that run but first we have to check if the branch calorimeter already exist in the file
507     // if it exists we are reprocessing data and we must delete that entries, if not we must create it.
508     //
509 mocchiut 1.8 if ( run == 0 ) reproc = true;
510 mocchiut 1.1 //
511     //
512     //
513     if ( !file->IsOpen() ){
514     if ( verbose ) printf(" CALORIMETER - ERROR: cannot open file for writing\n");
515     throw -101;
516     };
517     //
518 mocchiut 1.53 if ( delct ){
519     TTree *T = (TTree*)file->Get("Calorimeter");
520     if ( T ){
521 mocchiut 1.54 if ( verbose ) printf(" CALORIMETER - REMOVING Calorimeter TTree \n");
522 mocchiut 1.53 T->Delete("all");
523     }
524     }
525     //
526 mocchiut 1.55 TClonesArray *tcNucleiTrk = NULL;
527     TClonesArray *tcExtNucleiTrk = NULL;
528     TClonesArray *tcExtTrk = NULL;
529     TClonesArray *tcaloNucleiTrk = NULL;
530     TClonesArray *tcaloExtNucleiTrk = NULL;
531     TClonesArray *tcaloExtTrk = NULL;
532     Bool_t hasNucleiTrk = false;
533     Bool_t hasExtNucleiTrk = false;
534     Bool_t hasExtTrk = false;
535 mocchiut 1.11 if ( withtrk ){
536     //
537     // Does it contain the Tracker tree?
538     //
539     tracker = (TTree*)file->Get("Tracker");
540     if ( !tracker ) {
541     if ( verbose ) printf(" CALORIMETER - ERROR: no tracker tree\n");
542     code = -102;
543     goto closeandexit;
544 mocchiut 1.55 }
545 mocchiut 1.57 // tracker->SetMaxVirtualSize(2500000000LL); // EM residual Tracker-new tree in level2 files when NEVENTS is big
546 mocchiut 1.11 tracker->SetBranchAddress("TrkLevel2",&trk);
547     nevtrkl2 = tracker->GetEntries();
548 mocchiut 1.55 //
549     // Look for extended tracking algorithm
550     //
551     if ( verbose ) printf("Look for extended and nuclei tracking algorithms\n");
552     // Nuclei tracking algorithm
553     Int_t checkAlgo = 0;
554     tcNucleiTrk = new TClonesArray("TrkTrack");
555     checkAlgo = tracker->SetBranchAddress("TrackNuclei",&tcNucleiTrk);
556     if ( !checkAlgo ){
557     if ( verbose ) printf(" Nuclei tracking algorithm branch found! :D \n");
558     hasNucleiTrk = true;
559     } else {
560     if ( verbose ) printf(" Nuclei tracking algorithm branch not found :( !\n");
561     printf(" ok, this is not a problem (it depends on tracker settings) \n");
562     delete tcNucleiTrk;
563 mocchiutti 1.60 tcNucleiTrk=NULL; // 10RED reprocessing bug
564 mocchiut 1.55 }
565     // Nuclei tracking algorithm using calorimeter points
566     tcExtNucleiTrk = new TClonesArray("ExtTrack");
567     checkAlgo = tracker->SetBranchAddress("RecoveredTrackNuclei",&tcExtNucleiTrk);
568     if ( !checkAlgo ){
569     if ( verbose ) printf(" Recovered nuclei tracking algorithm branch found! :D \n");
570     hasExtNucleiTrk = true;
571     } else {
572     if ( verbose ) printf(" Recovered nuclei tracking algorithm branch not found :( !\n");
573     printf(" ok, this is not a problem (it depends on tracker settings) \n");
574     delete tcExtNucleiTrk;
575 mocchiutti 1.60 tcExtNucleiTrk=NULL; // 10RED reprocessing bug
576 mocchiut 1.55 }
577     // Tracking algorithm using calorimeter points
578     tcExtTrk = new TClonesArray("ExtTrack");
579     checkAlgo = tracker->SetBranchAddress("RecoveredTrack",&tcExtTrk);
580     if ( !checkAlgo ){
581     if ( verbose ) printf(" Recovered track algorithm branch found! :D \n");
582     hasExtTrk = true;
583     } else {
584     if ( verbose ) printf(" Recovered track algorithm branch not found :( !\n");
585     printf(" ok, this is not a problem (it depends on tracker settings) \n");
586     delete tcExtTrk;
587 mocchiutti 1.60 tcExtTrk=NULL; // 10RED reprocessing bug
588 mocchiut 1.55 }
589     }
590 mocchiut 1.1 //
591 mocchiut 1.30 // Is it a file from simulations?
592     //
593     if ( checksimu ){
594     TTree *h20 = (TTree*)file->Get("h20");
595     if ( h20 ){
596     //
597     // yes!
598     //
599     crosst = false;
600 mocchiut 1.34 ctground = true;
601 mocchiut 1.30 mechal = true;
602 mocchiut 1.47 mask18 = true;
603     if ( verbose ) printf("\n\n SIMULATION!!! Setting --use-default-alig and --no-crosstalk --mask-plane18X flags!\n\n");
604 mocchiut 1.30 h20->Delete();
605     };
606     };
607     //
608 mocchiut 1.1 // Call runinfo
609     //
610     sgnl = 0;
611     runinfo = new ItoRunInfo(file);
612     //
613     // update versioning information and retrieve informations about the run to be processed
614     //
615     caloversion = CaloInfo(false);
616     sgnl = runinfo->Update(run,"CALO",caloversion);
617     //
618     if ( sgnl ){
619     if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
620     code = sgnl;
621     goto closeandexit;
622     } else {
623     sgnl = 0;
624     };
625     //
626     // number of events in the file BEFORE the first event of our run
627     //
628     nobefrun = runinfo->GetFirstEntry();
629     //
630     // total number of events in the file
631     //
632     totfileentries = runinfo->GetFileEntries();
633     //
634     // first file entry AFTER the last event of our run
635     //
636     noaftrun = runinfo->GetLastEntry() + 1;
637     //
638     // number of run to be processed
639     //
640     numbofrun = runinfo->GetNoRun();
641     //
642     // number of runs in the file
643     //
644     totnorun = runinfo->GetRunEntries();
645     //
646     // Does it contain already a Calorimeter branch? if so we are reprocessing data, if not we must create it.
647     //
648     caloclone = (TTree*)file->Get("Calorimeter");
649     //
650     if ( !caloclone ){
651     reproc = false;
652 mocchiut 1.8 if ( run == 0 && verbose ) printf(" CALORIMETER - WARNING: you are reprocessing data but calorimeter tree does not exist!\n");
653     if ( runinfo->IsReprocessing() && run != 0 && verbose ) printf(" CALORIMETER - WARNING: it seems you are not reprocessing data but calorimeter\n versioning information already exists in RunInfo.\n");
654 mocchiut 1.1 //
655     } else {
656     //
657 mocchiut 1.57 // caloclone->SetMaxVirtualSize(2500000000LL); // EM residual Tracker-new tree in level2 files when NEVENTS is big
658 mocchiut 1.15 caloclone->SetAutoSave(900000000000000LL);
659 mocchiut 1.1 reproc = true;
660     //
661     if ( verbose ) printf("\n Preparing the pre-processing...\n");
662     //
663 mocchiut 1.28 if ( run == 0 || totnorun == 1 ){
664 mocchiut 1.1 //
665     // if we are reprocessing everything we don't need to copy any old event and we can just create a new branch in the clone tree and jump steps 4/5/7.
666     //
667     if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing all runs in the file\n");
668     reprocall = true;
669     //
670     } else {
671     //
672     // we are reprocessing a single run
673     //
674 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing run number %u \n",run);
675 mocchiut 1.1 reprocall = false;
676     //
677     //
678     //
679 mocchiut 1.30 gSystem->MakeDirectory(calofolder.str().c_str());
680     myfold = true;
681 mocchiut 1.1 tempfile = new TFile(tempname.str().c_str(),"RECREATE");
682     tempcalo = caloclone->CloneTree(-1,"fast");
683     tempcalo->SetName("Calorimeter-old");
684     tempfile->Write();
685     tempfile->Close();
686     };
687     //
688     // delete old tree
689     //
690 mocchiut 1.53 if ( debug ) printf(" Delete caloclone \n");
691     if ( caloclone ) caloclone->Delete("all");
692 mocchiut 1.1 //
693     if ( verbose ) printf("\n ...done!\n");
694     //
695 mocchiut 1.59 }
696 mocchiut 1.1 //
697     // create calorimeter tree calo
698     //
699     file->cd();
700     calo = new TTree("Calorimeter-new","PAMELA Level2 calorimeter data");
701 mocchiut 1.57 // calo->SetMaxVirtualSize(2500000000LL); // EM residual Tracker-new tree in level2 files when NEVENTS is big
702 mocchiut 1.15 calo->SetAutoSave(900000000000000LL);
703 pam-fi 1.16 ca->Set();
704 mocchiut 1.59 if ( !cl1only ) calo->Branch("CaloLevel2","CaloLevel2",&ca);
705     //calo->Branch("CaloLevel2","CaloLevel2",&ca);
706 mocchiut 1.13 if ( getl1 ) calo->Branch("CaloLevel1","CaloLevel1",&c1);
707 mocchiut 1.1 //
708 mocchiut 1.55 // create new branches for new tracking algorithms
709     //
710     if ( hasNucleiTrk ){
711     tcaloNucleiTrk = new TClonesArray("CaloTrkVar",1);
712     calo->Branch("TrackNuclei",&tcaloNucleiTrk);
713     }
714     if ( hasExtNucleiTrk ){
715     tcaloExtNucleiTrk = new TClonesArray("CaloTrkVar",1);
716     calo->Branch("RecoveredTrackNuclei",&tcaloExtNucleiTrk);
717     }
718     if ( hasExtTrk ){
719     tcaloExtTrk = new TClonesArray("CaloTrkVar",1);
720     calo->Branch("RecoveredTrack",&tcaloExtTrk);
721     }
722    
723     //
724 mocchiut 1.1 if ( reproc && !reprocall ){
725     //
726     //
727     //
728 mocchiut 1.59 tempfile = TFile::Open(tempname.str().c_str(),"READ");
729 mocchiut 1.1 caloclone = (TTree*)tempfile->Get("Calorimeter-old");
730 mocchiut 1.57 // caloclone->SetMaxVirtualSize(2500000000LL); // EM residual Tracker-new tree in level2 files when NEVENTS is big
731 mocchiut 1.15 caloclone->SetAutoSave(900000000000000LL);
732 mocchiut 1.58 caloclone->SetBranchAddress("CaloLevel2",&caclone);// EM reprocessing bug fixed
733 mocchiut 1.13 TBranch *havel1 = caloclone->GetBranch("CaloLevel1");
734     if ( !havel1 && getl1 ) throw -117;
735     if ( havel1 && !getl1 ) throw -118;
736     if ( getl1 ) caloclone->SetBranchAddress("CaloLevel1",&c1clone);
737 mocchiut 1.1 //
738 mocchiut 1.53 if ( nobefrun > 0 ){
739 mocchiut 1.1 if ( verbose ) printf("\n Pre-processing: copying events from the old tree before the processed run\n");
740 mocchiut 1.8 if ( verbose ) printf(" Copying %u events in the file which are before the beginning of the run %u \n",nobefrun,run);
741 mocchiut 1.1 if ( verbose ) printf(" Start copying at event number 0, end copying at event number %u \n",nobefrun);
742     for (UInt_t j = 0; j < nobefrun; j++){
743     //
744 mocchiut 1.48 if ( caloclone->GetEntry(j) <= 0 ) throw -36;
745 mocchiut 1.1 //
746     // copy caclone to ca
747     //
748 mocchiut 1.13 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
749 mocchiut 1.58 memcpy(&ca,&caclone,sizeof(caclone));// EM reprocessing bug fixed
750 mocchiut 1.1 //
751     // Fill entry in the new tree
752     //
753 mocchiut 1.58 calo->Fill();// EM reprocessing bug fixed
754 mocchiut 1.1 //
755 mocchiut 1.13 if ( getl1 ) c1->Clear();
756 mocchiut 1.6 ca->Clear();
757     //
758 mocchiut 1.59 }
759 mocchiut 1.1 if ( verbose ) printf(" Finished successful copying!\n");
760 mocchiut 1.59 }
761 mocchiut 1.53 }
762 mocchiut 1.1 //
763     // Get the list of run to be processed
764     //
765     runlist = runinfo->GetRunList();
766     //
767     // Loop over the run to be processed
768     //
769     for (UInt_t irun=0; irun < numbofrun; irun++){
770     //
771     badevent = 0;
772     //
773     idRun = runlist->At(irun);
774     if ( verbose ) printf("\n\n\n ####################################################################### \n");
775 mocchiut 1.8 if ( verbose ) printf(" PROCESSING RUN NUMBER %u \n",idRun);
776 mocchiut 1.1 if ( verbose ) printf(" ####################################################################### \n\n\n");
777     //
778     sgnl = runinfo->GetRunInfo(idRun);
779     if ( sgnl ){
780     if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
781     code = sgnl;
782     goto closeandexit;
783     } else {
784     sgnl = 0;
785     };
786 mocchiut 1.8 id_reg_run = runinfo->ID_ROOT_L0;
787 mocchiut 1.1 runheadtime = runinfo->RUNHEADER_TIME;
788     runtrailtime = runinfo->RUNTRAILER_TIME;
789 mocchiut 1.8 evfrom = runinfo->EV_FROM;
790     evto = runinfo->EV_TO;
791 mocchiut 1.1 //
792     if ( id_reg_run == -1 ){
793 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - ERROR: no run with ID_RUN = %u \n\n Exiting... \n\n",idRun);
794 mocchiut 1.1 code = -5;
795     goto closeandexit;
796     };
797     //
798 mocchiut 1.8 // prepare the timesync for the db
799     //
800 mocchiut 1.32 TString host = glt->CGetHost();
801     TString user = glt->CGetUser();
802     TString psw = glt->CGetPsw();
803     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
804     if ( !dbc->IsConnected() ) throw -116;
805 mocchiut 1.33 stringstream myquery;
806     myquery.str("");
807     myquery << "SET time_zone='+0:00'";
808 mocchiut 1.52 delete dbc->Query(myquery.str().c_str());
809 mocchiut 1.32 //
810 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
811 mocchiut 1.8 dbtime = new GL_TIMESYNC(runinfo->ID_ROOT_L0,"ID",dbc);
812     //
813 mocchiut 1.1 // Search in the DB the path and name of the LEVEL0 file to be processed.
814     //
815 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
816 mocchiut 1.8 glroot->Query_GL_ROOT(runinfo->ID_ROOT_L0,dbc);
817 mocchiut 1.1 //
818     ftmpname.str("");
819     ftmpname << glroot->PATH.Data() << "/";
820     ftmpname << glroot->NAME.Data();
821     fname = ftmpname.str().c_str();
822     //
823     // print out informations
824     //
825 mocchiut 1.8 totevent = runinfo->NEVENTS;
826 mocchiut 1.1 if ( verbose ) printf("\n LEVEL0 data file: %s \n",fname.Data());
827 mocchiut 1.8 if ( verbose ) printf(" RUN HEADER absolute time is: %u \n",runheadtime);
828     if ( verbose ) printf(" RUN TRAILER absolute time is: %u \n",runtrailtime);
829     if ( verbose ) printf(" %i events to be processed for run %u: from %i to %i (reg entries)\n\n",totevent,idRun,evfrom,evfrom+totevent);
830 mocchiut 1.1 //
831 mocchiut 1.43 // if ( !totevent ) goto closeandexit;
832 mocchiut 1.42 //
833 mocchiut 1.1 // Open Level0 file
834     //
835 mocchiut 1.52 if ( l0File ) l0File->Close();
836 mocchiut 1.1 l0File = new TFile(fname.Data());
837     if ( !l0File ) {
838     if ( verbose ) printf(" CALORIMETER - ERROR: problems opening Level0 file\n");
839     code = -6;
840     goto closeandexit;
841     };
842     l0tr = (TTree*)l0File->Get("Physics");
843     if ( !l0tr ) {
844     if ( verbose ) printf(" CALORIMETER - ERROR: no Physics tree in Level0 file\n");
845     l0File->Close();
846     code = -7;
847     goto closeandexit;
848     };
849     l0head = l0tr->GetBranch("Header");
850     if ( !l0head ) {
851     if ( verbose ) printf(" CALORIMETER - ERROR: no Header branch in Level0 tree\n");
852     l0File->Close();
853     code = -8;
854     goto closeandexit;
855     };
856     l0calo = l0tr->GetBranch("Calorimeter");
857     if ( !l0calo ) {
858     if ( verbose ) printf(" CALORIMETER - ERROR: no Calorimeter branch in Level0 tree\n");
859     l0File->Close();
860     code = -103;
861     goto closeandexit;
862     };
863     l0trig = l0tr->GetBranch("Trigger");
864     if ( !l0trig ) {
865     if ( verbose ) printf(" CALORIMETER - ERROR: no Trigger branch in Level0 tree\n");
866     l0File->Close();
867     code = -104;
868     goto closeandexit;
869     };
870     //
871     l0tr->SetBranchAddress("Trigger", &trig);
872 mocchiut 1.8 l0tr->SetBranchAddress("Header", &eh);
873 mocchiut 1.1 //
874 mocchiut 1.24 softinfo = (TTree*)l0File->Get("SoftInfo");
875     if ( softinfo ){
876     softinfo->SetBranchAddress("SoftInfo",&yodaver);
877 mocchiut 1.48 if ( softinfo->GetEntry(0) <= 0 ) throw -36;
878 mocchiut 1.24 };
879     if ( debug ) printf(" LEVEL0 FILE GENERATED WITH YODA VERSION %i \n",yodaver);
880     //
881 mocchiut 1.38 if ( withtrk ){
882     if ( trkpar1 || ( tttrkpar1 != 0 && tttrkpar1 < runheadtime ) ){
883     trkpar1 = false;
884     Int_t glpar = q4->Query_GL_PARAM(runinfo->RUNHEADER_TIME,1,dbc);
885     if ( glpar < 0 ){
886     code = glpar;
887     goto closeandexit;
888     };
889     tttrkpar1 = q4->TO_TIME;
890     // ----------------------------
891     // Read the magnetic field
892     // ----------------------------
893     if ( verbose ) printf(" Reading magnetic field maps at %s\n",(q4->PATH+q4->NAME).Data());
894     trk->LoadField(q4->PATH+q4->NAME);
895     if ( verbose ) printf("\n");
896     };
897     };
898     //
899     // Close DB connection
900     //
901     if ( dbc ){
902     dbc->Close();
903     delete dbc;
904 pam-fi 1.50 dbc = 0;
905 mocchiut 1.38 };
906     //
907     //
908 mocchiut 1.1 // Construct the event object, look for the calibration which include the first header
909     //
910     sgnl = 0;
911     if ( verbose ) printf(" Check for calorimeter calibrations and initialize event object \n");
912 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
913 mocchiut 1.32 event->ProcessingInit(glt,runheadtime,sgnl,l0tr,debug,verbose);
914 mocchiut 1.47 if ( mask18 ){
915     event->MaskPlane18X();
916     } else {
917     event->UnMaskPlane18X();
918     };
919 mocchiut 1.1 if ( verbose ) printf("\n");
920     if ( sgnl == 100 ) {
921     code = sgnl;
922     if ( verbose ) printf(" CALORIMETER - WARNING: run header not included in any calibration interval\n");
923     sgnl = 0;
924     };
925     if ( sgnl ){
926     l0File->Close();
927     code = sgnl;
928     goto closeandexit;
929     };
930     //
931     qy.str("");
932     //
933 mocchiut 1.8 nevents = l0head->GetEntries();
934 mocchiut 1.1 caloevents = l0calo->GetEntries();
935     //
936 mocchiut 1.43 if ( nevents < 1 && totevent ) {
937 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: Level0 file is empty\n\n");
938     l0File->Close();
939     code = -11;
940     goto closeandexit;
941     };
942     //
943 mocchiut 1.43 if ( evto > nevents-1 && totevent ) {
944 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: too few entries in the registry tree\n");
945     l0File->Close();
946     code = -12;
947     goto closeandexit;
948     };
949     //
950     // Check if we have to load parameter files
951     //
952     sgnl = 0;
953 mocchiut 1.32 sgnl = event->ChkParam(glt,runheadtime,mechal); // calorimeter parameter files
954 mocchiut 1.1 if ( sgnl < 0 ){
955     code = sgnl;
956     l0File->Close();
957     goto closeandexit;
958     };
959     //
960 mocchiut 1.34 // Calculate the cross talk corrections needed for this run using calibration informations
961     //
962     if ( crosst && !ctground ){
963     sgnl = 0;
964 mocchiut 1.41 sgnl = event->CalcCrossTalkCorr(glt,runheadtime,usetable);
965 mocchiut 1.34 if ( sgnl < 0 ){
966     code = sgnl;
967     l0File->Close();
968     goto closeandexit;
969     };
970     };
971     //
972 mocchiut 1.1 // run over all the events of the run
973     //
974     if ( verbose ) printf("\n Ready to start! \n\n Processed events: \n\n");
975     //
976 mocchiut 1.19 for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+runinfo->NEVENTS); re++){
977 mocchiut 1.44 // for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+5000); re++){
978 mocchiut 1.49 //for ( re = 92012+3400; re < 92012+3700; re++){
979     //for ( re = 4451+2833; re < 4451+2836 ; re++){
980 mocchiut 1.18 // printf(" i = %i \n",re-200241);
981 mocchiut 1.1 //
982     if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000);
983     //
984 mocchiut 1.39 if ( debug ) printf("\n\n\n EVENT NUMBER %i \n",procev);
985     //
986 mocchiut 1.48 if ( l0head->GetEntry(re) <= 0 ) throw -36;
987 mocchiut 1.1 //
988     // absolute time of this event
989     //
990 mocchiut 1.8 ph = eh->GetPscuHeader();
991     atime = dbtime->DBabsTime(ph->GetOrbitalTime());
992 mocchiut 1.1 //
993     //
994     //
995 mocchiut 1.8 if ( re > caloevents-1 ){
996 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: no physics events with entry = %i in Level0 file\n",i);
997     l0File->Close();
998     code = -112;
999     goto closeandexit;
1000     };
1001     //
1002 mocchiut 1.45 if ( atime > (runtrailtime+1) || atime < (runheadtime-1) ) {
1003 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - WARNING: event at time outside the run time window, skipping it\n");
1004 mocchiut 1.9 jumped++;
1005 mocchiut 1.1 goto jumpev;
1006     };
1007     //
1008     // retrieve tracker informations
1009     //
1010     if ( !reprocall ){
1011 mocchiut 1.19 itr = nobefrun + (re - evfrom - jumped);
1012     //itr = re-(46438+200241);
1013 mocchiut 1.1 } else {
1014 mocchiut 1.9 itr = runinfo->GetFirstEntry() + (re - evfrom - jumped);
1015 mocchiut 1.1 };
1016 mocchiut 1.11 //
1017     if ( withtrk ){
1018     if ( itr > nevtrkl2 ){
1019     if ( verbose ) printf(" CALORIMETER - ERROR: no tracker events with entry = %i in Level2 file\n",itr);
1020     l0File->Close();
1021     code = -113;
1022     goto closeandexit;
1023 mocchiut 1.55 }
1024 mocchiut 1.11 //
1025     trk->Clear();
1026     //
1027 mocchiut 1.55 //
1028     // Clones array must be cleared before going on
1029     //
1030     if ( hasNucleiTrk ){
1031     tcNucleiTrk->Delete();
1032     tcaloNucleiTrk->Delete();
1033     }
1034     if ( hasExtNucleiTrk ){
1035     tcExtNucleiTrk->Delete();
1036     tcaloExtNucleiTrk->Delete();
1037     }
1038     if ( hasExtTrk ){
1039     tcExtTrk->Delete();
1040     tcaloExtTrk->Delete();
1041     }
1042     //
1043 mocchiut 1.48 if ( tracker->GetEntry(itr) <= 0 ) throw -36;
1044 mocchiut 1.55 //
1045     }
1046 mocchiut 1.10 //
1047 mocchiut 1.1 procev++;
1048     //
1049     // start processing
1050     //
1051 mocchiut 1.13 if ( getl1 ) c1->Clear();
1052 mocchiut 1.6 ca->Clear();
1053 mocchiut 1.1 //
1054     // determine who generate the trigger for this event (TOF, S4/PULSER, CALO)
1055     //
1056 mocchiut 1.48 if ( l0trig->GetEntry(re) <= 0 ) throw -36;
1057 mocchiut 1.47 // S3 = 0;
1058     // S2 = 0;
1059     // S12 = 0;
1060     // S11 = 0;
1061     // S3 = trig->patterntrig[2];
1062     // S2 = trig->patterntrig[3];
1063     // S12 = trig->patterntrig[4];
1064     // S11 = trig->patterntrig[5];
1065     // if ( trig->patterntrig[1] & (1<<0) ) tmptrigty = 1.;
1066     // if ( trig->patterntrig[0] ) tmptrigty = 2.;
1067     // if ( S3 || S2 || S12 || S11 ) tmptrigty = 0.;
1068     // if ( !(trig->patterntrig[1] & (1<<0)) && !trig->patterntrig[0] && !S3 && !S2 && !S12 && !S11 ) tmptrigty = 1.;
1069     // event->clevel2->trigty = tmptrigty;
1070     event->clevel2->trigty = (Float_t)IsCaloSelfTrigger(trig->trigconf,(UInt_t*)trig->patterntrig);
1071     if ( debug ) printf(" starting with trigty %f \n",event->clevel2->trigty);
1072 mocchiut 1.1 //
1073     // check if the calibration we are using is still good, if not load another calibration
1074     //
1075     sgnl = 0;
1076 mocchiut 1.32 sgnl = event->ChkCalib(glt,atime);
1077 mocchiut 1.1 if ( sgnl < 0 ){
1078     code = sgnl;
1079     goto closeandexit;
1080     };
1081     if ( sgnl == 100 ){
1082     code = sgnl;
1083     if ( verbose ) printf(" CALORIMETER - WARNING: data not associated to any calibration interval\n");
1084     badevent++;
1085     sgnl = 0;
1086     };
1087     //
1088 mocchiut 1.11 // do we have at least one track from the tracker? this check has been disabled
1089 mocchiut 1.1 //
1090     event->clevel1->good2 = 1;
1091 mocchiut 1.8 //
1092 mocchiut 1.29 // use the whole calorimeter, 22 W planes
1093     //
1094     event->clevel1->npla = 22;
1095     //
1096     // Use the standard silicon planes shifting
1097     //
1098     event->clevel1->reverse = 0;
1099     //
1100 mocchiut 1.1 //
1101 mocchiut 1.8 // Calibrate calorimeter event "re" and store output in the two structures that will be passed to fortran routine
1102 mocchiut 1.1 //
1103 mocchiut 1.11 event->Calibrate(re);
1104 mocchiut 1.1 //
1105     // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract topological variables.
1106     //
1107     //
1108     // Calculate variables common to all tracks (qtot, nstrip, etc.)
1109     //
1110 mocchiut 1.53 if ( !cl1only ) event->GetCommonVar();
1111 mocchiut 1.1 //
1112     // Fill common variables
1113     //
1114 mocchiut 1.53 if ( cl1only ){
1115     event->FillCommonVar(c1);
1116     } else {
1117     event->FillCommonVar(c1,ca);
1118     }
1119 mocchiut 1.1 //
1120     // Calculate variables related to tracks only if we have at least one track (from selftrigger and/or tracker)
1121     //
1122     ntrkentry = 0;
1123     //
1124     filled = false;
1125 mocchiut 1.11 //
1126 mocchiut 1.53 // skip track-related if processing level1 only
1127 mocchiut 1.1 //
1128 mocchiut 1.53 if ( !cl1only ){
1129     //
1130     // Run over tracks (tracker or calorimeter )
1131     //
1132     if ( withtrk ){
1133     for(Int_t nt=0; nt < trk->ntrk(); nt++){
1134     //
1135     event->clevel1->good2 = 1;
1136     //
1137     TrkTrack *ptt = trk->GetStoredTrack(nt);
1138     //
1139     event->clevel1->trkchi2 = 0;
1140     //
1141     // Copy the alpha vector in the input structure
1142     //
1143     for (Int_t e = 0; e < 5 ; e++){
1144     event->clevel1->al_p[e][0] = ptt->al[e];
1145     };
1146     //
1147     // Get tracker related variables for this track
1148     //
1149     event->GetTrkVar();
1150     //
1151     // Save tracker track sequence number
1152     //
1153     event->trkseqno = nt;
1154     //
1155     // Copy values in the class ca from the structure clevel2
1156     //
1157     event->FillTrkVar(ca,ntrkentry);
1158     ntrkentry++;
1159     filled = true;
1160     //
1161 mocchiut 1.55 } // loop on all the tracks
1162     //
1163     // Code for extended tracking algorithm:
1164     //
1165     if ( hasNucleiTrk ){
1166     Int_t ttentry = 0;
1167     for(Int_t nt=0; nt < tcNucleiTrk->GetEntries() ; nt++){
1168     //
1169     event->clevel1->good2 = 1;
1170     //
1171     TrkTrack *ptt = (TrkTrack*)(tcNucleiTrk->At(nt));
1172     //
1173     event->clevel1->trkchi2 = 0;
1174     //
1175     // Copy the alpha vector in the input structure
1176     //
1177     for (Int_t e = 0; e < 5 ; e++){
1178     event->clevel1->al_p[e][0] = ptt->al[e];
1179     };
1180     //
1181     // Get tracker related variables for this track
1182     //
1183     event->GetTrkVar();
1184     //
1185     // Save tracker track sequence number
1186     //
1187     event->trkseqno = nt;
1188     //
1189     // Copy values in the class ca from the structure clevel2
1190     //
1191     event->FillTrkVar(tcaloNucleiTrk,ttentry);
1192     ttentry++;
1193     //
1194     } // loop on all the tracks
1195     }
1196    
1197     if ( hasExtNucleiTrk ){
1198     Int_t ttentry = 0;
1199     for(Int_t nt=0; nt < tcExtNucleiTrk->GetEntries() ; nt++){
1200     //
1201     event->clevel1->good2 = 1;
1202     //
1203     ExtTrack *ptt = (ExtTrack*)(tcExtNucleiTrk->At(nt));
1204     //
1205     event->clevel1->trkchi2 = 0;
1206     //
1207     // Copy the alpha vector in the input structure
1208     //
1209     for (Int_t e = 0; e < 5 ; e++){
1210     event->clevel1->al_p[e][0] = ptt->al[e];
1211     };
1212     //
1213     // Get tracker related variables for this track
1214     //
1215     event->GetTrkVar();
1216     //
1217     // Save tracker track sequence number
1218     //
1219     event->trkseqno = nt;
1220     //
1221     // Copy values in the class ca from the structure clevel2
1222     //
1223     event->FillTrkVar(tcaloExtNucleiTrk,ttentry);
1224     ttentry++;
1225     //
1226     } // loop on all the tracks
1227     }
1228    
1229     if ( hasExtTrk ){
1230     Int_t ttentry = 0;
1231     for(Int_t nt=0; nt < tcExtTrk->GetEntries() ; nt++){
1232     //
1233     event->clevel1->good2 = 1;
1234     //
1235     ExtTrack *ptt = (ExtTrack*)(tcExtTrk->At(nt));
1236     //
1237     event->clevel1->trkchi2 = 0;
1238     //
1239     // Copy the alpha vector in the input structure
1240     //
1241     for (Int_t e = 0; e < 5 ; e++){
1242     event->clevel1->al_p[e][0] = ptt->al[e];
1243     };
1244     //
1245     // Get tracker related variables for this track
1246     //
1247     event->GetTrkVar();
1248     //
1249     // Save tracker track sequence number
1250     //
1251     event->trkseqno = nt;
1252     //
1253     // Copy values in the class ca from the structure clevel2
1254     //
1255     event->FillTrkVar(tcaloExtTrk,ttentry);
1256     ttentry++;
1257     //
1258     } // loop on all the tracks
1259     }
1260    
1261     }
1262 mocchiut 1.53 //
1263     // if no tracks found but there is the possibility to have a good track we should try to calculate anyway the track related variables using the calorimeter
1264     // fit of the track (to be used for example when TRK is off due to any reason like IPM3/5 off).
1265     // here we make an event selection so it must be done very carefully...
1266     //
1267     // conditions are: 0) no track from the tracker 1) we have a track fit both in x and y 2) no problems with calo for this event 3) no selftrigger event
1268     //
1269     if ( trackanyway && !filled && event->clevel2->npcfit[0] >= 2 && event->clevel2->npcfit[1] >= 2 && event->clevel2->good != 0 && event->clevel2->trigty < 2. ){
1270     if ( debug ) printf(" Event with a track not fitted by the tracker at entry %i \n",itr);
1271     //
1272     // Disable "track mode" in the fortran routine
1273     //
1274     event->clevel1->good2 = 0;
1275     event->clevel1->riginput = rigdefault;
1276     if ( debug ) printf(" Using as default rigidity: %f \n",event->clevel1->riginput);
1277     //
1278     // We have a selftrigger event to analyze.
1279     //
1280     for (Int_t e = 0; e < 5 ; e++){
1281     event->clevel1->al_p[e][0] = 0.;
1282     event->clevel1->al_p[e][1] = 0.;
1283     }
1284     event->clevel1->trkchi2 = 0;
1285     //
1286     event->GetTrkVar();
1287     //
1288     // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
1289     //
1290     if ( event->clevel1->good2 == 0 ) {
1291     //
1292     // In selftrigger mode the trkentry variable is set to -1
1293     //
1294     event->trkseqno = -3;
1295     //
1296     // Copy values in the class ca from the structure clevel2
1297     //
1298     if ( debug ) printf(" -3 begin \n");
1299     event->FillTrkVar(ca,ntrkentry);
1300     if ( debug ) printf(" -3 end \n");
1301     ntrkentry++;
1302     filled = true;
1303     //
1304     } else {
1305     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1306     }
1307     //
1308     }
1309     //
1310     // Call high energy nuclei routine
1311     //
1312     if ( hZn && event->clevel2->trigty >= 2. ){
1313     if ( debug ) printf(" Calling selftrigger high energy nuclei routine for entry %i \n",itr);
1314     //
1315     // Disable "track mode" in the fortran routine
1316     //
1317     event->clevel1->good2 = 0;
1318     //
1319     // Set high energy nuclei flag to one
1320     //
1321     event->clevel1->hzn = 1;
1322     event->clevel1->riginput = rigdefault;
1323     //
1324     // We have a selftrigger event to analyze.
1325     //
1326     for (Int_t e = 0; e < 5 ; e++){
1327     event->clevel1->al_p[e][0] = 0.;
1328     event->clevel1->al_p[e][1] = 0.;
1329     }
1330     event->clevel1->trkchi2 = 0;
1331     //
1332     event->GetTrkVar();
1333     //
1334     // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
1335     //
1336     if ( event->clevel1->good2 == 0 ) {
1337     //
1338     // In selftrigger mode the trkentry variable is set to -1
1339     //
1340     event->trkseqno = -2;
1341     //
1342     // Copy values in the class ca from the structure clevel2
1343     //
1344     event->FillTrkVar(ca,ntrkentry);
1345     ntrkentry++;
1346     filled = true;
1347     //
1348     } else {
1349     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1350     }
1351     //
1352     }
1353     //
1354     // self trigger event
1355     //
1356     if ( st && event->clevel2->trigty >= 2. ){
1357     if ( verbose ) printf(" Selftrigger event at entry %i \n",itr);
1358     //
1359     // Disable "track mode" in the fortran routine
1360     //
1361     event->clevel1->good2 = 0;
1362     //
1363     // disable high enery nuclei flag;
1364     //
1365     event->clevel1->hzn = 0;
1366     //
1367     // We have a selftrigger event to analyze.
1368     //
1369     for (Int_t e = 0; e < 5 ; e++){
1370     event->clevel1->al_p[e][0] = 0.;
1371     event->clevel1->al_p[e][1] = 0.;
1372     };
1373     event->clevel1->trkchi2 = 0;
1374     //
1375     event->GetTrkVar();
1376     //
1377     // if we had no problem (clevel2->good = 0, NOTICE zero, not one in selftrigger mode!), fill and go on
1378     //
1379     if ( event->clevel1->good2 == 0 ) {
1380     //
1381     // In selftrigger mode the trkentry variable is set to -1
1382     //
1383     event->trkseqno = -1;
1384     //
1385     // Copy values in the class ca from the structure clevel2
1386     //
1387     event->FillTrkVar(ca,ntrkentry);
1388     ntrkentry++;
1389     filled = true;
1390     //
1391     } else {
1392     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1393     }
1394     }
1395     //
1396     if ( !filled ) badevent++;
1397     }
1398 mocchiut 1.1 //
1399     // Clear structures used to communicate with fortran
1400     //
1401     event->ClearStructs();
1402     //
1403     // Fill the rootple
1404     //
1405     calo->Fill();
1406     //
1407 mocchiut 1.6 // delete ca;
1408 mocchiut 1.1 //
1409     jumpev:
1410 mocchiut 1.23 if ( !debug ) debug = false;
1411 mocchiut 1.1 //
1412     };
1413     //
1414     if ( verbose ) printf("\n SUMMARY:\n");
1415     if ( verbose ) printf(" Total number of events: %i \n",totevent);
1416     if ( verbose ) printf(" Events with at least one track: %i \n",totevent-badevent);
1417     if ( verbose ) printf(" Events without tracks: %i \n",badevent);
1418     //
1419     if ( badevent == totevent ){
1420 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - WARNING no tracks or good events in run %u \n",idRun);
1421 mocchiut 1.1 code = 101;
1422     };
1423     //
1424 mocchiut 1.8 delete dbtime;
1425 mocchiut 1.1 //
1426     // Clear variables before processing another run (needed to have always the same result when reprocessing data with the same software).
1427     //
1428     event->RunClose();
1429     if ( l0File ) l0File->Close();
1430     //
1431     }; // process all the runs
1432     //
1433     if ( verbose ) printf("\n Finished processing data \n");
1434     //
1435     closeandexit:
1436     //
1437     if ( !reprocall && reproc && code >= 0 ){
1438 mocchiut 1.59 //
1439 mocchiut 1.1 if ( totfileentries > noaftrun ){
1440     if ( verbose ) printf("\n Post-processing: copying events from the old tree after the processed run\n");
1441 mocchiut 1.8 if ( verbose ) printf(" Copying %u events in the file which are after the end of the run %u \n",(totfileentries-noaftrun),run);
1442     if ( verbose ) printf(" Start copying at event number %u end copying at event number %u \n",noaftrun,totfileentries);
1443 mocchiut 1.1 for (UInt_t j = noaftrun; j < totfileentries; j++ ){
1444     //
1445     // Get entry from old tree
1446     //
1447 mocchiut 1.48 if ( caloclone->GetEntry(j) <= 0 ) throw -36;
1448 mocchiut 1.1 //
1449     // copy caclone to ca
1450 mocchiut 1.59 //
1451     //EE if ( getl1 ) c1->Clear();
1452     //ca->Clear();
1453 mocchiut 1.1 //
1454 mocchiut 1.59 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
1455 mocchiut 1.58 memcpy(&ca,&caclone,sizeof(caclone));// EM reprocessing bug fixed
1456 mocchiut 1.59 if ( debug ) printf(" CA %f CACL %f \n",ca->qtot,caclone->qtot);
1457 mocchiut 1.1 //
1458     // Fill entry in the new tree
1459     //
1460 mocchiut 1.59 calo->Fill();// EM reprocessing bug fixed
1461 mocchiut 1.1 //
1462 mocchiut 1.59 }
1463 mocchiut 1.1 if ( verbose ) printf(" Finished successful copying!\n");
1464 mocchiut 1.59 }
1465     }
1466 mocchiut 1.1 //
1467     // Case of no errors: close files, delete old tree(s), write and close level2 file
1468     //
1469     if ( l0File ) l0File->Close();
1470     if ( tempfile ) tempfile->Close();
1471 mocchiut 1.30 if ( myfold ) gSystem->Unlink(tempname.str().c_str());
1472 mocchiut 1.1 //
1473     if ( code < 0 ) printf("\n CALORIMETER - ERROR: an error occurred, try to save anyway...\n");
1474     if ( verbose ) printf("\n Writing and closing rootple\n");
1475 mocchiut 1.53 if ( debug ) printf("\n Setting name\n");
1476 mocchiut 1.1 if ( calo ) calo->SetName("Calorimeter");
1477 mocchiut 1.53 if ( debug ) printf("\n Writing calo tree\n");
1478 mocchiut 1.1 if ( file ){
1479     file->cd();
1480 mocchiut 1.53 if ( calo ) calo->Write("Calorimeter", TObject::kOverwrite); // 10RED CRASH when exiting with no tracker, calo is not defined
1481 mocchiut 1.1 };
1482     //
1483 mocchiut 1.53 if ( debug ) printf("\n Unlink folder\n");
1484 mocchiut 1.30 if ( myfold ) gSystem->Unlink(calofolder.str().c_str());
1485 mocchiut 1.1 //
1486 mocchiut 1.51 // if ( ca ) delete ca;
1487     // if ( caclone ) delete caclone;
1488     // if ( c1 ) delete c1;
1489     // if ( c1clone ) delete c1clone;
1490     // if ( trk ) delete trk; // cannot delete anything since now all these objects are owned by TFile file!
1491    
1492 mocchiut 1.55 if ( tcNucleiTrk ){
1493     tcNucleiTrk->Delete();
1494     delete tcNucleiTrk;
1495     tcNucleiTrk = NULL;
1496     }
1497     if ( tcExtNucleiTrk ){
1498     tcExtNucleiTrk->Delete();
1499     delete tcExtNucleiTrk;
1500     tcExtNucleiTrk = NULL;
1501     }
1502     if ( tcExtTrk ){
1503     tcExtTrk->Delete();
1504     delete tcExtTrk;
1505     tcExtTrk = NULL;
1506     }
1507    
1508 mocchiut 1.53 if ( debug ) printf("\n Delete q4\n");
1509 mocchiut 1.7 if ( q4 ) delete q4;
1510 mocchiut 1.53 if ( debug ) printf("\n Delete glroot\n");
1511 mocchiut 1.7 if ( glroot ) delete glroot;
1512 mocchiut 1.53 if ( debug ) printf("\n Close runinfo\n");
1513 mocchiut 1.51 if ( runinfo ) runinfo->Close();
1514 mocchiut 1.53 if ( debug ) printf("\n Delete runinfo\n");
1515 mocchiut 1.7 if ( runinfo ) delete runinfo;
1516     //
1517 mocchiut 1.1 // the end
1518     //
1519     if ( verbose ) printf("\n Exiting...\n");
1520     if ( code < 0 ) throw code;
1521     return(code);
1522     }
1523    
1524    
1525    

  ViewVC Help
Powered by ViewVC 1.1.23