1 |
// |
// |
2 |
|
// Given a calibration and a data file this program create an ntuple with LEVEL2 S4 variables |
3 |
|
// |
4 |
|
// |
5 |
// C/C++ headers |
// C/C++ headers |
6 |
// |
// |
7 |
#include <fstream> |
#include <fstream> |
32 |
#include <stdlib.h> |
#include <stdlib.h> |
33 |
#include <math.h> |
#include <math.h> |
34 |
// |
// |
35 |
|
// RunInfo header |
36 |
|
// |
37 |
|
#include <RunInfo.h> |
38 |
|
#include <GLTables.h> |
39 |
|
// |
40 |
// YODA headers |
// YODA headers |
41 |
// |
// |
42 |
#include <PamelaRun.h> |
#include <PamelaRun.h> |
45 |
#include <EventHeader.h> |
#include <EventHeader.h> |
46 |
#include <CalibS4Event.h> |
#include <CalibS4Event.h> |
47 |
#include <physics/S4/S4Event.h> |
#include <physics/S4/S4Event.h> |
|
//#include <yodaUtility.h> |
|
|
// |
|
|
// RunInfo header |
|
|
// |
|
|
#include <RunInfo.h> |
|
|
#include <GLTables.h> |
|
48 |
// |
// |
49 |
// This program headers |
// This program headers |
50 |
// |
// |
51 |
#include <S4Level2.h> |
#include <S4Level2.h> |
52 |
#include <S4Core.h> |
#include <S4Core.h> |
53 |
#include <S4Verl2.h> |
#include <S4Verl2.h> |
|
|
|
|
using namespace std; |
|
|
|
|
|
// |
|
54 |
// |
// |
55 |
|
using namespace std; |
56 |
// |
// |
57 |
/* |
/* |
58 |
* Fit function Received from Valeria Malvezzi 06/02/2006 |
* Fit function Received from Valeria Malvezzi 06/02/2006 |
71 |
// |
// |
72 |
GL_S4_CALIB *glS4calib = new GL_S4_CALIB(); |
GL_S4_CALIB *glS4calib = new GL_S4_CALIB(); |
73 |
// |
// |
74 |
glS4calib->Query_GL_S4_CALIB(atime, dbc); |
// if ( !dbc->IsConnected() ) throw -504; |
75 |
|
Int_t s4sig = glS4calib->Query_GL_S4_CALIB(atime, dbc); |
76 |
|
if ( s4sig != -57 && s4sig < 0 ){ |
77 |
|
parametri->AddAt(0.,0); |
78 |
|
parametri->AddAt(0.,1); |
79 |
|
return parametri; |
80 |
|
}; |
81 |
|
if ( s4sig < 0 ) throw s4sig; |
82 |
// |
// |
83 |
GL_ROOT *glroot = new GL_ROOT(); |
GL_ROOT *glroot = new GL_ROOT(); |
84 |
|
// if ( !dbc->IsConnected() ) throw -504; |
85 |
glroot->Query_GL_ROOT(glS4calib->ID_ROOT_L0,dbc); |
glroot->Query_GL_ROOT(glS4calib->ID_ROOT_L0,dbc); |
86 |
// |
// |
87 |
stringstream ftmpname; |
stringstream ftmpname; |
107 |
// |
// |
108 |
tr->GetEntry(glS4calib->EV_ROOT); |
tr->GetEntry(glS4calib->EV_ROOT); |
109 |
// |
// |
110 |
//----------- variable initialization ------------------------------------------------- |
// Variables initialization |
111 |
// |
// |
112 |
Double_t mip[3]={1, 30, 300}; |
Double_t mip[3]={1, 30, 300}; |
113 |
Double_t adc[3] = {0.,0.,0.}; |
Double_t adc[3] = {0.,0.,0.}; |
114 |
// |
// |
115 |
//------------ Fit calibrations and find parameters to calibrate data ------------------ |
// Fit calibrations and find parameters to calibrate data |
116 |
// |
// |
117 |
pamela::S4::S4Event *s4Record = 0; |
pamela::S4::S4Event *s4Record = 0; |
118 |
// |
// |
141 |
adc[2]=adc[2]/128; |
adc[2]=adc[2]/128; |
142 |
// if ( IsDebug() ) printf(" adc1 = %g adc2 = %g adc3 = %g\n",adc[0],adc[1],adc[2]); |
// if ( IsDebug() ) printf(" adc1 = %g adc2 = %g adc3 = %g\n",adc[0],adc[1],adc[2]); |
143 |
TGraph *fitpar = new TGraph (3, adc, mip); |
TGraph *fitpar = new TGraph (3, adc, mip); |
144 |
TF1 *func = new TF1("fitf", fitf, -1., 400., 2); // definizione della funzione, 2 = num. parametri |
TF1 *func = new TF1("fitf", fitf, -1., 400., 2); // function definition with 2 parameters |
|
// TF1 *func = new TF1("fitf", fitf, -0., 1000., 2); // definizione della funzione, 2 = num. parametri |
|
145 |
// |
// |
146 |
func->SetParameters(0.3,1.); //inizializzazione dei parametri a 1 |
func->SetParameters(0.3,1.); // parameters initialization to 1 |
147 |
func->SetParNames("m","q"); //definisce il nome dei parametri |
func->SetParNames("m","q"); //parameter's name |
148 |
fitpar->Fit(func,"qr"); //fitta fitpar con la funzione func nel range definito nella funzione |
fitpar->Fit(func,"qr"); //function fit |
|
//fitpar->Fit(func,"r"); //fitta fitpar con la funzione func nel range definito nella funzione |
|
149 |
parametri->AddAt(func->GetParameter(0),0); |
parametri->AddAt(func->GetParameter(0),0); |
150 |
parametri->AddAt(func->GetParameter(1),1); |
parametri->AddAt(func->GetParameter(1),1); |
151 |
// if ( parametri[0] < 0. || parametri[0] > 0.5 || parametri[1] < 0.7 || parametri[1] > 1.1 ) valid = 0; |
// if ( parametri[0] < 0. || parametri[0] > 0.5 || parametri[1] < 0.7 || parametri[1] > 1.1 ) valid = 0; |
163 |
// CORE ROUTINE |
// CORE ROUTINE |
164 |
// |
// |
165 |
// |
// |
166 |
int S4Core(UInt_t run, TFile *file, TSQLServer *dbc, Int_t S4argc, char *S4argv[]){ |
int S4Core(UInt_t run, TFile *file,GL_TABLES *glt, Int_t S4argc, char *S4argv[]){ |
|
Int_t i = 0; |
|
167 |
// |
// |
168 |
TString processFolder = "S4Folder"; |
// Set these to true to have a verbose output. |
169 |
// |
// |
170 |
Bool_t debug = false; |
Bool_t debug = false; |
|
// |
|
171 |
Bool_t verbose = false; |
Bool_t verbose = false; |
172 |
// |
// |
173 |
|
// |
174 |
|
// Output directory is the working directoy. |
175 |
|
// |
176 |
|
const char* outDir = gSystem->DirName(gSystem->DirName(file->GetPath())); |
177 |
|
// |
178 |
|
Int_t i = 0; |
179 |
|
TString processFolder = Form("S4Folder_%u",run); |
180 |
if ( S4argc > 0 ){ |
if ( S4argc > 0 ){ |
181 |
i = 0; |
i = 0; |
182 |
while ( i < S4argc ){ |
while ( i < S4argc ){ |
189 |
}; |
}; |
190 |
if ( (!strcmp(S4argv[i],"--debug")) || (!strcmp(S4argv[i],"-g"))) { |
if ( (!strcmp(S4argv[i],"--debug")) || (!strcmp(S4argv[i],"-g"))) { |
191 |
verbose = true; |
verbose = true; |
192 |
|
debug = true; |
193 |
}; |
}; |
194 |
if ( (!strcmp(S4argv[i],"--verbose")) || (!strcmp(S4argv[i],"-v"))) { |
if ( (!strcmp(S4argv[i],"--verbose")) || (!strcmp(S4argv[i],"-v"))) { |
195 |
verbose = true; |
verbose = true; |
197 |
i++; |
i++; |
198 |
}; |
}; |
199 |
}; |
}; |
|
// |
|
|
const char* outDir = gSystem->DirName(gSystem->DirName(file->GetPath())); |
|
|
// |
|
200 |
// Variables for level2 |
// Variables for level2 |
201 |
// |
// |
202 |
TTree *S4tr = 0; |
TTree *S4tr = 0; |
203 |
UInt_t nevents = 0; |
UInt_t nevents = 0; |
204 |
// |
// |
205 |
// variables needed to reprocess data |
// Variables needed to reprocess data |
206 |
|
// |
207 |
|
Long64_t maxsize = 10000000000LL; |
208 |
|
TTree::SetMaxTreeSize(maxsize); |
209 |
// |
// |
210 |
TString S4version; |
TString S4version; |
211 |
ItoRunInfo *runinfo = 0; |
ItoRunInfo *runinfo = 0; |
267 |
TTree *tempS4 = 0; |
TTree *tempS4 = 0; |
268 |
stringstream tempname; |
stringstream tempname; |
269 |
stringstream S4folder; |
stringstream S4folder; |
270 |
|
Bool_t myfold = false; |
271 |
tempname.str(""); |
tempname.str(""); |
272 |
tempname << outDir; |
tempname << outDir; |
273 |
tempname << "/" << processFolder.Data(); |
tempname << "/" << processFolder.Data(); |
274 |
S4folder.str(""); |
S4folder.str(""); |
275 |
S4folder << tempname.str().c_str(); |
S4folder << tempname.str().c_str(); |
|
gSystem->MakeDirectory(S4folder.str().c_str()); |
|
276 |
tempname << "/S4tree_run"; |
tempname << "/S4tree_run"; |
277 |
tempname << run << ".root"; |
tempname << run << ".root"; |
278 |
// |
// |
332 |
// number of run to be processed |
// number of run to be processed |
333 |
// |
// |
334 |
numbofrun = runinfo->GetNoRun(); |
numbofrun = runinfo->GetNoRun(); |
335 |
|
UInt_t totnorun = runinfo->GetRunEntries(); |
336 |
// |
// |
337 |
// Try to access the S4 tree in the file, if it exists we are reprocessing data if not we are processing a new run |
// Try to access the S4 tree in the file, if it exists we are reprocessing data if not we are processing a new run |
338 |
// |
// |
353 |
// |
// |
354 |
// tree exists, we are reprocessing data. Are we reprocessing a single run or all the file? |
// tree exists, we are reprocessing data. Are we reprocessing a single run or all the file? |
355 |
// |
// |
356 |
|
S4trclone->SetAutoSave(900000000000000LL); |
357 |
reproc = true; |
reproc = true; |
358 |
// |
// |
359 |
// update versioning information |
// update versioning information |
360 |
// |
// |
361 |
if (verbose) printf("\n Preparing the pre-processing...\n"); |
if (verbose) printf("\n Preparing the pre-processing...\n"); |
362 |
// |
// |
363 |
if ( run == 0 ){ |
if ( run == 0 || totnorun == 1 ){ |
364 |
// |
// |
365 |
// we are reprocessing all the file |
// we are reprocessing all the file |
366 |
// 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 |
// 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 |
379 |
// |
// |
380 |
// copying old tree to a new file |
// copying old tree to a new file |
381 |
// |
// |
382 |
|
gSystem->MakeDirectory(S4folder.str().c_str()); |
383 |
|
myfold = true; |
384 |
tempfile = new TFile(tempname.str().c_str(),"RECREATE"); |
tempfile = new TFile(tempname.str().c_str(),"RECREATE"); |
385 |
tempS4 = S4trclone->CloneTree(-1,"fast"); |
tempS4 = S4trclone->CloneTree(-1,"fast"); |
386 |
tempS4->SetName("S4-old"); |
tempS4->SetName("S4-old"); |
400 |
// |
// |
401 |
file->cd(); |
file->cd(); |
402 |
S4tr = new TTree("S4-new","PAMELA Level2 S4 data"); |
S4tr = new TTree("S4-new","PAMELA Level2 S4 data"); |
403 |
|
S4tr->SetAutoSave(900000000000000LL); |
404 |
S4tr->Branch("S4Level2","S4Level2",&s4); |
S4tr->Branch("S4Level2","S4Level2",&s4); |
405 |
// |
// |
406 |
if ( reproc && !reprocall ){ |
if ( reproc && !reprocall ){ |
409 |
// |
// |
410 |
tempfile = new TFile(tempname.str().c_str(),"READ"); |
tempfile = new TFile(tempname.str().c_str(),"READ"); |
411 |
S4trclone = (TTree*)tempfile->Get("S4-old"); |
S4trclone = (TTree*)tempfile->Get("S4-old"); |
412 |
|
S4trclone->SetAutoSave(900000000000000LL); |
413 |
S4trclone->SetBranchAddress("S4Level2",&s4clone); |
S4trclone->SetBranchAddress("S4Level2",&s4clone); |
414 |
// |
// |
415 |
if ( nobefrun > 0 ){ |
if ( nobefrun > 0 ){ |
478 |
// |
// |
479 |
// prepare the timesync for the db |
// prepare the timesync for the db |
480 |
// |
// |
481 |
|
TString host = glt->CGetHost(); |
482 |
|
TString user = glt->CGetUser(); |
483 |
|
TString psw = glt->CGetPsw(); |
484 |
|
TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data()); |
485 |
|
if ( !dbc->IsConnected() ) throw -504; |
486 |
|
stringstream myquery; |
487 |
|
myquery.str(""); |
488 |
|
myquery << "SET time_zone='+0:00'"; |
489 |
|
dbc->Query(myquery.str().c_str()); |
490 |
dbtime = new GL_TIMESYNC(runinfo->ID_ROOT_L0,"ID",dbc); |
dbtime = new GL_TIMESYNC(runinfo->ID_ROOT_L0,"ID",dbc); |
491 |
// |
// |
492 |
// Search in the DB the path and name of the LEVEL0 file to be processed. |
// Search in the DB the path and name of the LEVEL0 file to be processed. |
493 |
// |
// |
494 |
|
// if ( !dbc->IsConnected() ) throw -504; |
495 |
glroot->Query_GL_ROOT(runinfo->ID_ROOT_L0,dbc); |
glroot->Query_GL_ROOT(runinfo->ID_ROOT_L0,dbc); |
496 |
// |
// |
497 |
ftmpname.str(""); |
ftmpname.str(""); |
509 |
printf(" %i events to be processed for run %u: from %u to %u \n\n",totevent,idRun,runinfo->EV_FROM,(runinfo->EV_FROM+totevent)); |
printf(" %i events to be processed for run %u: from %u to %u \n\n",totevent,idRun,runinfo->EV_FROM,(runinfo->EV_FROM+totevent)); |
510 |
} |
} |
511 |
// |
// |
512 |
|
// if ( !totevent ) goto closeandexit; |
513 |
|
// |
514 |
// Open Level0 file |
// Open Level0 file |
515 |
// |
// |
516 |
l0File = new TFile(fname.Data()); |
l0File = new TFile(fname.Data()); |
547 |
// |
// |
548 |
nevents = l0S4->GetEntries(); |
nevents = l0S4->GetEntries(); |
549 |
// |
// |
550 |
if ( nevents < 1 ) { |
if ( nevents < 1 && totevent ) { |
551 |
if ( debug ) printf(" S4 - ERROR: Level0 file is empty\n\n"); |
if ( debug ) printf(" S4 - ERROR: Level0 file is empty\n\n"); |
552 |
l0File->Close(); |
l0File->Close(); |
553 |
code = -11; |
code = -11; |
554 |
goto closeandexit; |
goto closeandexit; |
555 |
}; |
}; |
556 |
// |
// |
557 |
if ( runinfo->EV_TO > nevents-1 ) { |
if ( runinfo->EV_TO > nevents-1 && totevent ) { |
558 |
if ( debug ) printf(" S4 - ERROR: too few entries in the S4 tree\n"); |
if ( debug ) printf(" S4 - ERROR: too few entries in the S4 tree\n"); |
559 |
l0File->Close(); |
l0File->Close(); |
560 |
code = -12; |
code = -12; |
578 |
// |
// |
579 |
if (verbose) printf("\n Ready to start! \n\n Processed events: \n\n"); |
if (verbose) printf("\n Ready to start! \n\n Processed events: \n\n"); |
580 |
// |
// |
581 |
|
if ( dbc ){ |
582 |
|
dbc->Close(); |
583 |
|
// delete dbc; |
584 |
|
}; |
585 |
|
// |
586 |
for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+runinfo->NEVENTS); re++){ |
for ( re = runinfo->EV_FROM; re < (runinfo->EV_FROM+runinfo->NEVENTS); re++){ |
587 |
// |
// |
588 |
if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000); |
if ( procev%1000 == 0 && procev > 0 && verbose ) printf(" %iK \n",procev/1000); |
596 |
// |
// |
597 |
// paranoid check |
// paranoid check |
598 |
// |
// |
599 |
if ( (atime > runinfo->RUNTRAILER_TIME) || (atime < runinfo->RUNHEADER_TIME) ) { |
if ( (atime > (runinfo->RUNTRAILER_TIME+1)) || (atime < (runinfo->RUNHEADER_TIME-1)) ) { |
600 |
if (verbose) printf(" S4 - WARNING: event at time outside the run time window, skipping it\n"); |
if (verbose) printf(" S4 - WARNING: event at time outside the run time window, skipping it\n"); |
601 |
goto jumpev; |
goto jumpev; |
602 |
}; |
}; |
671 |
// |
// |
672 |
if ( l0File ) l0File->Close(); |
if ( l0File ) l0File->Close(); |
673 |
if ( tempfile ) tempfile->Close(); |
if ( tempfile ) tempfile->Close(); |
674 |
gSystem->Unlink(tempname.str().c_str()); |
if ( myfold) gSystem->Unlink(tempname.str().c_str()); |
675 |
// |
// |
676 |
if ( runinfo ) runinfo->Close(); |
if ( runinfo ) runinfo->Close(); |
677 |
if ( S4tr ) S4tr->SetName("S4"); |
if ( S4tr ) S4tr->SetName("S4"); |
680 |
file->Write(); |
file->Write(); |
681 |
}; |
}; |
682 |
// |
// |
683 |
gSystem->Unlink(S4folder.str().c_str()); |
if ( myfold ) gSystem->Unlink(S4folder.str().c_str()); |
684 |
// |
// |
685 |
// the end |
// the end |
686 |
// |
// |