/[PAMELA software]/DarthVader/S4Level2/src/S4Core.cpp
ViewVC logotype

Annotation of /DarthVader/S4Level2/src/S4Core.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (hide annotations) (download)
Fri Aug 4 10:31:31 2006 UTC (18 years, 4 months ago) by mocchiut
Branch: MAIN
Changes since 1.2: +6 -0 lines
Small memory leaks fixed

1 mocchiut 1.1 //
2     // C/C++ headers
3     //
4     #include <fstream>
5     #include <string.h>
6     #include <iostream>
7     #include <cstring>
8     #include <stdio.h>
9     //
10     // ROOT headers
11     //
12     #include <TTree.h>
13     #include <TClassEdit.h>
14     #include <TObject.h>
15     #include <TList.h>
16     #include <TArrayL.h>
17     #include <TSystem.h>
18     #include <TSystemDirectory.h>
19     #include <TString.h>
20     #include <TFile.h>
21     #include <TClass.h>
22     #include <TSQLServer.h>
23     #include <TSQLRow.h>
24     #include <TSQLResult.h>
25     #include <TClonesArray.h>
26     #include <stdlib.h>
27     #include <math.h>
28     //
29     // YODA headers
30     //
31     #include <PamelaRun.h>
32     #include <RegistryEvent.h>
33     #include <PscuHeader.h>
34     #include <PscuEvent.h>
35     #include <EventHeader.h>
36     #include <RegistryEvent.h>
37     #include <CalibS4Event.h>
38     #include <physics/S4/S4Event.h>
39     //#include <yodaUtility.h>
40     //
41     // RunInfo header
42     //
43     #include <RunInfo.h>
44     #include <GLTables.h>
45     //
46     // This program headers
47     //
48     #include <S4Level2.h>
49     #include <S4Core.h>
50     #include <S4Verl2.h>
51    
52     using namespace std;
53    
54     //
55     // CORE ROUTINE
56     //
57     //
58    
59     int S4Core(ULong64_t run, TFile *file, TSQLServer *dbc, Int_t S4argc, char *S4argv[]){
60     Int_t i = 0;
61     //
62     TString processFolder = "S4Folder";
63     //
64     // Set these to true to have a very verbose output.
65     Bool_t debug = false;
66     //
67     Bool_t verbose = false;
68     //
69     if ( S4argc > 0 ){
70     i = 0;
71     while ( i < S4argc ){
72     if ( !strcmp(S4argv[i],"-processFolder") ) {
73     if ( S4argc < i+1 ){
74     throw -3;
75     };
76     processFolder = (TString)S4argv[i+1];
77     i++;
78     };
79     if ( (!strcmp(S4argv[i],"--debug")) || (!strcmp(S4argv[i],"-g"))) {
80     verbose = true;
81     };
82     if ( (!strcmp(S4argv[i],"--verbose")) || (!strcmp(S4argv[i],"-v"))) {
83     verbose = true;
84     };
85     i++;
86     };
87     };
88    
89     const char* outDir = gSystem->DirName(gSystem->DirName(file->GetPath()));
90     // Variables for level2
91     //
92     TTree *S4tr = 0;
93     Long64_t nevents = 0LL;
94     //
95     // variables needed to reprocess data
96     //
97     TString S4version;
98     ItoRunInfo *runinfo = 0;
99     TArrayL *runlist = 0;
100     TTree *S4trclone = 0;
101     Bool_t reproc = false;
102     Bool_t reprocall = false;
103     UInt_t nobefrun = 0;
104     UInt_t noaftrun = 0;
105     UInt_t numbofrun = 0;
106     stringstream ftmpname;
107     TString fname;
108     Long64_t totfileentries = 0ULL;
109     Long64_t idRun = 0LL;
110     Double_t ParamFit0 = 0ULL;
111     Double_t ParamFit1 = 0ULL;
112     //
113     // variables needed to handle error signals
114     //
115     Int_t code = 0;
116     Int_t sgnl;
117     //
118     // S4 level2 classes
119     //
120     S4Level2 *s4 = new S4Level2();
121     S4Level2 *s4clone = new S4Level2();
122     //
123     // define variables for opening and reading level0 file
124     //
125     TFile *l0File = 0;
126     TTree *l0tr = 0;
127     TBranch *l0head = 0;
128     TBranch *l0registry = 0;
129     TBranch *l0S4 =0;
130     pamela::RegistryEvent *l0reg=0;
131     pamela::S4::S4Event *l0s4e = 0;
132     pamela::EventHeader *eh = 0;
133     //
134     // Define other basic variables
135     //
136     UInt_t procev = 0;
137     stringstream file2;
138     stringstream file3;
139     stringstream qy;
140     Int_t totevent = 0;
141     ULong64_t atime = 0ULL;
142     Int_t ei = 0;
143     Int_t re = 0;
144     //
145     // Working filename
146     //
147     TString outputfile;
148     stringstream name;
149     name.str("");
150     name << outDir << "/";
151     //
152     // temporary file and folder
153     //
154     TFile *tempfile = 0;
155     TTree *tempS4 = 0;
156     stringstream tempname;
157     stringstream S4folder;
158     tempname.str("");
159     tempname << outDir;
160     tempname << "/" << processFolder.Data();
161     S4folder.str("");
162     S4folder << tempname.str().c_str();
163     gSystem->MakeDirectory(S4folder.str().c_str());
164     tempname << "/S4tree_run";
165     tempname << run << ".root";
166     //
167     // DB classes
168     //
169     GL_ROOT *glroot = new GL_ROOT();
170     GL_S4_CALIB *glS4calib = new GL_S4_CALIB();
171     //
172     // Let's start!
173     //
174     // As a first thing we must check what we have to do: if run = 0 we must process all events in the file has been passed
175     // if run != 0 we must process only that run but first we have to check if the tree MyDetector2 already exist in the file
176     // if it exists we are reprocessing data and we must delete that entries, if not we must create it.
177     //
178     if ( run == 0ULL ) reproc = true;
179     //
180     //
181     // Output file is "outputfile"
182     //
183     if ( !file->IsOpen() ){
184     //printf(" S4 - ERROR: cannot open file for writing\n");
185     throw -501;
186     };
187     //
188     // Retrieve GL_RUN variables from the level2 file
189     //
190     S4version = S4Info(false); // we should decide how to handle versioning system
191     //
192     // create an interface to RunInfo called "runinfo"
193     //
194     runinfo = new ItoRunInfo(file);
195     //
196     // open "Run" tree in level2 file, if not existing return an error (sngl != 0)
197     //
198     sgnl = 0;
199     sgnl = runinfo->Update(run,"S4",S4version);
200     if ( sgnl ){
201     //printf(" S4 - ERROR: RunInfo exited with non-zero status\n");
202     code = sgnl;
203     goto closeandexit;
204     } else {
205     sgnl = 0;
206     };
207     //
208     // number of events in the file BEFORE the first event of our run
209     //
210     nobefrun = runinfo->GetFirstEntry();
211     //
212     // total number of events in the file
213     //
214     totfileentries = runinfo->GetFileEntries();
215     //
216     // first file entry AFTER the last event of our run
217     //
218     noaftrun = runinfo->GetLastEntry() + 1;
219     //
220     // number of run to be processed
221     //
222     numbofrun = runinfo->GetNoRun();
223     //
224     // Try to access the S4 tree in the file, if it exists we are reprocessing data if not we are processing a new run
225     //
226     S4trclone = (TTree*)file->Get("S4");
227     //
228     if ( !S4trclone ){
229     //
230     // tree does not exist, we are not reprocessing
231     //
232     reproc = false;
233     if ( run == 0ULL ){
234     if (verbose) printf(" S4 - WARNING: you are reprocessing data but S4 tree does not exist!\n");
235     }
236     if ( runinfo->IsReprocessing() && run != 0ULL ) {
237     if (verbose) printf(" S4 - WARNING: it seems you are not reprocessing data but S4\n versioning information already exists in RunInfo.\n");
238     }
239     } else {
240     //
241     // tree exists, we are reprocessing data. Are we reprocessing a single run or all the file?
242     //
243     reproc = true;
244     //
245     // update versioning information
246     //
247     //sgnl = runinfo->Update(run, "S4",S4version);
248     if (verbose) printf("\n Preparing the pre-processing...\n");
249     //
250     if ( run == 0ULL ){
251     //
252     // we are reprocessing all the file
253     // if we are reprocessing everything we don't need to copy any old event and we can just work with the new tree and delete the old one immediately
254     //
255     reprocall = true;
256     //
257     if (verbose) printf("\n S4 - WARNING: Reprocessing all runs\n");
258     //
259     } else {
260     //
261     // we are reprocessing a single run, we must copy to the new tree the events in the file which preceed the first event of the run
262     //
263     reprocall = false;
264     //
265     if (verbose) printf("\n S4 - WARNING: Reprocessing run number %llu \n",run);
266     //
267     // copying old tree to a new file
268     //
269     tempfile = new TFile(tempname.str().c_str(),"RECREATE");
270     tempS4 = S4trclone->CloneTree(-1,"fast");
271     tempS4->SetName("S4-old");
272     tempfile->Write();
273     tempfile->Close();
274     }
275     //
276     // Delete the old tree from old file and memory
277     //
278     S4trclone->Delete("all");
279     //
280     if (verbose) printf(" ...done!\n");
281     //
282     };
283     //
284     // create mydetector tree mydect
285     //
286     file->cd();
287     S4tr = new TTree("S4-new","PAMELA Level2 S4 data");
288     S4tr->Branch("S4Level2","S4Level2",&s4);
289     //
290     if ( reproc && !reprocall ){
291     //
292     // open new file and retrieve also tree informations
293     //
294     tempfile = new TFile(tempname.str().c_str(),"READ");
295     S4trclone = (TTree*)tempfile->Get("S4-old");
296     S4trclone->SetBranchAddress("S4Level2",&s4clone);
297     //
298     if ( nobefrun > 0 ){
299     if (verbose){
300     printf("\n Pre-processing: copying events from the old tree before the processed run\n");
301     printf(" Copying %u events in the file which are before the beginning of the run %llu \n",nobefrun,run);
302     printf(" Start copying at event number 0, end copying at event number %u \n",nobefrun);
303     }
304     for (UInt_t j = 0; j < nobefrun; j++){
305     //
306     S4trclone->GetEntry(j);
307     //
308     // copy s4clone to mydec
309     //
310 mocchiut 1.2 // s4 = new S4Level2();
311     s4->Clear();
312 mocchiut 1.1 memcpy(&s4,&s4clone,sizeof(s4clone));
313     //
314     // Fill entry in the new tree
315     //
316     S4tr->Fill();
317     //
318     };
319     if (verbose) printf(" Finished successful copying!\n");
320     };
321     };
322     //
323     // Get the list of run to be processed, if only one run has to be processed the list will contain one entry only.
324     //
325     runlist = runinfo->GetRunList();
326     //
327     // Loop over the run to be processed
328     //
329     for (UInt_t irun=0; irun < numbofrun; irun++){
330     //
331     // retrieve the first run ID to be processed using the RunInfo list
332     //
333     idRun = runlist->At(irun);
334     if (verbose){
335     printf("\n\n\n ####################################################################### \n");
336     printf(" PROCESSING RUN NUMBER %i \n",(int)idRun);
337     printf(" ####################################################################### \n\n\n");
338     }
339     //
340     runinfo->ID_REG_RUN = 0ULL;
341     //
342     // store in the runinfo class the GL_RUN variables for our run
343     //
344     sgnl = 0;
345     sgnl = runinfo->GetRunInfo(idRun);
346     if ( sgnl ){
347     //printf(" S4 - ERROR: RunInfo exited with non-zero status\n");
348     code = sgnl;
349     goto closeandexit;
350     } else {
351     sgnl = 0;
352     };
353     //
354     // now you can access that variables using the RunInfo class this way runinfo->ID_REG_RUN
355     //
356     if ( runinfo->ID_REG_RUN == 0 ){
357     //printf("\n S4 - ERROR: no run with ID_RUN = %i \n\n Exiting... \n\n",(int)idRun);
358     code = -5;
359     goto closeandexit;
360     };
361     //
362     // Search in the DB the path and name of the LEVEL0 file to be processed.
363     //
364     glroot->Query_GL_ROOT(runinfo->ID_REG_RUN,dbc);
365     //
366     ftmpname.str("");
367     ftmpname << glroot->PATH.Data() << "/";
368     ftmpname << glroot->NAME.Data();
369     fname = ftmpname.str().c_str();
370     //
371     // print out informations
372     //
373     totevent = runinfo->EV_REG_PHYS_TO - runinfo->EV_REG_PHYS_FROM + 1;
374     if (verbose){
375     printf("\n LEVEL0 data file: %s \n",fname.Data());
376     printf(" RUN HEADER absolute time is: %llu \n",runinfo->RUNHEADER_TIME);
377     printf(" RUN TRAILER absolute time is: %llu \n",runinfo->RUNTRAILER_TIME);
378     printf(" %i events to be processed for run %llu: from %i to %i (reg entries)\n\n",totevent,idRun,runinfo->EV_REG_PHYS_FROM,runinfo->EV_REG_PHYS_TO);
379     }
380     //
381     // Open Level0 file
382     //
383     l0File = new TFile(fname.Data());
384     if ( !l0File ) {
385     //printf(" S4 - ERROR: problems opening Level0 file\n");
386     code = -6;
387     goto closeandexit;
388     };
389     l0tr = (TTree*)l0File->Get("Physics");
390     if ( !l0tr ) {
391     //printf(" S4 - ERROR: no Physics tree in Level0 file\n");
392     l0File->Close();
393     code = -7;
394     goto closeandexit;
395     };
396     l0head = l0tr->GetBranch("Header");
397     if ( !l0head ) {
398     //printf(" S4 - ERROR: no Header branch in Level0 tree\n");
399     l0File->Close();
400     code = -8;
401     goto closeandexit;
402     };
403     l0registry = l0tr->GetBranch("Registry");
404     if ( !l0registry ) {
405     //printf(" S4 - ERROR: no Registry branch in Level0 tree\n");
406     l0File->Close();
407     code = -9;
408     goto closeandexit;
409     };
410     l0S4 = l0tr->GetBranch("S4");
411     if ( !l0S4 ) {
412     //printf(" S4 - ERROR: no S4 branch in Level0 tree\n");
413     l0File->Close();
414     code = -503;
415     goto closeandexit;
416     };
417    
418     //
419     l0tr->SetBranchAddress("Registry", &l0reg);
420     l0tr->SetBranchAddress("S4", &l0s4e);
421     l0tr->SetBranchAddress("Header", &eh);
422     //
423     nevents = l0registry->GetEntries();
424     //
425     if ( nevents < 1 ) {
426     //printf(" S4 - ERROR: Level0 file is empty\n\n");
427     l0File->Close();
428     code = -11;
429     goto closeandexit;
430     };
431     //
432     if ( runinfo->EV_REG_PHYS_TO > nevents-1 ) {
433     //printf(" S4 - ERROR: too few entries in the registry tree\n");
434     l0File->Close();
435     code = -12;
436     goto closeandexit;
437     };
438     //
439     // Check if we have to load parameter files (or calibration associated to runs and not to events)
440     // second query: which is the value of paramfit relative to the calibration?
441     //
442     glS4calib->Query_GL_S4_CALIB(runinfo->RUNHEADER_TIME, dbc);
443     ParamFit0 = glS4calib->PARAM_FIT0;
444     ParamFit1 = glS4calib->PARAM_FIT1;
445     //
446     // run over all the events of the run
447     //
448     if (verbose) printf("\n Ready to start! \n\n Processed events: \n\n");
449     //
450     for ( re = runinfo->EV_REG_PHYS_FROM; re <= runinfo->EV_REG_PHYS_TO; re++){
451     //
452     if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000);
453     //
454     l0registry->GetEntry(re);
455     //
456     // absolute time of this event
457     //
458     atime = l0reg->absTime;
459     //
460     // physics events is at entry number ei where
461     //
462     ei = l0reg->event;
463     //
464     // paranoid check
465     //
466     if ( (atime > runinfo->RUNTRAILER_TIME) || (atime < runinfo->RUNHEADER_TIME) ) {
467     if (verbose) printf(" S4 - WARNING: event at time outside the run time window, skipping it\n");
468     goto jumpev;
469     };
470     //
471     procev++;
472     //
473     // start processing
474     //
475 mocchiut 1.2 // s4 = new S4Level2();
476     s4->Clear();
477 mocchiut 1.1 l0S4->GetEntry(ei);
478     if (l0s4e->unpackError == 1) continue;
479     s4->S4adc = l0s4e->S4_DATA;
480     if ((l0s4e->S4_DATA) > 31 ){
481     s4->S4calibrated = ParamFit0*((l0s4e->S4_DATA)-32)+ParamFit1;
482     }else{
483     s4->S4calibrated = 0;
484     }
485     S4tr->Fill();
486     //
487     //
488     jumpev:
489     debug = false;
490     //
491     };
492     //
493     // Here you may want to clear some variables before processing another run
494     //
495     ei = 0;
496     }; // process all the runs
497     //
498     if (verbose) printf("\n Finished processing data \n");
499     //
500     closeandexit:
501     //
502     // we have finished processing the run(s). If we processed a single run now we must copy all the events after our run from the old tree to the new one and delete the old tree.
503     //
504     if ( !reprocall && reproc && code >= 0 ){
505     if ( totfileentries > noaftrun ){
506     if (verbose){
507     printf("\n Post-processing: copying events from the old tree after the processed run\n");
508     printf(" Copying %i events in the file which are after the end of the run %i \n",(int)(totfileentries-noaftrun),(int)run);
509     printf(" Start copying at event number %i end copying at event number %i \n",(int)noaftrun,(int)totfileentries);
510     }
511     for (UInt_t j = noaftrun; j < totfileentries; j++ ){
512     //
513     // Get entry from old tree
514     //
515     S4trclone->GetEntry(j);
516     //
517     // copy s4clone to s4
518     //
519 mocchiut 1.2 // s4 = new S4Level2();
520     s4->Clear();
521 mocchiut 1.1 memcpy(&s4,&s4clone,sizeof(s4clone));
522     //
523     // Fill entry in the new tree
524     //
525     S4tr->Fill();
526     };
527     if (verbose) printf(" Finished successful copying!\n");
528     };
529     };
530     //
531     // Close files, delete old tree(s), write and close level2 file
532     //
533     if ( l0File ) l0File->Close();
534     if ( tempfile ) tempfile->Close();
535     gSystem->Unlink(tempname.str().c_str());
536    
537     //if ( code < 0 ) printf("\n S4 - ERROR: an error occurred, try to save anyway...\n");
538     //printf("\n Writing and closing rootple\n");
539     if ( runinfo ) runinfo->Close();
540     if ( S4tr ) S4tr->SetName("S4");
541     if ( file ){
542     file->cd();
543     file->Write();
544     };
545     //
546     gSystem->Unlink(S4folder.str().c_str());
547     //
548     // the end
549     //
550     if (verbose) printf("\n Exiting...\n");
551     if(S4tr)S4tr->Delete();
552 mocchiut 1.3 //
553     if ( s4 ) delete s4;
554     if ( s4clone ) delete s4clone;
555     if ( glroot ) delete glroot;
556     if ( runinfo ) delete runinfo;
557     //
558 mocchiut 1.1 if(code < 0) throw code;
559     return(code);
560     }
561    
562    
563    

  ViewVC Help
Powered by ViewVC 1.1.23