/[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.55 - (hide annotations) (download)
Fri Jun 6 14:18:14 2014 UTC (10 years, 10 months ago) by mocchiut
Branch: MAIN
Changes since 1.54: +202 -7 lines
New extended algorithm track-related added to CALO, TOF and ORB

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.11 tracker->SetBranchAddress("TrkLevel2",&trk);
546     nevtrkl2 = tracker->GetEntries();
547 mocchiut 1.55 //
548     // Look for extended tracking algorithm
549     //
550     if ( verbose ) printf("Look for extended and nuclei tracking algorithms\n");
551     // Nuclei tracking algorithm
552     Int_t checkAlgo = 0;
553     tcNucleiTrk = new TClonesArray("TrkTrack");
554     checkAlgo = tracker->SetBranchAddress("TrackNuclei",&tcNucleiTrk);
555     if ( !checkAlgo ){
556     if ( verbose ) printf(" Nuclei tracking algorithm branch found! :D \n");
557     hasNucleiTrk = true;
558     } else {
559     if ( verbose ) printf(" Nuclei tracking algorithm branch not found :( !\n");
560     printf(" ok, this is not a problem (it depends on tracker settings) \n");
561     delete tcNucleiTrk;
562     }
563     // Nuclei tracking algorithm using calorimeter points
564     tcExtNucleiTrk = new TClonesArray("ExtTrack");
565     checkAlgo = tracker->SetBranchAddress("RecoveredTrackNuclei",&tcExtNucleiTrk);
566     if ( !checkAlgo ){
567     if ( verbose ) printf(" Recovered nuclei tracking algorithm branch found! :D \n");
568     hasExtNucleiTrk = true;
569     } else {
570     if ( verbose ) printf(" Recovered nuclei tracking algorithm branch not found :( !\n");
571     printf(" ok, this is not a problem (it depends on tracker settings) \n");
572     delete tcExtNucleiTrk;
573     }
574     // Tracking algorithm using calorimeter points
575     tcExtTrk = new TClonesArray("ExtTrack");
576     checkAlgo = tracker->SetBranchAddress("RecoveredTrack",&tcExtTrk);
577     if ( !checkAlgo ){
578     if ( verbose ) printf(" Recovered track algorithm branch found! :D \n");
579     hasExtTrk = true;
580     } else {
581     if ( verbose ) printf(" Recovered track algorithm branch not found :( !\n");
582     printf(" ok, this is not a problem (it depends on tracker settings) \n");
583     delete tcExtTrk;
584     }
585     }
586 mocchiut 1.1 //
587 mocchiut 1.30 // Is it a file from simulations?
588     //
589     if ( checksimu ){
590     TTree *h20 = (TTree*)file->Get("h20");
591     if ( h20 ){
592     //
593     // yes!
594     //
595     crosst = false;
596 mocchiut 1.34 ctground = true;
597 mocchiut 1.30 mechal = true;
598 mocchiut 1.47 mask18 = true;
599     if ( verbose ) printf("\n\n SIMULATION!!! Setting --use-default-alig and --no-crosstalk --mask-plane18X flags!\n\n");
600 mocchiut 1.30 h20->Delete();
601     };
602     };
603     //
604 mocchiut 1.1 // Call runinfo
605     //
606     sgnl = 0;
607     runinfo = new ItoRunInfo(file);
608     //
609     // update versioning information and retrieve informations about the run to be processed
610     //
611     caloversion = CaloInfo(false);
612     sgnl = runinfo->Update(run,"CALO",caloversion);
613     //
614     if ( sgnl ){
615     if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
616     code = sgnl;
617     goto closeandexit;
618     } else {
619     sgnl = 0;
620     };
621     //
622     // number of events in the file BEFORE the first event of our run
623     //
624     nobefrun = runinfo->GetFirstEntry();
625     //
626     // total number of events in the file
627     //
628     totfileentries = runinfo->GetFileEntries();
629     //
630     // first file entry AFTER the last event of our run
631     //
632     noaftrun = runinfo->GetLastEntry() + 1;
633     //
634     // number of run to be processed
635     //
636     numbofrun = runinfo->GetNoRun();
637     //
638     // number of runs in the file
639     //
640     totnorun = runinfo->GetRunEntries();
641     //
642     // Does it contain already a Calorimeter branch? if so we are reprocessing data, if not we must create it.
643     //
644     caloclone = (TTree*)file->Get("Calorimeter");
645     //
646     if ( !caloclone ){
647     reproc = false;
648 mocchiut 1.8 if ( run == 0 && verbose ) printf(" CALORIMETER - WARNING: you are reprocessing data but calorimeter tree does not exist!\n");
649     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");
650 mocchiut 1.1 //
651     } else {
652     //
653 mocchiut 1.15 caloclone->SetAutoSave(900000000000000LL);
654 mocchiut 1.1 reproc = true;
655     //
656     if ( verbose ) printf("\n Preparing the pre-processing...\n");
657     //
658 mocchiut 1.28 if ( run == 0 || totnorun == 1 ){
659 mocchiut 1.1 //
660     // 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.
661     //
662     if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing all runs in the file\n");
663     reprocall = true;
664     //
665     } else {
666     //
667     // we are reprocessing a single run
668     //
669 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing run number %u \n",run);
670 mocchiut 1.1 reprocall = false;
671     //
672     //
673     //
674 mocchiut 1.30 gSystem->MakeDirectory(calofolder.str().c_str());
675     myfold = true;
676 mocchiut 1.1 tempfile = new TFile(tempname.str().c_str(),"RECREATE");
677     tempcalo = caloclone->CloneTree(-1,"fast");
678     tempcalo->SetName("Calorimeter-old");
679     tempfile->Write();
680     tempfile->Close();
681     };
682     //
683     // delete old tree
684     //
685 mocchiut 1.53 if ( debug ) printf(" Delete caloclone \n");
686     if ( caloclone ) caloclone->Delete("all");
687 mocchiut 1.1 //
688     if ( verbose ) printf("\n ...done!\n");
689     //
690     };
691     //
692     // create calorimeter tree calo
693     //
694     file->cd();
695     calo = new TTree("Calorimeter-new","PAMELA Level2 calorimeter data");
696 mocchiut 1.15 calo->SetAutoSave(900000000000000LL);
697 pam-fi 1.16 ca->Set();
698 mocchiut 1.53 if ( !cl1only ) calo->Branch("CaloLevel2","CaloLevel2",&ca);
699 mocchiut 1.13 if ( getl1 ) calo->Branch("CaloLevel1","CaloLevel1",&c1);
700 mocchiut 1.1 //
701 mocchiut 1.55 // create new branches for new tracking algorithms
702     //
703     if ( hasNucleiTrk ){
704     tcaloNucleiTrk = new TClonesArray("CaloTrkVar",1);
705     calo->Branch("TrackNuclei",&tcaloNucleiTrk);
706     }
707     if ( hasExtNucleiTrk ){
708     tcaloExtNucleiTrk = new TClonesArray("CaloTrkVar",1);
709     calo->Branch("RecoveredTrackNuclei",&tcaloExtNucleiTrk);
710     }
711     if ( hasExtTrk ){
712     tcaloExtTrk = new TClonesArray("CaloTrkVar",1);
713     calo->Branch("RecoveredTrack",&tcaloExtTrk);
714     }
715    
716     //
717 mocchiut 1.1 if ( reproc && !reprocall ){
718     //
719     //
720     //
721     tempfile = new TFile(tempname.str().c_str(),"READ");
722     caloclone = (TTree*)tempfile->Get("Calorimeter-old");
723 mocchiut 1.15 caloclone->SetAutoSave(900000000000000LL);
724 mocchiut 1.53 if ( !cl1only ) caloclone->SetBranchAddress("CaloLevel2",&caclone);
725 mocchiut 1.13 TBranch *havel1 = caloclone->GetBranch("CaloLevel1");
726     if ( !havel1 && getl1 ) throw -117;
727     if ( havel1 && !getl1 ) throw -118;
728     if ( getl1 ) caloclone->SetBranchAddress("CaloLevel1",&c1clone);
729 mocchiut 1.1 //
730 mocchiut 1.53 if ( nobefrun > 0 ){
731 mocchiut 1.1 if ( verbose ) printf("\n Pre-processing: copying events from the old tree before the processed run\n");
732 mocchiut 1.8 if ( verbose ) printf(" Copying %u events in the file which are before the beginning of the run %u \n",nobefrun,run);
733 mocchiut 1.1 if ( verbose ) printf(" Start copying at event number 0, end copying at event number %u \n",nobefrun);
734     for (UInt_t j = 0; j < nobefrun; j++){
735     //
736 mocchiut 1.48 if ( caloclone->GetEntry(j) <= 0 ) throw -36;
737 mocchiut 1.1 //
738     // copy caclone to ca
739     //
740 mocchiut 1.13 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
741 mocchiut 1.53 if ( !cl1only ) memcpy(&ca,&caclone,sizeof(caclone));
742 mocchiut 1.1 //
743     // Fill entry in the new tree
744     //
745 mocchiut 1.53 if ( getl1 || !cl1only ) calo->Fill();
746 mocchiut 1.1 //
747 mocchiut 1.13 if ( getl1 ) c1->Clear();
748 mocchiut 1.6 ca->Clear();
749     //
750 mocchiut 1.1 };
751     if ( verbose ) printf(" Finished successful copying!\n");
752 mocchiut 1.53 }
753     }
754 mocchiut 1.1 //
755     // Get the list of run to be processed
756     //
757     runlist = runinfo->GetRunList();
758     //
759     // Loop over the run to be processed
760     //
761     for (UInt_t irun=0; irun < numbofrun; irun++){
762     //
763     badevent = 0;
764     //
765     idRun = runlist->At(irun);
766     if ( verbose ) printf("\n\n\n ####################################################################### \n");
767 mocchiut 1.8 if ( verbose ) printf(" PROCESSING RUN NUMBER %u \n",idRun);
768 mocchiut 1.1 if ( verbose ) printf(" ####################################################################### \n\n\n");
769     //
770     sgnl = runinfo->GetRunInfo(idRun);
771     if ( sgnl ){
772     if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
773     code = sgnl;
774     goto closeandexit;
775     } else {
776     sgnl = 0;
777     };
778 mocchiut 1.8 id_reg_run = runinfo->ID_ROOT_L0;
779 mocchiut 1.1 runheadtime = runinfo->RUNHEADER_TIME;
780     runtrailtime = runinfo->RUNTRAILER_TIME;
781 mocchiut 1.8 evfrom = runinfo->EV_FROM;
782     evto = runinfo->EV_TO;
783 mocchiut 1.1 //
784     if ( id_reg_run == -1 ){
785 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - ERROR: no run with ID_RUN = %u \n\n Exiting... \n\n",idRun);
786 mocchiut 1.1 code = -5;
787     goto closeandexit;
788     };
789     //
790 mocchiut 1.8 // prepare the timesync for the db
791     //
792 mocchiut 1.32 TString host = glt->CGetHost();
793     TString user = glt->CGetUser();
794     TString psw = glt->CGetPsw();
795     TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
796     if ( !dbc->IsConnected() ) throw -116;
797 mocchiut 1.33 stringstream myquery;
798     myquery.str("");
799     myquery << "SET time_zone='+0:00'";
800 mocchiut 1.52 delete dbc->Query(myquery.str().c_str());
801 mocchiut 1.32 //
802 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
803 mocchiut 1.8 dbtime = new GL_TIMESYNC(runinfo->ID_ROOT_L0,"ID",dbc);
804     //
805 mocchiut 1.1 // Search in the DB the path and name of the LEVEL0 file to be processed.
806     //
807 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
808 mocchiut 1.8 glroot->Query_GL_ROOT(runinfo->ID_ROOT_L0,dbc);
809 mocchiut 1.1 //
810     ftmpname.str("");
811     ftmpname << glroot->PATH.Data() << "/";
812     ftmpname << glroot->NAME.Data();
813     fname = ftmpname.str().c_str();
814     //
815     // print out informations
816     //
817 mocchiut 1.8 totevent = runinfo->NEVENTS;
818 mocchiut 1.1 if ( verbose ) printf("\n LEVEL0 data file: %s \n",fname.Data());
819 mocchiut 1.8 if ( verbose ) printf(" RUN HEADER absolute time is: %u \n",runheadtime);
820     if ( verbose ) printf(" RUN TRAILER absolute time is: %u \n",runtrailtime);
821     if ( verbose ) printf(" %i events to be processed for run %u: from %i to %i (reg entries)\n\n",totevent,idRun,evfrom,evfrom+totevent);
822 mocchiut 1.1 //
823 mocchiut 1.43 // if ( !totevent ) goto closeandexit;
824 mocchiut 1.42 //
825 mocchiut 1.1 // Open Level0 file
826     //
827 mocchiut 1.52 if ( l0File ) l0File->Close();
828 mocchiut 1.1 l0File = new TFile(fname.Data());
829     if ( !l0File ) {
830     if ( verbose ) printf(" CALORIMETER - ERROR: problems opening Level0 file\n");
831     code = -6;
832     goto closeandexit;
833     };
834     l0tr = (TTree*)l0File->Get("Physics");
835     if ( !l0tr ) {
836     if ( verbose ) printf(" CALORIMETER - ERROR: no Physics tree in Level0 file\n");
837     l0File->Close();
838     code = -7;
839     goto closeandexit;
840     };
841     l0head = l0tr->GetBranch("Header");
842     if ( !l0head ) {
843     if ( verbose ) printf(" CALORIMETER - ERROR: no Header branch in Level0 tree\n");
844     l0File->Close();
845     code = -8;
846     goto closeandexit;
847     };
848     l0calo = l0tr->GetBranch("Calorimeter");
849     if ( !l0calo ) {
850     if ( verbose ) printf(" CALORIMETER - ERROR: no Calorimeter branch in Level0 tree\n");
851     l0File->Close();
852     code = -103;
853     goto closeandexit;
854     };
855     l0trig = l0tr->GetBranch("Trigger");
856     if ( !l0trig ) {
857     if ( verbose ) printf(" CALORIMETER - ERROR: no Trigger branch in Level0 tree\n");
858     l0File->Close();
859     code = -104;
860     goto closeandexit;
861     };
862     //
863     l0tr->SetBranchAddress("Trigger", &trig);
864 mocchiut 1.8 l0tr->SetBranchAddress("Header", &eh);
865 mocchiut 1.1 //
866 mocchiut 1.24 softinfo = (TTree*)l0File->Get("SoftInfo");
867     if ( softinfo ){
868     softinfo->SetBranchAddress("SoftInfo",&yodaver);
869 mocchiut 1.48 if ( softinfo->GetEntry(0) <= 0 ) throw -36;
870 mocchiut 1.24 };
871     if ( debug ) printf(" LEVEL0 FILE GENERATED WITH YODA VERSION %i \n",yodaver);
872     //
873 mocchiut 1.38 if ( withtrk ){
874     if ( trkpar1 || ( tttrkpar1 != 0 && tttrkpar1 < runheadtime ) ){
875     trkpar1 = false;
876     Int_t glpar = q4->Query_GL_PARAM(runinfo->RUNHEADER_TIME,1,dbc);
877     if ( glpar < 0 ){
878     code = glpar;
879     goto closeandexit;
880     };
881     tttrkpar1 = q4->TO_TIME;
882     // ----------------------------
883     // Read the magnetic field
884     // ----------------------------
885     if ( verbose ) printf(" Reading magnetic field maps at %s\n",(q4->PATH+q4->NAME).Data());
886     trk->LoadField(q4->PATH+q4->NAME);
887     if ( verbose ) printf("\n");
888     };
889     };
890     //
891     // Close DB connection
892     //
893     if ( dbc ){
894     dbc->Close();
895     delete dbc;
896 pam-fi 1.50 dbc = 0;
897 mocchiut 1.38 };
898     //
899     //
900 mocchiut 1.1 // Construct the event object, look for the calibration which include the first header
901     //
902     sgnl = 0;
903     if ( verbose ) printf(" Check for calorimeter calibrations and initialize event object \n");
904 mocchiut 1.31 // if ( !dbc->IsConnected() ) throw -116;
905 mocchiut 1.32 event->ProcessingInit(glt,runheadtime,sgnl,l0tr,debug,verbose);
906 mocchiut 1.47 if ( mask18 ){
907     event->MaskPlane18X();
908     } else {
909     event->UnMaskPlane18X();
910     };
911 mocchiut 1.1 if ( verbose ) printf("\n");
912     if ( sgnl == 100 ) {
913     code = sgnl;
914     if ( verbose ) printf(" CALORIMETER - WARNING: run header not included in any calibration interval\n");
915     sgnl = 0;
916     };
917     if ( sgnl ){
918     l0File->Close();
919     code = sgnl;
920     goto closeandexit;
921     };
922     //
923     qy.str("");
924     //
925 mocchiut 1.8 nevents = l0head->GetEntries();
926 mocchiut 1.1 caloevents = l0calo->GetEntries();
927     //
928 mocchiut 1.43 if ( nevents < 1 && totevent ) {
929 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: Level0 file is empty\n\n");
930     l0File->Close();
931     code = -11;
932     goto closeandexit;
933     };
934     //
935 mocchiut 1.43 if ( evto > nevents-1 && totevent ) {
936 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: too few entries in the registry tree\n");
937     l0File->Close();
938     code = -12;
939     goto closeandexit;
940     };
941     //
942     // Check if we have to load parameter files
943     //
944     sgnl = 0;
945 mocchiut 1.32 sgnl = event->ChkParam(glt,runheadtime,mechal); // calorimeter parameter files
946 mocchiut 1.1 if ( sgnl < 0 ){
947     code = sgnl;
948     l0File->Close();
949     goto closeandexit;
950     };
951     //
952 mocchiut 1.34 // Calculate the cross talk corrections needed for this run using calibration informations
953     //
954     if ( crosst && !ctground ){
955     sgnl = 0;
956 mocchiut 1.41 sgnl = event->CalcCrossTalkCorr(glt,runheadtime,usetable);
957 mocchiut 1.34 if ( sgnl < 0 ){
958     code = sgnl;
959     l0File->Close();
960     goto closeandexit;
961     };
962     };
963     //
964 mocchiut 1.1 // run over all the events of the run
965     //
966     if ( verbose ) printf("\n Ready to start! \n\n Processed events: \n\n");
967     //
968 mocchiut 1.19 for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+runinfo->NEVENTS); re++){
969 mocchiut 1.44 // for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+5000); re++){
970 mocchiut 1.49 //for ( re = 92012+3400; re < 92012+3700; re++){
971     //for ( re = 4451+2833; re < 4451+2836 ; re++){
972 mocchiut 1.18 // printf(" i = %i \n",re-200241);
973 mocchiut 1.1 //
974     if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000);
975     //
976 mocchiut 1.39 if ( debug ) printf("\n\n\n EVENT NUMBER %i \n",procev);
977     //
978 mocchiut 1.48 if ( l0head->GetEntry(re) <= 0 ) throw -36;
979 mocchiut 1.1 //
980     // absolute time of this event
981     //
982 mocchiut 1.8 ph = eh->GetPscuHeader();
983     atime = dbtime->DBabsTime(ph->GetOrbitalTime());
984 mocchiut 1.1 //
985     //
986     //
987 mocchiut 1.8 if ( re > caloevents-1 ){
988 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - ERROR: no physics events with entry = %i in Level0 file\n",i);
989     l0File->Close();
990     code = -112;
991     goto closeandexit;
992     };
993     //
994 mocchiut 1.45 if ( atime > (runtrailtime+1) || atime < (runheadtime-1) ) {
995 mocchiut 1.1 if ( verbose ) printf(" CALORIMETER - WARNING: event at time outside the run time window, skipping it\n");
996 mocchiut 1.9 jumped++;
997 mocchiut 1.1 goto jumpev;
998     };
999     //
1000     // retrieve tracker informations
1001     //
1002     if ( !reprocall ){
1003 mocchiut 1.19 itr = nobefrun + (re - evfrom - jumped);
1004     //itr = re-(46438+200241);
1005 mocchiut 1.1 } else {
1006 mocchiut 1.9 itr = runinfo->GetFirstEntry() + (re - evfrom - jumped);
1007 mocchiut 1.1 };
1008 mocchiut 1.11 //
1009     if ( withtrk ){
1010     if ( itr > nevtrkl2 ){
1011     if ( verbose ) printf(" CALORIMETER - ERROR: no tracker events with entry = %i in Level2 file\n",itr);
1012     l0File->Close();
1013     code = -113;
1014     goto closeandexit;
1015 mocchiut 1.55 }
1016 mocchiut 1.11 //
1017     trk->Clear();
1018     //
1019 mocchiut 1.55 //
1020     // Clones array must be cleared before going on
1021     //
1022     if ( hasNucleiTrk ){
1023     tcNucleiTrk->Delete();
1024     tcaloNucleiTrk->Delete();
1025     }
1026     if ( hasExtNucleiTrk ){
1027     tcExtNucleiTrk->Delete();
1028     tcaloExtNucleiTrk->Delete();
1029     }
1030     if ( hasExtTrk ){
1031     tcExtTrk->Delete();
1032     tcaloExtTrk->Delete();
1033     }
1034     //
1035 mocchiut 1.48 if ( tracker->GetEntry(itr) <= 0 ) throw -36;
1036 mocchiut 1.55 //
1037     }
1038 mocchiut 1.10 //
1039 mocchiut 1.1 procev++;
1040     //
1041     // start processing
1042     //
1043 mocchiut 1.13 if ( getl1 ) c1->Clear();
1044 mocchiut 1.6 ca->Clear();
1045 mocchiut 1.1 //
1046     // determine who generate the trigger for this event (TOF, S4/PULSER, CALO)
1047     //
1048 mocchiut 1.48 if ( l0trig->GetEntry(re) <= 0 ) throw -36;
1049 mocchiut 1.47 // S3 = 0;
1050     // S2 = 0;
1051     // S12 = 0;
1052     // S11 = 0;
1053     // S3 = trig->patterntrig[2];
1054     // S2 = trig->patterntrig[3];
1055     // S12 = trig->patterntrig[4];
1056     // S11 = trig->patterntrig[5];
1057     // if ( trig->patterntrig[1] & (1<<0) ) tmptrigty = 1.;
1058     // if ( trig->patterntrig[0] ) tmptrigty = 2.;
1059     // if ( S3 || S2 || S12 || S11 ) tmptrigty = 0.;
1060     // if ( !(trig->patterntrig[1] & (1<<0)) && !trig->patterntrig[0] && !S3 && !S2 && !S12 && !S11 ) tmptrigty = 1.;
1061     // event->clevel2->trigty = tmptrigty;
1062     event->clevel2->trigty = (Float_t)IsCaloSelfTrigger(trig->trigconf,(UInt_t*)trig->patterntrig);
1063     if ( debug ) printf(" starting with trigty %f \n",event->clevel2->trigty);
1064 mocchiut 1.1 //
1065     // check if the calibration we are using is still good, if not load another calibration
1066     //
1067     sgnl = 0;
1068 mocchiut 1.32 sgnl = event->ChkCalib(glt,atime);
1069 mocchiut 1.1 if ( sgnl < 0 ){
1070     code = sgnl;
1071     goto closeandexit;
1072     };
1073     if ( sgnl == 100 ){
1074     code = sgnl;
1075     if ( verbose ) printf(" CALORIMETER - WARNING: data not associated to any calibration interval\n");
1076     badevent++;
1077     sgnl = 0;
1078     };
1079     //
1080 mocchiut 1.11 // do we have at least one track from the tracker? this check has been disabled
1081 mocchiut 1.1 //
1082     event->clevel1->good2 = 1;
1083 mocchiut 1.8 //
1084 mocchiut 1.29 // use the whole calorimeter, 22 W planes
1085     //
1086     event->clevel1->npla = 22;
1087     //
1088     // Use the standard silicon planes shifting
1089     //
1090     event->clevel1->reverse = 0;
1091     //
1092 mocchiut 1.1 //
1093 mocchiut 1.8 // Calibrate calorimeter event "re" and store output in the two structures that will be passed to fortran routine
1094 mocchiut 1.1 //
1095 mocchiut 1.11 event->Calibrate(re);
1096 mocchiut 1.1 //
1097     // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract topological variables.
1098     //
1099     //
1100     // Calculate variables common to all tracks (qtot, nstrip, etc.)
1101     //
1102 mocchiut 1.53 if ( !cl1only ) event->GetCommonVar();
1103 mocchiut 1.1 //
1104     // Fill common variables
1105     //
1106 mocchiut 1.53 if ( cl1only ){
1107     event->FillCommonVar(c1);
1108     } else {
1109     event->FillCommonVar(c1,ca);
1110     }
1111 mocchiut 1.1 //
1112     // Calculate variables related to tracks only if we have at least one track (from selftrigger and/or tracker)
1113     //
1114     ntrkentry = 0;
1115     //
1116     filled = false;
1117 mocchiut 1.11 //
1118 mocchiut 1.53 // skip track-related if processing level1 only
1119 mocchiut 1.1 //
1120 mocchiut 1.53 if ( !cl1only ){
1121     //
1122     // Run over tracks (tracker or calorimeter )
1123     //
1124     if ( withtrk ){
1125     for(Int_t nt=0; nt < trk->ntrk(); nt++){
1126     //
1127     event->clevel1->good2 = 1;
1128     //
1129     TrkTrack *ptt = trk->GetStoredTrack(nt);
1130     //
1131     event->clevel1->trkchi2 = 0;
1132     //
1133     // Copy the alpha vector in the input structure
1134     //
1135     for (Int_t e = 0; e < 5 ; e++){
1136     event->clevel1->al_p[e][0] = ptt->al[e];
1137     };
1138     //
1139     // Get tracker related variables for this track
1140     //
1141     event->GetTrkVar();
1142     //
1143     // Save tracker track sequence number
1144     //
1145     event->trkseqno = nt;
1146     //
1147     // Copy values in the class ca from the structure clevel2
1148     //
1149     event->FillTrkVar(ca,ntrkentry);
1150     ntrkentry++;
1151     filled = true;
1152     //
1153 mocchiut 1.55 } // loop on all the tracks
1154     //
1155     // Code for extended tracking algorithm:
1156     //
1157     if ( hasNucleiTrk ){
1158     Int_t ttentry = 0;
1159     for(Int_t nt=0; nt < tcNucleiTrk->GetEntries() ; nt++){
1160     //
1161     event->clevel1->good2 = 1;
1162     //
1163     TrkTrack *ptt = (TrkTrack*)(tcNucleiTrk->At(nt));
1164     //
1165     event->clevel1->trkchi2 = 0;
1166     //
1167     // Copy the alpha vector in the input structure
1168     //
1169     for (Int_t e = 0; e < 5 ; e++){
1170     event->clevel1->al_p[e][0] = ptt->al[e];
1171     };
1172     //
1173     // Get tracker related variables for this track
1174     //
1175     event->GetTrkVar();
1176     //
1177     // Save tracker track sequence number
1178     //
1179     event->trkseqno = nt;
1180     //
1181     // Copy values in the class ca from the structure clevel2
1182     //
1183     event->FillTrkVar(tcaloNucleiTrk,ttentry);
1184     ttentry++;
1185     //
1186     } // loop on all the tracks
1187     }
1188    
1189     if ( hasExtNucleiTrk ){
1190     Int_t ttentry = 0;
1191     for(Int_t nt=0; nt < tcExtNucleiTrk->GetEntries() ; nt++){
1192     //
1193     event->clevel1->good2 = 1;
1194     //
1195     ExtTrack *ptt = (ExtTrack*)(tcExtNucleiTrk->At(nt));
1196     //
1197     event->clevel1->trkchi2 = 0;
1198     //
1199     // Copy the alpha vector in the input structure
1200     //
1201     for (Int_t e = 0; e < 5 ; e++){
1202     event->clevel1->al_p[e][0] = ptt->al[e];
1203     };
1204     //
1205     // Get tracker related variables for this track
1206     //
1207     event->GetTrkVar();
1208     //
1209     // Save tracker track sequence number
1210     //
1211     event->trkseqno = nt;
1212     //
1213     // Copy values in the class ca from the structure clevel2
1214     //
1215     event->FillTrkVar(tcaloExtNucleiTrk,ttentry);
1216     ttentry++;
1217     //
1218     } // loop on all the tracks
1219     }
1220    
1221     if ( hasExtTrk ){
1222     Int_t ttentry = 0;
1223     for(Int_t nt=0; nt < tcExtTrk->GetEntries() ; nt++){
1224     //
1225     event->clevel1->good2 = 1;
1226     //
1227     ExtTrack *ptt = (ExtTrack*)(tcExtTrk->At(nt));
1228     //
1229     event->clevel1->trkchi2 = 0;
1230     //
1231     // Copy the alpha vector in the input structure
1232     //
1233     for (Int_t e = 0; e < 5 ; e++){
1234     event->clevel1->al_p[e][0] = ptt->al[e];
1235     };
1236     //
1237     // Get tracker related variables for this track
1238     //
1239     event->GetTrkVar();
1240     //
1241     // Save tracker track sequence number
1242     //
1243     event->trkseqno = nt;
1244     //
1245     // Copy values in the class ca from the structure clevel2
1246     //
1247     event->FillTrkVar(tcaloExtTrk,ttentry);
1248     ttentry++;
1249     //
1250     } // loop on all the tracks
1251     }
1252    
1253     }
1254 mocchiut 1.53 //
1255     // 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
1256     // fit of the track (to be used for example when TRK is off due to any reason like IPM3/5 off).
1257     // here we make an event selection so it must be done very carefully...
1258     //
1259     // 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
1260     //
1261     if ( trackanyway && !filled && event->clevel2->npcfit[0] >= 2 && event->clevel2->npcfit[1] >= 2 && event->clevel2->good != 0 && event->clevel2->trigty < 2. ){
1262     if ( debug ) printf(" Event with a track not fitted by the tracker at entry %i \n",itr);
1263     //
1264     // Disable "track mode" in the fortran routine
1265     //
1266     event->clevel1->good2 = 0;
1267     event->clevel1->riginput = rigdefault;
1268     if ( debug ) printf(" Using as default rigidity: %f \n",event->clevel1->riginput);
1269     //
1270     // We have a selftrigger event to analyze.
1271     //
1272     for (Int_t e = 0; e < 5 ; e++){
1273     event->clevel1->al_p[e][0] = 0.;
1274     event->clevel1->al_p[e][1] = 0.;
1275     }
1276     event->clevel1->trkchi2 = 0;
1277     //
1278     event->GetTrkVar();
1279     //
1280     // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
1281     //
1282     if ( event->clevel1->good2 == 0 ) {
1283     //
1284     // In selftrigger mode the trkentry variable is set to -1
1285     //
1286     event->trkseqno = -3;
1287     //
1288     // Copy values in the class ca from the structure clevel2
1289     //
1290     if ( debug ) printf(" -3 begin \n");
1291     event->FillTrkVar(ca,ntrkentry);
1292     if ( debug ) printf(" -3 end \n");
1293     ntrkentry++;
1294     filled = true;
1295     //
1296     } else {
1297     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1298     }
1299     //
1300     }
1301     //
1302     // Call high energy nuclei routine
1303     //
1304     if ( hZn && event->clevel2->trigty >= 2. ){
1305     if ( debug ) printf(" Calling selftrigger high energy nuclei routine for entry %i \n",itr);
1306     //
1307     // Disable "track mode" in the fortran routine
1308     //
1309     event->clevel1->good2 = 0;
1310     //
1311     // Set high energy nuclei flag to one
1312     //
1313     event->clevel1->hzn = 1;
1314     event->clevel1->riginput = rigdefault;
1315     //
1316     // We have a selftrigger event to analyze.
1317     //
1318     for (Int_t e = 0; e < 5 ; e++){
1319     event->clevel1->al_p[e][0] = 0.;
1320     event->clevel1->al_p[e][1] = 0.;
1321     }
1322     event->clevel1->trkchi2 = 0;
1323     //
1324     event->GetTrkVar();
1325     //
1326     // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
1327     //
1328     if ( event->clevel1->good2 == 0 ) {
1329     //
1330     // In selftrigger mode the trkentry variable is set to -1
1331     //
1332     event->trkseqno = -2;
1333     //
1334     // Copy values in the class ca from the structure clevel2
1335     //
1336     event->FillTrkVar(ca,ntrkentry);
1337     ntrkentry++;
1338     filled = true;
1339     //
1340     } else {
1341     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1342     }
1343     //
1344     }
1345     //
1346     // self trigger event
1347     //
1348     if ( st && event->clevel2->trigty >= 2. ){
1349     if ( verbose ) printf(" Selftrigger event at entry %i \n",itr);
1350     //
1351     // Disable "track mode" in the fortran routine
1352     //
1353     event->clevel1->good2 = 0;
1354     //
1355     // disable high enery nuclei flag;
1356     //
1357     event->clevel1->hzn = 0;
1358     //
1359     // We have a selftrigger event to analyze.
1360     //
1361     for (Int_t e = 0; e < 5 ; e++){
1362     event->clevel1->al_p[e][0] = 0.;
1363     event->clevel1->al_p[e][1] = 0.;
1364     };
1365     event->clevel1->trkchi2 = 0;
1366     //
1367     event->GetTrkVar();
1368     //
1369     // if we had no problem (clevel2->good = 0, NOTICE zero, not one in selftrigger mode!), fill and go on
1370     //
1371     if ( event->clevel1->good2 == 0 ) {
1372     //
1373     // In selftrigger mode the trkentry variable is set to -1
1374     //
1375     event->trkseqno = -1;
1376     //
1377     // Copy values in the class ca from the structure clevel2
1378     //
1379     event->FillTrkVar(ca,ntrkentry);
1380     ntrkentry++;
1381     filled = true;
1382     //
1383     } else {
1384     if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1385     }
1386     }
1387     //
1388     if ( !filled ) badevent++;
1389     }
1390 mocchiut 1.1 //
1391     // Clear structures used to communicate with fortran
1392     //
1393     event->ClearStructs();
1394     //
1395     // Fill the rootple
1396     //
1397     calo->Fill();
1398     //
1399 mocchiut 1.6 // delete ca;
1400 mocchiut 1.1 //
1401     jumpev:
1402 mocchiut 1.23 if ( !debug ) debug = false;
1403 mocchiut 1.1 //
1404     };
1405     //
1406     if ( verbose ) printf("\n SUMMARY:\n");
1407     if ( verbose ) printf(" Total number of events: %i \n",totevent);
1408     if ( verbose ) printf(" Events with at least one track: %i \n",totevent-badevent);
1409     if ( verbose ) printf(" Events without tracks: %i \n",badevent);
1410     //
1411     if ( badevent == totevent ){
1412 mocchiut 1.8 if ( verbose ) printf("\n CALORIMETER - WARNING no tracks or good events in run %u \n",idRun);
1413 mocchiut 1.1 code = 101;
1414     };
1415     //
1416 mocchiut 1.8 delete dbtime;
1417 mocchiut 1.1 //
1418     // Clear variables before processing another run (needed to have always the same result when reprocessing data with the same software).
1419     //
1420     event->RunClose();
1421     if ( l0File ) l0File->Close();
1422     //
1423     }; // process all the runs
1424     //
1425     if ( verbose ) printf("\n Finished processing data \n");
1426     //
1427     closeandexit:
1428     //
1429     if ( !reprocall && reproc && code >= 0 ){
1430     if ( totfileentries > noaftrun ){
1431     if ( verbose ) printf("\n Post-processing: copying events from the old tree after the processed run\n");
1432 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);
1433     if ( verbose ) printf(" Start copying at event number %u end copying at event number %u \n",noaftrun,totfileentries);
1434 mocchiut 1.1 for (UInt_t j = noaftrun; j < totfileentries; j++ ){
1435     //
1436     // Get entry from old tree
1437     //
1438 mocchiut 1.48 if ( caloclone->GetEntry(j) <= 0 ) throw -36;
1439 mocchiut 1.1 //
1440     // copy caclone to ca
1441     //
1442 mocchiut 1.13 if ( getl1 ) c1->Clear();
1443 mocchiut 1.6 ca->Clear();
1444 mocchiut 1.8 //
1445 mocchiut 1.13 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
1446 mocchiut 1.54 if ( !cl1only ) memcpy(&ca,&caclone,sizeof(caclone));
1447 mocchiut 1.1 //
1448     // Fill entry in the new tree
1449     //
1450 mocchiut 1.53 if ( getl1 || !cl1only ) calo->Fill();
1451 mocchiut 1.1 //
1452     };
1453     if ( verbose ) printf(" Finished successful copying!\n");
1454     };
1455     };
1456     //
1457     // Case of no errors: close files, delete old tree(s), write and close level2 file
1458     //
1459     if ( l0File ) l0File->Close();
1460     if ( tempfile ) tempfile->Close();
1461 mocchiut 1.30 if ( myfold ) gSystem->Unlink(tempname.str().c_str());
1462 mocchiut 1.1 //
1463     if ( code < 0 ) printf("\n CALORIMETER - ERROR: an error occurred, try to save anyway...\n");
1464     if ( verbose ) printf("\n Writing and closing rootple\n");
1465 mocchiut 1.53 if ( debug ) printf("\n Setting name\n");
1466 mocchiut 1.1 if ( calo ) calo->SetName("Calorimeter");
1467 mocchiut 1.53 if ( debug ) printf("\n Writing calo tree\n");
1468 mocchiut 1.1 if ( file ){
1469     file->cd();
1470 mocchiut 1.53 if ( calo ) calo->Write("Calorimeter", TObject::kOverwrite); // 10RED CRASH when exiting with no tracker, calo is not defined
1471 mocchiut 1.1 };
1472     //
1473 mocchiut 1.53 if ( debug ) printf("\n Unlink folder\n");
1474 mocchiut 1.30 if ( myfold ) gSystem->Unlink(calofolder.str().c_str());
1475 mocchiut 1.1 //
1476 mocchiut 1.51 // if ( ca ) delete ca;
1477     // if ( caclone ) delete caclone;
1478     // if ( c1 ) delete c1;
1479     // if ( c1clone ) delete c1clone;
1480     // if ( trk ) delete trk; // cannot delete anything since now all these objects are owned by TFile file!
1481    
1482 mocchiut 1.55 if ( tcNucleiTrk ){
1483     tcNucleiTrk->Delete();
1484     delete tcNucleiTrk;
1485     tcNucleiTrk = NULL;
1486     }
1487     if ( tcExtNucleiTrk ){
1488     tcExtNucleiTrk->Delete();
1489     delete tcExtNucleiTrk;
1490     tcExtNucleiTrk = NULL;
1491     }
1492     if ( tcExtTrk ){
1493     tcExtTrk->Delete();
1494     delete tcExtTrk;
1495     tcExtTrk = NULL;
1496     }
1497    
1498 mocchiut 1.53 if ( debug ) printf("\n Delete q4\n");
1499 mocchiut 1.7 if ( q4 ) delete q4;
1500 mocchiut 1.53 if ( debug ) printf("\n Delete glroot\n");
1501 mocchiut 1.7 if ( glroot ) delete glroot;
1502 mocchiut 1.53 if ( debug ) printf("\n Close runinfo\n");
1503 mocchiut 1.51 if ( runinfo ) runinfo->Close();
1504 mocchiut 1.53 if ( debug ) printf("\n Delete runinfo\n");
1505 mocchiut 1.7 if ( runinfo ) delete runinfo;
1506     //
1507 mocchiut 1.1 // the end
1508     //
1509     if ( verbose ) printf("\n Exiting...\n");
1510     if ( code < 0 ) throw code;
1511     return(code);
1512     }
1513    
1514    
1515    

  ViewVC Help
Powered by ViewVC 1.1.23