/[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.56 - (hide annotations) (download)
Mon Jul 21 13:01:49 2014 UTC (10 years, 7 months ago) by mocchiut
Branch: MAIN
Changes since 1.55: +4 -0 lines
TTree replicas in high NEVENTS run bug (hopefully) fixed

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

  ViewVC Help
Powered by ViewVC 1.1.23