/[PAMELA software]/calo/flight/FCaloREPROC/src/CaloCore.cpp
ViewVC logotype

Contents of /calo/flight/FCaloREPROC/src/CaloCore.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1.1.1 - (show annotations) (download) (vendor branch)
Fri Oct 26 11:41:32 2007 UTC (17 years, 1 month ago) by mocchiut
Branch: MAIN, FCaloREPROC
CVS Tags: start, v0r00, HEAD
Changes since 1.1: +0 -0 lines
Error occurred while calculating annotation data.
Imported sources

1 //
2 // Given a calibration and a data file this program create an ntuple with LEVEL2 calorimeter variables - Emiliano Mocchiutti
3 //
4 // CaloCore.cxx version 3.05 (2006-05-30)
5 //
6 // The only input needed is the path to the directory created by YODA for the data file you want to analyze.
7 //
8 // Changelog:
9 //
10 // 3.08 (2006-11-13): Added high energy nuclei capability and "process all events" capability.
11 //
12 // 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.
13 //
14 // 3.03 - 3.04 (2006-05-23): Forgot to put impx and impy in the PAMELA reference system, fixed.
15 //
16 // 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.
17 //
18 // 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.
19 // Added variable planemax[2], plane of maximum energy release (x and y) in final output. Use ItoRunInfo instead of RunInfo.
20 //
21 // 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
22 // small bug in fortran routines
23 //
24 // 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
25 // working class "CaloProcessing", linked to the preliminary tracker flight software v0r00, reads YODA unique output files,
26 // reduced the number of installed libraries, F77 programs splitted depending on the function contained, introduced the readout and
27 // processing of self-trigger events, if the tracker provides more than one track all calorimeter track-related variables are saved
28 // as many times as the number of tracks (via the TClonesArray object in the level2 rootple) and many other small changes.
29 //
30 // 2.00 - 2.01 (2006-01-26): bug: wrong calculation of baselines in some cases, fixed.
31 //
32 // 1.00 - 2.00 (2006-01-11): use TSQL ROOT classes instead of directly calling MySQL.
33 //
34 // 0.00 - 1.00 (2005-09-14): seems working.
35 //
36 // 0.00 (2005-09-09): clone of CaloLEVEL2.c .
37 //
38 // C/C++ headers
39 //
40 #include <fstream>
41 #include <string.h>
42 //
43 // ROOT headers
44 //
45 #include <TTree.h>
46 #include <TClassEdit.h>
47 #include <TObject.h>
48 #include <TList.h>
49 #include <TArrayI.h>
50 #include <TSystem.h>
51 #include <TSystemDirectory.h>
52 #include <TString.h>
53 #include <TFile.h>
54 #include <TClass.h>
55 #include <TCanvas.h>
56 #include <TH1.h>
57 #include <TH1F.h>
58 #include <TH2D.h>
59 #include <TLatex.h>
60 #include <TPad.h>
61 #include <TSQLServer.h>
62 #include <TSQLRow.h>
63 #include <TSQLResult.h>
64 #include <TClonesArray.h>
65 #include <TStreamerInfo.h>
66 //
67 // This program headers
68 //
69 #include <RunInfo.h>
70 //
71 // YODA headers
72 //
73 #include <PamelaRun.h>
74 #include <physics/trigger/TriggerEvent.h>
75 //
76 // This program headers
77 //
78 #include <CaloCore.h>
79 #include <CaloLevel1.h>
80 #include <CaloLevel2.h>
81 #include <CaloLevel0.h>
82 #include <CaloVerl2.h>
83 //
84 // Tracker classes headers and definitions
85 //
86 #include <TrkLevel2.h>
87 //
88 using namespace std;
89 //
90 // CORE ROUTINE
91 //
92 int CaloCore(UInt_t run, TFile *file, GL_TABLES *glt, Int_t calargc, char *calargv[]){
93 //
94 // Set these to true to have a very verbose output.
95 //
96 Bool_t verbose = false;
97 Bool_t debug = false;
98 //
99 Bool_t crosst = true;
100 Bool_t ctground = true;
101 //
102 Bool_t trackanyway = true;
103 //
104 Float_t rigdefault = 50.;
105 //
106 Bool_t hZn = true;
107 //
108 Bool_t withtrk = true;
109 //
110 Bool_t st = true;
111 //
112 Bool_t getl1 = true;
113 //
114 Bool_t mechal = false;
115 //
116 Bool_t checksimu = true;
117 //
118 // Output directory is the working directoy.
119 //
120 const char* outdir = gSystem->DirName(gSystem->DirName(file->GetPath()));
121 //
122 Int_t ri = 0;
123 TString processFolder = Form("calorimeterFolder_%u",run);
124 if ( calargc > 0 ){
125 ri = 0;
126 while ( ri < calargc ){
127 if ( !strcmp(calargv[ri],"-processFolder") ) {
128 if ( calargc < ri+1 ){
129 throw -3;
130 };
131 processFolder = (TString)calargv[ri+1];
132 ri++;
133 };
134 if ( !strcmp(calargv[ri],"-v") || !strcmp(calargv[ri],"--verbose") ) {
135 verbose = true;
136 };
137 if ( !strcmp(calargv[ri],"-g") || !strcmp(calargv[ri],"--debug") ) {
138 verbose = true;
139 debug = true;
140 };
141 if ( !strcmp(calargv[ri],"--alltracks") ) {
142 trackanyway = true;
143 };
144 if ( !strcmp(calargv[ri],"--use-default-alig") ) {
145 mechal = true;
146 };
147 if ( !strcmp(calargv[ri],"--no-crosstalk") ) {
148 crosst = false;
149 };
150 if ( !strcmp(calargv[ri],"--flight-crosstalk") ) {
151 ctground = false;
152 };
153 if ( !strcmp(calargv[ri],"--ground-crosstalk") ) {
154 ctground = true;
155 };
156 if ( !strcmp(calargv[ri],"--no-tracker") ) {
157 withtrk = false;
158 };
159 if ( !strcmp(calargv[ri],"--with-tracker") ) {
160 withtrk = true;
161 };
162 if ( !strcmp(calargv[ri],"--defrig") ) {
163 if ( calargc < ri+1 ){
164 throw -3;
165 };
166 rigdefault = atof(calargv[ri+1]);
167 ri++;
168 };
169 if ( !strcmp(calargv[ri],"--no-alltracks") ) {
170 trackanyway = false;
171 };
172 if ( !strcmp(calargv[ri],"--highZnuclei") ) {
173 hZn = true;
174 };
175 if ( !strcmp(calargv[ri],"--no-highZnuclei") ) {
176 hZn = false;
177 };
178 if ( !strcmp(calargv[ri],"--selftrigger") ) {
179 st = true;
180 };
181 if ( !strcmp(calargv[ri],"--no-selftrigger") ) {
182 st = false;
183 };
184 if ( !strcmp(calargv[ri],"--no-level1") ) {
185 getl1 = false;
186 };
187 if ( !strcmp(calargv[ri],"--ignore-h20") ) {
188 checksimu = false;
189 };
190 if ( !strcmp(calargv[ri],"--help") ) {
191 printf("\n\n CALORIMETER HELP CALLED\n\n");
192 printf(" CaloCore options: \n");
193 printf(" -v | --verbose be verbose\n");
194 printf(" -g | --debug be really verbose\n");
195 printf(" --defrig rig rig is the default rigidity in GV to be used to\n");
196 printf(" obtain calorimeter variables in the routines\n");
197 printf(" \"alltracks\" and \"higZnuclei\" [default = 50]\n");
198 printf(" --alltracks fill the track related variables even in the case\n");
199 printf(" of no tracks from tracker and no selftrigger event\n");
200 printf(" when we have a calorimeter fit for both views [default]\n");
201 printf(" --no-alltracks fill the track related variables only in the case\n");
202 printf(" of a good track from tracker or selftrigger\n");
203 printf(" --highZnuclei call the routine to analyze high Z nuclei\n");
204 printf(" selftrigger events [default]\n");
205 printf(" --no-highZnuclei do not call the routine to analyze high Z nuclei\n");
206 printf(" selftrigger events\n");
207 printf(" --no-tracker do not use tracker level2\n");
208 printf(" --with-tracker use tracker level2 [default]\n");
209 printf(" --no-crosstalk do not apply crosstalk corrections\n");
210 printf(" --ground-crosstalk apply ground crosstalk corrections [default]\n");
211 printf(" --flight-crosstalk apply flight crosstalk corrections\n");
212 printf(" --selftrigger process selftrigger events [default]\n");
213 printf(" --no-selftrigger skip selftrigger events\n");
214 printf(" --no-level1 do not save Level1 TBranch\n");
215 printf(" --ignore-h20 do not check if it is a file from simulations\n");
216 printf(" --use-default-alig use default mechanical alignement (defined in CaloLevel1.h)\n");
217 throw -114;
218 };
219 ri++;
220 };
221 };
222 //
223 if ( verbose ){
224 printf("\n");
225 if ( getl1 ) printf(" Saving calorimeter level1 data \n");
226 if ( st ) printf(" Calling selftrigger subroutine \n");
227 if ( hZn ) printf(" Calling high energy nuclei subroutine \n");
228 if ( trackanyway ) printf(" Filling track related variables for all the possible tracks \n");
229 if ( hZn || trackanyway ) printf(" => default assumed rigidity %f \n",rigdefault);
230 if ( withtrk ) printf(" Using tracker level2 data \n");
231 if ( mechal ) printf(" Using default mechanical alignement (defined in CaloLevel1.h)\n");
232 if ( crosst ){
233 printf(" Applying cross-talk corrections \n");
234 if ( ctground ){
235 printf(" => Using ground cross-talk coefficients \n");
236 } else {
237 printf(" => Using flight cross-talk coefficients \n");
238 };
239 } else {
240 printf(" Do not applying cross-talk corrections \n");
241 };
242 if ( !checksimu ) printf(" Check on h20 TTree in level2 file skipped \n");
243 printf("\n");
244 };
245 //
246 // Working filename
247 //
248 TString outputfile;
249 stringstream name;
250 name.str("");
251 name << outdir << "/";
252 //
253 // Variables.
254 //
255 TTree *tracker = 0;
256 TTree *calo = 0;
257 TTree *caloclone = 0;
258 Bool_t reproc = false;
259 Bool_t reprocall = false;
260 UInt_t nevents = 0;
261 UInt_t nobefrun = 0;
262 UInt_t noaftrun = 0;
263 UInt_t numbofrun = 0;
264 UInt_t totnorun = 0;
265 //
266 Int_t code = 0;
267 Int_t sgnl;
268 //
269 // calorimeter level2 classes
270 //
271 CaloLevel1 *c1 = 0;
272 CaloLevel1 *c1clone = 0;
273 if ( getl1 ){
274 c1 = new CaloLevel1();
275 c1clone = new CaloLevel1();
276 };
277 //
278 // calorimeter level2 classes
279 //
280 CaloLevel2 *ca = new CaloLevel2();
281 CaloLevel2 *caclone = new CaloLevel2();
282 //
283 TrkLevel2 *trk = new TrkLevel2();
284 Int_t nevtrkl2 = 0;
285 //
286 UInt_t procev = 0;
287 //
288 // 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).
289 //
290 UInt_t runheadtime = 0;
291 UInt_t runtrailtime = 0;
292 UInt_t evfrom = 0;
293 UInt_t evto = 0;
294 UInt_t totfileentries = 0;
295 UInt_t idRun = 0;
296 Int_t id_reg_run=-1;
297 stringstream ftmpname;
298 TString fname;
299 //
300 // define variables for opening and reading level0 file
301 //
302 TFile *l0File = 0;
303 TTree *l0tr = 0;
304 TTree *softinfo = 0;
305 TBranch *l0head = 0;
306 TBranch *l0calo = 0;
307 TBranch *l0trig = 0;
308 pamela::EventHeader *eh = 0;
309 pamela::PscuHeader *ph = 0;
310 pamela::trigger::TriggerEvent *trig = 0;
311 Int_t yodaver = 0;
312 //
313 // Define some basic variables
314 //
315 CaloLevel0 *event = new CaloLevel0(); // NOTICE: very important to call here the constructor!
316 event->SetCrossTalk(crosst);
317 event->SetCrossTalkType(ctground);
318 stringstream file2;
319 stringstream file3;
320 stringstream qy;
321 // Bool_t imtrack = false;
322 Bool_t filled = false;
323 //
324 UInt_t caloevents = 0;
325 stringstream calfile;
326 stringstream aligfile;
327 //
328 Int_t i = -1;
329 Int_t itr = -1;
330 Int_t badevent = 0;
331 Int_t totevent = 0;
332 //
333 UInt_t atime = 0;
334 //
335 Int_t S3 = 0;
336 Int_t S2 = 0;
337 Int_t S12 = 0;
338 Int_t S11 = 0;
339 UInt_t re = 0;
340 UInt_t jumped = 0;
341 //
342 TString caloversion;
343 ItoRunInfo *runinfo = 0;
344 TArrayI *runlist = 0;
345 //
346 Float_t tmptrigty = -1.;
347 Int_t ntrkentry = 0;
348 GL_PARAM *q4 = new GL_PARAM();
349 UInt_t tttrkpar1 = 0;
350 Bool_t trkpar1 = true;
351 GL_ROOT *glroot = new GL_ROOT();
352 GL_TIMESYNC *dbtime = 0;
353 //
354 Long64_t maxsize = 10000000000LL;
355 TTree::SetMaxTreeSize(maxsize);
356 //
357 //
358 TFile *tempfile = 0;
359 TTree *tempcalo = 0;
360 Bool_t myfold = false;
361 stringstream tempname;
362 stringstream calofolder;
363 tempname.str("");
364 tempname << outdir;
365 tempname << "/" << processFolder.Data();
366 calofolder.str("");
367 calofolder << tempname.str().c_str();
368 tempname << "/calotree_run";
369 tempname << run << ".root";
370 //
371 // 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
372 // if run != -1 we must process only that run but first we have to check if the branch calorimeter already exist in the file
373 // if it exists we are reprocessing data and we must delete that entries, if not we must create it.
374 //
375 if ( run == 0 ) reproc = true;
376 //
377 //
378 //
379 if ( !file->IsOpen() ){
380 if ( verbose ) printf(" CALORIMETER - ERROR: cannot open file for writing\n");
381 throw -101;
382 };
383 //
384 if ( withtrk ){
385 //
386 // Does it contain the Tracker tree?
387 //
388 tracker = (TTree*)file->Get("Tracker");
389 if ( !tracker ) {
390 if ( verbose ) printf(" CALORIMETER - ERROR: no tracker tree\n");
391 code = -102;
392 goto closeandexit;
393 };
394 tracker->SetBranchAddress("TrkLevel2",&trk);
395 nevtrkl2 = tracker->GetEntries();
396 };
397 //
398 // Is it a file from simulations?
399 //
400 if ( checksimu ){
401 TTree *h20 = (TTree*)file->Get("h20");
402 if ( h20 ){
403 //
404 // yes!
405 //
406 crosst = false;
407 ctground = true;
408 mechal = true;
409 if ( verbose ) printf("\n\n SIMULATION!!! Setting --use-default-alig and --no-crosstalk flags!\n\n");
410 h20->Delete();
411 };
412 };
413 //
414 // Call runinfo
415 //
416 sgnl = 0;
417 runinfo = new ItoRunInfo(file);
418 //
419 // update versioning information and retrieve informations about the run to be processed
420 //
421 caloversion = CaloInfo(false);
422 sgnl = runinfo->Update(run,"CALO",caloversion);
423 //
424 if ( sgnl ){
425 if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
426 code = sgnl;
427 goto closeandexit;
428 } else {
429 sgnl = 0;
430 };
431 //
432 // number of events in the file BEFORE the first event of our run
433 //
434 nobefrun = runinfo->GetFirstEntry();
435 //
436 // total number of events in the file
437 //
438 totfileentries = runinfo->GetFileEntries();
439 //
440 // first file entry AFTER the last event of our run
441 //
442 noaftrun = runinfo->GetLastEntry() + 1;
443 //
444 // number of run to be processed
445 //
446 numbofrun = runinfo->GetNoRun();
447 //
448 // number of runs in the file
449 //
450 totnorun = runinfo->GetRunEntries();
451 //
452 // Does it contain already a Calorimeter branch? if so we are reprocessing data, if not we must create it.
453 //
454 caloclone = (TTree*)file->Get("Calorimeter");
455 //
456 if ( !caloclone ){
457 reproc = false;
458 if ( run == 0 && verbose ) printf(" CALORIMETER - WARNING: you are reprocessing data but calorimeter tree does not exist!\n");
459 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");
460 //
461 } else {
462 //
463 caloclone->SetAutoSave(900000000000000LL);
464 reproc = true;
465 //
466 if ( verbose ) printf("\n Preparing the pre-processing...\n");
467 //
468 if ( run == 0 || totnorun == 1 ){
469 //
470 // 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.
471 //
472 if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing all runs in the file\n");
473 reprocall = true;
474 //
475 } else {
476 //
477 // we are reprocessing a single run
478 //
479 if ( verbose ) printf("\n CALORIMETER - WARNING: Reprocessing run number %u \n",run);
480 reprocall = false;
481 //
482 //
483 //
484 gSystem->MakeDirectory(calofolder.str().c_str());
485 myfold = true;
486 tempfile = new TFile(tempname.str().c_str(),"RECREATE");
487 tempcalo = caloclone->CloneTree(-1,"fast");
488 tempcalo->SetName("Calorimeter-old");
489 tempfile->Write();
490 tempfile->Close();
491 };
492 //
493 // delete old tree
494 //
495 caloclone->Delete("all");
496 //
497 if ( verbose ) printf("\n ...done!\n");
498 //
499 };
500 //
501 // create calorimeter tree calo
502 //
503 file->cd();
504 calo = new TTree("Calorimeter-new","PAMELA Level2 calorimeter data");
505 calo->SetAutoSave(900000000000000LL);
506 ca->Set();
507 calo->Branch("CaloLevel2","CaloLevel2",&ca);
508 if ( getl1 ) calo->Branch("CaloLevel1","CaloLevel1",&c1);
509 //
510 if ( reproc && !reprocall ){
511 //
512 //
513 //
514 tempfile = new TFile(tempname.str().c_str(),"READ");
515 caloclone = (TTree*)tempfile->Get("Calorimeter-old");
516 caloclone->SetAutoSave(900000000000000LL);
517 caloclone->SetBranchAddress("CaloLevel2",&caclone);
518 TBranch *havel1 = caloclone->GetBranch("CaloLevel1");
519 if ( !havel1 && getl1 ) throw -117;
520 if ( havel1 && !getl1 ) throw -118;
521 if ( getl1 ) caloclone->SetBranchAddress("CaloLevel1",&c1clone);
522 //
523 if ( nobefrun > 0 ){
524 if ( verbose ) printf("\n Pre-processing: copying events from the old tree before the processed run\n");
525 if ( verbose ) printf(" Copying %u events in the file which are before the beginning of the run %u \n",nobefrun,run);
526 if ( verbose ) printf(" Start copying at event number 0, end copying at event number %u \n",nobefrun);
527 for (UInt_t j = 0; j < nobefrun; j++){
528 //
529 caloclone->GetEntry(j);
530 //
531 // copy caclone to ca
532 //
533 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
534 memcpy(&ca,&caclone,sizeof(caclone));
535 //
536 // Fill entry in the new tree
537 //
538 calo->Fill();
539 //
540 if ( getl1 ) c1->Clear();
541 ca->Clear();
542 //
543 };
544 if ( verbose ) printf(" Finished successful copying!\n");
545 };
546 };
547 //
548 // Get the list of run to be processed
549 //
550 runlist = runinfo->GetRunList();
551 //
552 // Loop over the run to be processed
553 //
554 for (UInt_t irun=0; irun < numbofrun; irun++){
555 //
556 badevent = 0;
557 //
558 idRun = runlist->At(irun);
559 if ( verbose ) printf("\n\n\n ####################################################################### \n");
560 if ( verbose ) printf(" PROCESSING RUN NUMBER %u \n",idRun);
561 if ( verbose ) printf(" ####################################################################### \n\n\n");
562 //
563 sgnl = runinfo->GetRunInfo(idRun);
564 if ( sgnl ){
565 if ( verbose ) printf(" CALORIMETER - ERROR: RunInfo exited with non-zero status\n");
566 code = sgnl;
567 goto closeandexit;
568 } else {
569 sgnl = 0;
570 };
571 id_reg_run = runinfo->ID_ROOT_L0;
572 runheadtime = runinfo->RUNHEADER_TIME;
573 runtrailtime = runinfo->RUNTRAILER_TIME;
574 evfrom = runinfo->EV_FROM;
575 evto = runinfo->EV_TO;
576 //
577 if ( id_reg_run == -1 ){
578 if ( verbose ) printf("\n CALORIMETER - ERROR: no run with ID_RUN = %u \n\n Exiting... \n\n",idRun);
579 code = -5;
580 goto closeandexit;
581 };
582 //
583 // prepare the timesync for the db
584 //
585 TString host = glt->CGetHost();
586 TString user = glt->CGetUser();
587 TString psw = glt->CGetPsw();
588 TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
589 if ( !dbc->IsConnected() ) throw -116;
590 stringstream myquery;
591 myquery.str("");
592 myquery << "SET time_zone='+0:00'";
593 dbc->Query(myquery.str().c_str());
594 //
595 // if ( !dbc->IsConnected() ) throw -116;
596 dbtime = new GL_TIMESYNC(runinfo->ID_ROOT_L0,"ID",dbc);
597 //
598 // Search in the DB the path and name of the LEVEL0 file to be processed.
599 //
600 // if ( !dbc->IsConnected() ) throw -116;
601 glroot->Query_GL_ROOT(runinfo->ID_ROOT_L0,dbc);
602 //
603 ftmpname.str("");
604 ftmpname << glroot->PATH.Data() << "/";
605 ftmpname << glroot->NAME.Data();
606 fname = ftmpname.str().c_str();
607 //
608 // print out informations
609 //
610 totevent = runinfo->NEVENTS;
611 if ( verbose ) printf("\n LEVEL0 data file: %s \n",fname.Data());
612 if ( verbose ) printf(" RUN HEADER absolute time is: %u \n",runheadtime);
613 if ( verbose ) printf(" RUN TRAILER absolute time is: %u \n",runtrailtime);
614 if ( verbose ) printf(" %i events to be processed for run %u: from %i to %i (reg entries)\n\n",totevent,idRun,evfrom,evfrom+totevent);
615 //
616 // Open Level0 file
617 //
618 l0File = new TFile(fname.Data());
619 if ( !l0File ) {
620 if ( verbose ) printf(" CALORIMETER - ERROR: problems opening Level0 file\n");
621 code = -6;
622 goto closeandexit;
623 };
624 l0tr = (TTree*)l0File->Get("Physics");
625 if ( !l0tr ) {
626 if ( verbose ) printf(" CALORIMETER - ERROR: no Physics tree in Level0 file\n");
627 l0File->Close();
628 code = -7;
629 goto closeandexit;
630 };
631 l0head = l0tr->GetBranch("Header");
632 if ( !l0head ) {
633 if ( verbose ) printf(" CALORIMETER - ERROR: no Header branch in Level0 tree\n");
634 l0File->Close();
635 code = -8;
636 goto closeandexit;
637 };
638 l0calo = l0tr->GetBranch("Calorimeter");
639 if ( !l0calo ) {
640 if ( verbose ) printf(" CALORIMETER - ERROR: no Calorimeter branch in Level0 tree\n");
641 l0File->Close();
642 code = -103;
643 goto closeandexit;
644 };
645 l0trig = l0tr->GetBranch("Trigger");
646 if ( !l0trig ) {
647 if ( verbose ) printf(" CALORIMETER - ERROR: no Trigger branch in Level0 tree\n");
648 l0File->Close();
649 code = -104;
650 goto closeandexit;
651 };
652 //
653 l0tr->SetBranchAddress("Trigger", &trig);
654 l0tr->SetBranchAddress("Header", &eh);
655 //
656 softinfo = (TTree*)l0File->Get("SoftInfo");
657 if ( softinfo ){
658 softinfo->SetBranchAddress("SoftInfo",&yodaver);
659 softinfo->GetEntry(0);
660 };
661 if ( debug ) printf(" LEVEL0 FILE GENERATED WITH YODA VERSION %i \n",yodaver);
662 //
663 // Construct the event object, look for the calibration which include the first header
664 //
665 sgnl = 0;
666 if ( verbose ) printf(" Check for calorimeter calibrations and initialize event object \n");
667 // if ( !dbc->IsConnected() ) throw -116;
668 event->ProcessingInit(glt,runheadtime,sgnl,l0tr,debug,verbose);
669 if ( verbose ) printf("\n");
670 if ( sgnl == 100 ) {
671 code = sgnl;
672 if ( verbose ) printf(" CALORIMETER - WARNING: run header not included in any calibration interval\n");
673 sgnl = 0;
674 };
675 if ( sgnl ){
676 l0File->Close();
677 code = sgnl;
678 goto closeandexit;
679 };
680 //
681 qy.str("");
682 //
683 nevents = l0head->GetEntries();
684 caloevents = l0calo->GetEntries();
685 //
686 if ( nevents < 1 ) {
687 if ( verbose ) printf(" CALORIMETER - ERROR: Level0 file is empty\n\n");
688 l0File->Close();
689 code = -11;
690 goto closeandexit;
691 };
692 //
693 if ( evto > nevents-1 ) {
694 if ( verbose ) printf(" CALORIMETER - ERROR: too few entries in the registry tree\n");
695 l0File->Close();
696 code = -12;
697 goto closeandexit;
698 };
699 //
700 // Check if we have to load parameter files
701 //
702 sgnl = 0;
703 // if ( !dbc->IsConnected() ) throw -116;
704 sgnl = event->ChkParam(glt,runheadtime,mechal); // calorimeter parameter files
705 if ( sgnl < 0 ){
706 code = sgnl;
707 l0File->Close();
708 goto closeandexit;
709 };
710 //
711 // Calculate the cross talk corrections needed for this run using calibration informations
712 //
713 if ( crosst && !ctground ){
714 sgnl = 0;
715 sgnl = event->CalcCrossTalkCorr(glt,runheadtime);
716 if ( sgnl < 0 ){
717 code = sgnl;
718 l0File->Close();
719 goto closeandexit;
720 };
721 };
722 //
723 //
724 if ( withtrk ){
725 if ( trkpar1 || ( tttrkpar1 != 0 && tttrkpar1 < runheadtime ) ){
726 trkpar1 = false;
727 // if ( !dbc->IsConnected() ) throw -116;
728 Int_t glpar = q4->Query_GL_PARAM(runinfo->RUNHEADER_TIME,1,dbc);
729 if ( glpar < 0 ){
730 code = glpar;
731 goto closeandexit;
732 };
733 tttrkpar1 = q4->TO_TIME;
734 // ----------------------------
735 // Read the magnetic field
736 // ----------------------------
737 if ( verbose ) printf(" Reading magnetic field maps at %s\n",(q4->PATH+q4->NAME).Data());
738 trk->LoadField(q4->PATH+q4->NAME);
739 if ( verbose ) printf("\n");
740 };
741 };
742 //
743 // run over all the events of the run
744 //
745 if ( verbose ) printf("\n Ready to start! \n\n Processed events: \n\n");
746 //
747 if ( dbc ){
748 dbc->Close();
749 delete dbc;
750 };
751 //
752 for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+runinfo->NEVENTS); re++){
753 //for ( re = 46438+200241; re < 46441+200241; re++){
754 // printf(" i = %i \n",re-200241);
755 //
756 if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000);
757 //
758 l0head->GetEntry(re);
759 //
760 // absolute time of this event
761 //
762 ph = eh->GetPscuHeader();
763 atime = dbtime->DBabsTime(ph->GetOrbitalTime());
764 //
765 //
766 //
767 if ( re > caloevents-1 ){
768 if ( verbose ) printf(" CALORIMETER - ERROR: no physics events with entry = %i in Level0 file\n",i);
769 l0File->Close();
770 code = -112;
771 goto closeandexit;
772 };
773 //
774 if ( atime > runtrailtime || atime < runheadtime ) {
775 if ( verbose ) printf(" CALORIMETER - WARNING: event at time outside the run time window, skipping it\n");
776 jumped++;
777 goto jumpev;
778 };
779 //
780 // retrieve tracker informations
781 //
782 if ( !reprocall ){
783 itr = nobefrun + (re - evfrom - jumped);
784 //itr = re-(46438+200241);
785 } else {
786 itr = runinfo->GetFirstEntry() + (re - evfrom - jumped);
787 };
788 //
789 if ( withtrk ){
790 if ( itr > nevtrkl2 ){
791 if ( verbose ) printf(" CALORIMETER - ERROR: no tracker events with entry = %i in Level2 file\n",itr);
792 l0File->Close();
793 code = -113;
794 goto closeandexit;
795 };
796 //
797 trk->Clear();
798 //
799 tracker->GetEntry(itr);
800 //
801 };
802 //
803 procev++;
804 //
805 // start processing
806 //
807 if ( getl1 ) c1->Clear();
808 ca->Clear();
809 //
810 // determine who generate the trigger for this event (TOF, S4/PULSER, CALO)
811 //
812 l0trig->GetEntry(re);
813 S3 = 0;
814 S2 = 0;
815 S12 = 0;
816 S11 = 0;
817 S3 = trig->patterntrig[2];
818 S2 = trig->patterntrig[3];
819 S12 = trig->patterntrig[4];
820 S11 = trig->patterntrig[5];
821 if ( trig->patterntrig[1] & (1<<0) ) tmptrigty = 1.;
822 if ( trig->patterntrig[0] ) tmptrigty = 2.;
823 if ( S3 || S2 || S12 || S11 ) tmptrigty = 0.;
824 if ( !(trig->patterntrig[1] & (1<<0)) && !trig->patterntrig[0] && !S3 && !S2 && !S12 && !S11 ) tmptrigty = 1.;
825 event->clevel2->trigty = tmptrigty;
826 //
827 // check if the calibration we are using is still good, if not load another calibration
828 //
829 sgnl = 0;
830 sgnl = event->ChkCalib(glt,atime);
831 if ( sgnl < 0 ){
832 code = sgnl;
833 goto closeandexit;
834 };
835 if ( sgnl == 100 ){
836 code = sgnl;
837 if ( verbose ) printf(" CALORIMETER - WARNING: data not associated to any calibration interval\n");
838 badevent++;
839 sgnl = 0;
840 };
841 //
842 // do we have at least one track from the tracker? this check has been disabled
843 //
844 event->clevel1->good2 = 1;
845 //
846 // use the whole calorimeter, 22 W planes
847 //
848 event->clevel1->npla = 22;
849 //
850 // Use the standard silicon planes shifting
851 //
852 event->clevel1->reverse = 0;
853 //
854 //
855 // Calibrate calorimeter event "re" and store output in the two structures that will be passed to fortran routine
856 //
857 event->Calibrate(re);
858 //
859 // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract topological variables.
860 //
861 //
862 // Calculate variables common to all tracks (qtot, nstrip, etc.)
863 //
864 event->GetCommonVar();
865 //
866 // Fill common variables
867 //
868 event->FillCommonVar(c1,ca);
869 //
870 // Calculate variables related to tracks only if we have at least one track (from selftrigger and/or tracker)
871 //
872 ntrkentry = 0;
873 //
874 filled = false;
875 //
876 // Run over tracks (tracker or calorimeter )
877 //
878 if ( withtrk ){
879 for(Int_t nt=0; nt < trk->ntrk(); nt++){
880 //
881 event->clevel1->good2 = 1;
882 //
883 TrkTrack *ptt = trk->GetStoredTrack(nt);
884 //
885 event->clevel1->trkchi2 = 0;
886 //
887 // Copy the alpha vector in the input structure
888 //
889 for (Int_t e = 0; e < 5 ; e++){
890 event->clevel1->al_p[e][0] = ptt->al[e];
891 };
892 //
893 // Get tracker related variables for this track
894 //
895 event->GetTrkVar();
896 //
897 // Save tracker track sequence number
898 //
899 event->trkseqno = nt;
900 //
901 // Copy values in the class ca from the structure clevel2
902 //
903 event->FillTrkVar(ca,ntrkentry);
904 ntrkentry++;
905 filled = true;
906 //
907 }; // loop on all the tracks
908 };
909 //
910 // 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
911 // fit of the track (to be used for example when TRK is off due to any reason like IPM3/5 off).
912 // here we make an event selection so it must be done very carefully...
913 //
914 // 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
915 //
916 if ( trackanyway && !filled && event->clevel2->npcfit[0] >= 2 && event->clevel2->npcfit[1] >= 2 && event->clevel2->good != 0 && event->clevel2->trigty < 2. ){
917 if ( debug ) printf(" Event with a track not fitted by the tracker at entry %i \n",itr);
918 //
919 // Disable "track mode" in the fortran routine
920 //
921 event->clevel1->good2 = 0;
922 event->clevel1->riginput = rigdefault;
923 if ( debug ) printf(" Using as default rigidity: %f \n",event->clevel1->riginput);
924 //
925 // We have a selftrigger event to analyze.
926 //
927 for (Int_t e = 0; e < 5 ; e++){
928 event->clevel1->al_p[e][0] = 0.;
929 event->clevel1->al_p[e][1] = 0.;
930 };
931 event->clevel1->trkchi2 = 0;
932 //
933 event->GetTrkVar();
934 //
935 // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
936 //
937 if ( event->clevel1->good2 == 0 ) {
938 //
939 // In selftrigger mode the trkentry variable is set to -1
940 //
941 event->trkseqno = -3;
942 //
943 // Copy values in the class ca from the structure clevel2
944 //
945 event->FillTrkVar(ca,ntrkentry);
946 ntrkentry++;
947 filled = true;
948 //
949 } else {
950 if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
951 };
952 //
953 };
954 //
955 // Call high energy nuclei routine
956 //
957 if ( hZn && event->clevel2->trigty >= 2. ){
958 if ( debug ) printf(" Calling selftrigger high energy nuclei routine for entry %i \n",itr);
959 //
960 // Disable "track mode" in the fortran routine
961 //
962 event->clevel1->good2 = 0;
963 //
964 // Set high energy nuclei flag to one
965 //
966 event->clevel1->hzn = 1;
967 event->clevel1->riginput = rigdefault;
968 //
969 // We have a selftrigger event to analyze.
970 //
971 for (Int_t e = 0; e < 5 ; e++){
972 event->clevel1->al_p[e][0] = 0.;
973 event->clevel1->al_p[e][1] = 0.;
974 };
975 event->clevel1->trkchi2 = 0;
976 //
977 event->GetTrkVar();
978 //
979 // if we had no problem (clevel1->good2 = 0, NOTICE zero, not one in this mode!), fill and go on
980 //
981 if ( event->clevel1->good2 == 0 ) {
982 //
983 // In selftrigger mode the trkentry variable is set to -1
984 //
985 event->trkseqno = -2;
986 //
987 // Copy values in the class ca from the structure clevel2
988 //
989 event->FillTrkVar(ca,ntrkentry);
990 ntrkentry++;
991 filled = true;
992 //
993 } else {
994 if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
995 };
996 //
997 };
998 //
999 // self trigger event
1000 //
1001 if ( st && event->clevel2->trigty >= 2. ){
1002 if ( verbose ) printf(" Selftrigger event at entry %i \n",itr);
1003 //
1004 // Disable "track mode" in the fortran routine
1005 //
1006 event->clevel1->good2 = 0;
1007 //
1008 // disable high enery nuclei flag;
1009 //
1010 event->clevel1->hzn = 0;
1011 //
1012 // We have a selftrigger event to analyze.
1013 //
1014 for (Int_t e = 0; e < 5 ; e++){
1015 event->clevel1->al_p[e][0] = 0.;
1016 event->clevel1->al_p[e][1] = 0.;
1017 };
1018 event->clevel1->trkchi2 = 0;
1019 //
1020 event->GetTrkVar();
1021 //
1022 // if we had no problem (clevel2->good = 0, NOTICE zero, not one in selftrigger mode!), fill and go on
1023 //
1024 if ( event->clevel1->good2 == 0 ) {
1025 //
1026 // In selftrigger mode the trkentry variable is set to -1
1027 //
1028 event->trkseqno = -1;
1029 //
1030 // Copy values in the class ca from the structure clevel2
1031 //
1032 event->FillTrkVar(ca,ntrkentry);
1033 ntrkentry++;
1034 filled = true;
1035 //
1036 } else {
1037 if ( verbose ) printf(" Selftrigger: problems with event at entry %i \n",itr);
1038 };
1039 };
1040 //
1041 if ( !filled ) badevent++;
1042 //
1043 // Clear structures used to communicate with fortran
1044 //
1045 event->ClearStructs();
1046 //
1047 // Fill the rootple
1048 //
1049 calo->Fill();
1050 //
1051 // delete ca;
1052 //
1053 jumpev:
1054 if ( !debug ) debug = false;
1055 //
1056 };
1057 //
1058 if ( verbose ) printf("\n SUMMARY:\n");
1059 if ( verbose ) printf(" Total number of events: %i \n",totevent);
1060 if ( verbose ) printf(" Events with at least one track: %i \n",totevent-badevent);
1061 if ( verbose ) printf(" Events without tracks: %i \n",badevent);
1062 //
1063 if ( badevent == totevent ){
1064 if ( verbose ) printf("\n CALORIMETER - WARNING no tracks or good events in run %u \n",idRun);
1065 code = 101;
1066 };
1067 //
1068 delete dbtime;
1069 //
1070 // Clear variables before processing another run (needed to have always the same result when reprocessing data with the same software).
1071 //
1072 event->RunClose();
1073 if ( l0File ) l0File->Close();
1074 //
1075 }; // process all the runs
1076 //
1077 if ( verbose ) printf("\n Finished processing data \n");
1078 //
1079 closeandexit:
1080 //
1081 if ( !reprocall && reproc && code >= 0 ){
1082 if ( totfileentries > noaftrun ){
1083 if ( verbose ) printf("\n Post-processing: copying events from the old tree after the processed run\n");
1084 if ( verbose ) printf(" Copying %u events in the file which are after the end of the run %u \n",(totfileentries-noaftrun),run);
1085 if ( verbose ) printf(" Start copying at event number %u end copying at event number %u \n",noaftrun,totfileentries);
1086 for (UInt_t j = noaftrun; j < totfileentries; j++ ){
1087 //
1088 // Get entry from old tree
1089 //
1090 caloclone->GetEntry(j);
1091 //
1092 // copy caclone to ca
1093 //
1094 if ( getl1 ) c1->Clear();
1095 ca->Clear();
1096 //
1097 if ( getl1 ) memcpy(&c1,&c1clone,sizeof(c1clone));
1098 memcpy(&ca,&caclone,sizeof(caclone));
1099 //
1100 // Fill entry in the new tree
1101 //
1102 calo->Fill();
1103 //
1104 };
1105 if ( verbose ) printf(" Finished successful copying!\n");
1106 };
1107 };
1108 //
1109 // Case of no errors: close files, delete old tree(s), write and close level2 file
1110 //
1111 if ( l0File ) l0File->Close();
1112 if ( tempfile ) tempfile->Close();
1113 if ( myfold ) gSystem->Unlink(tempname.str().c_str());
1114 if ( tracker ) tracker->Delete();
1115 //
1116 if ( code < 0 ) printf("\n CALORIMETER - ERROR: an error occurred, try to save anyway...\n");
1117 if ( verbose ) printf("\n Writing and closing rootple\n");
1118 if ( runinfo ) runinfo->Close();
1119 if ( calo ) calo->SetName("Calorimeter");
1120 if ( file ){
1121 file->cd();
1122 file->Write();
1123 };
1124 if ( calo ) calo->Delete();
1125 //
1126 if ( myfold ) gSystem->Unlink(calofolder.str().c_str());
1127 //
1128 if ( ca ) delete ca;
1129 if ( caclone ) delete caclone;
1130 if ( c1 ) delete c1;
1131 if ( c1clone ) delete c1clone;
1132 if ( trk ) delete trk;
1133 if ( q4 ) delete q4;
1134 if ( glroot ) delete glroot;
1135 if ( runinfo ) delete runinfo;
1136 //
1137 // the end
1138 //
1139 if ( verbose ) printf("\n Exiting...\n");
1140 if ( code < 0 ) throw code;
1141 return(code);
1142 }
1143
1144
1145

  ViewVC Help
Powered by ViewVC 1.1.23