]>
Commit | Line | Data |
---|---|---|
53f515d3 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | /* $Id$ */ | |
17 | ||
18 | /* 21/11/07 | |
3f637cfc | 19 | |
90293ba6 | 20 | Version for MUONTRKda MUON tracking |
21 | (A. Baldisseri, J.-L. Charvet & Ch. Finck) | |
3f637cfc | 22 | |
23 | ||
24 | Rem: AliMUON2DMap stores all channels, even those which are not connected | |
25 | if pedMean == -1, channel not connected to a pad | |
26 | ||
3f637cfc | 27 | |
28 | */ | |
29 | extern "C" { | |
30 | #include <daqDA.h> | |
31 | } | |
32 | ||
33 | #include "event.h" | |
34 | #include "monitor.h" | |
35 | ||
36 | #include <Riostream.h> | |
37 | #include <stdio.h> | |
38 | #include <stdlib.h> | |
39 | ||
40 | //AliRoot | |
41 | #include "AliMUONRawStreamTracker.h" | |
42 | #include "AliMUONDspHeader.h" | |
43 | #include "AliMUONBlockHeader.h" | |
44 | #include "AliMUONBusStruct.h" | |
45 | #include "AliMUONDDLTracker.h" | |
2bf54af2 | 46 | #include "AliMUONVStore.h" |
3f637cfc | 47 | #include "AliMUON2DMap.h" |
0aa03544 | 48 | #include "AliMUONCalibParamND.h" |
3f637cfc | 49 | #include "AliMpIntPair.h" |
2bf54af2 | 50 | #include "AliMpConstants.h" |
3f637cfc | 51 | #include "AliRawReaderDate.h" |
52 | ||
3f637cfc | 53 | //ROOT |
54 | #include "TFile.h" | |
53f515d3 | 55 | #include "TSystem.h" |
3f637cfc | 56 | #include "TTree.h" |
57 | #include "TH1F.h" | |
58 | #include "TString.h" | |
59 | #include "TStopwatch.h" | |
60 | #include "TMath.h" | |
61 | #include "TTimeStamp.h" | |
62 | #include "TGraphErrors.h" | |
63 | #include "TF1.h" | |
fd99693f | 64 | #include "TROOT.h" |
65 | #include "TPluginManager.h" | |
a76f513d | 66 | #include "TFitter.h" |
3f637cfc | 67 | |
a76f513d | 68 | #define NFITPARAMS 4 |
3f637cfc | 69 | |
70 | // global variables | |
90293ba6 | 71 | const Int_t gkNChannels = AliMpConstants::ManuNofChannels(); |
72 | const Int_t gkADCMax = 4095; | |
73 | ||
74 | AliMUONVStore* gPedestalStore = new AliMUON2DMap(kFALSE); | |
75 | ||
76 | Int_t gNManu = 0; | |
77 | Int_t gNChannel = 0; | |
78 | UInt_t gRunNumber = 0; | |
79 | Int_t gNEvents = 0; | |
53f515d3 | 80 | Int_t gGlitchErrors= 0; |
81 | Int_t gParityErrors= 0; | |
82 | Int_t gPaddingErrors= 0; | |
90293ba6 | 83 | Int_t gNDateEvents = 0; |
84 | Int_t gPrintLevel = 1; // global printout variable | |
85 | Int_t gPlotLevel = 0; // global plot variable | |
86 | ||
87 | TH1F* gPedMeanHisto = 0x0; | |
88 | TH1F* gPedSigmaHisto = 0x0; | |
89 | Char_t gHistoFileName[256]; | |
90 | ||
91 | // used by makegain | |
53f515d3 | 92 | Char_t gHistoFileName_gain[256]="MUONTRKda_gain_data.root"; |
93 | Char_t gRootFileName[256]; | |
94 | Char_t gOutFolder[256]="."; | |
95 | Char_t filename[256]; | |
96 | Char_t filenam[256]="MUONTRKda_gain"; | |
97 | // Char_t filenam[256]="MUONTRKda_gain.param"; // if gPrintLevel = 2 | |
98 | // Char_t filenam_bad[256]="MUONTRKda_gain.bad"; // if gPrintLevel = 2 | |
99 | // Char_t filenam_peak[256]="MUONTRKda_gain.peak"; // if gPrintLevel = 3 | |
100 | // Char_t filenam_out[256]="MUONTRKda_gain.out"; | |
101 | Char_t flatFile[256]; | |
102 | ||
103 | ofstream filcout; | |
104 | ||
105 | TString flatOutputFile; | |
106 | TString logOutputFile; | |
90293ba6 | 107 | TString gCommand("ped"); |
53f515d3 | 108 | TTimeStamp date; |
3f637cfc | 109 | |
110 | // funtions | |
111 | ||
90293ba6 | 112 | |
113 | //________________ | |
114 | Double_t funcLin(Double_t *x, Double_t *par) | |
115 | { | |
116 | return par[0] + par[1]*x[0]; | |
117 | } | |
118 | ||
3f637cfc | 119 | //________________ |
90293ba6 | 120 | Double_t funcParabolic(Double_t *x, Double_t *par) |
3f637cfc | 121 | { |
90293ba6 | 122 | return par[0]*x[0]*x[0]; |
123 | } | |
3f637cfc | 124 | |
90293ba6 | 125 | //________________ |
126 | Double_t funcCalib(Double_t *x, Double_t *par) | |
127 | { | |
128 | Double_t xLim= par[3]; | |
3f637cfc | 129 | |
90293ba6 | 130 | if(x[0] <= xLim) return par[0] + par[1]*x[0]; |
3f637cfc | 131 | |
90293ba6 | 132 | Double_t yLim = par[0]+ par[1]*xLim; |
133 | return yLim + par[1]*(x[0] - xLim) + par[2]*(x[0] - xLim)*(x[0] - xLim); | |
3f637cfc | 134 | } |
135 | ||
136 | ||
3f637cfc | 137 | //__________ |
138 | void MakePed(Int_t busPatchId, Int_t manuId, Int_t channelId, Int_t charge) | |
139 | { | |
140 | ||
141 | AliMUONVCalibParam* ped = | |
90293ba6 | 142 | static_cast<AliMUONVCalibParam*>(gPedestalStore->FindObject(busPatchId, manuId)); |
3f637cfc | 143 | |
144 | if (!ped) { | |
90293ba6 | 145 | gNManu++; |
146 | ped = new AliMUONCalibParamND(2, gkNChannels,busPatchId, manuId, -1.); // put default wise -1, not connected channel | |
147 | gPedestalStore->Add(ped); | |
3f637cfc | 148 | } |
149 | ||
53f515d3 | 150 | if (gNEvents == 0) { |
0aa03544 | 151 | ped->SetValueAsDouble(channelId, 0, 0.); |
152 | ped->SetValueAsDouble(channelId, 1, 0.); | |
3f637cfc | 153 | } |
154 | ||
0aa03544 | 155 | Double_t pedMean = ped->ValueAsDouble(channelId, 0) + charge; |
156 | Double_t pedSigma = ped->ValueAsDouble(channelId, 1) + charge*charge; | |
3f637cfc | 157 | |
0aa03544 | 158 | ped->SetValueAsDouble(channelId, 0, pedMean); |
159 | ped->SetValueAsDouble(channelId, 1, pedSigma); | |
3f637cfc | 160 | |
161 | } | |
162 | ||
163 | //________________ | |
164 | void MakePedStore(TString flatOutputFile = "") | |
165 | { | |
53f515d3 | 166 | // TTimeStamp date; |
0aa03544 | 167 | Double_t pedMean; |
168 | Double_t pedSigma; | |
3f637cfc | 169 | ofstream fileout; |
170 | Int_t busPatchId; | |
171 | Int_t manuId; | |
172 | Int_t channelId; | |
173 | ||
174 | // histo | |
53f515d3 | 175 | TFile* histoFile = 0; |
176 | TTree* tree = 0; | |
177 | if (gCommand.CompareTo("ped") == 0) | |
178 | { | |
179 | sprintf(gHistoFileName,"%s/MUONTRKda_ped_%d.root",gOutFolder,gRunNumber); | |
180 | histoFile = new TFile(gHistoFileName,"RECREATE","MUON Tracking pedestals"); | |
181 | // } | |
3f637cfc | 182 | |
183 | Char_t name[255]; | |
184 | Char_t title[255]; | |
185 | sprintf(name,"pedmean_allch"); | |
186 | sprintf(title,"Pedestal mean all channels"); | |
53f515d3 | 187 | Int_t nx = 4096; |
3f637cfc | 188 | Int_t xmin = 0; |
53f515d3 | 189 | Int_t xmax = 4095; |
90293ba6 | 190 | gPedMeanHisto = new TH1F(name,title,nx,xmin,xmax); |
191 | gPedMeanHisto->SetDirectory(histoFile); | |
3f637cfc | 192 | |
193 | sprintf(name,"pedsigma_allch"); | |
194 | sprintf(title,"Pedestal sigma all channels"); | |
53f515d3 | 195 | nx = 201; |
3f637cfc | 196 | xmin = 0; |
53f515d3 | 197 | xmax = 200; |
90293ba6 | 198 | gPedSigmaHisto = new TH1F(name,title,nx,xmin,xmax); |
199 | gPedSigmaHisto->SetDirectory(histoFile); | |
3f637cfc | 200 | |
53f515d3 | 201 | // TTree* tree = new TTree("t","Pedestal tree"); |
202 | tree = new TTree("t","Pedestal tree"); | |
3f637cfc | 203 | tree->Branch("bp",&busPatchId,"bp/I"); |
204 | tree->Branch("manu",&manuId,",manu/I"); | |
205 | tree->Branch("channel",&channelId,",channel/I"); | |
90293ba6 | 206 | tree->Branch("pedMean",&pedMean,",pedMean/D"); |
207 | tree->Branch("pedSigma",&pedSigma,",pedSigma/D"); | |
3f637cfc | 208 | |
53f515d3 | 209 | } |
210 | ||
3f637cfc | 211 | if (!flatOutputFile.IsNull()) { |
212 | fileout.open(flatOutputFile.Data()); | |
213 | fileout<<"//===========================================================================" << endl; | |
214 | fileout<<"// Pedestal file calculated by MUONTRKda"<<endl; | |
215 | fileout<<"//===========================================================================" << endl; | |
90293ba6 | 216 | fileout<<"// * Run : " << gRunNumber << endl; |
3f637cfc | 217 | fileout<<"// * Date : " << date.AsString("l") <<endl; |
90293ba6 | 218 | fileout<<"// * Statictics : " << gNEvents << endl; |
219 | fileout<<"// * # of MANUS : " << gNManu << endl; | |
220 | fileout<<"// * # of channels : " << gNChannel << endl; | |
3f637cfc | 221 | fileout<<"//"<<endl; |
222 | fileout<<"//---------------------------------------------------------------------------" << endl; | |
223 | fileout<<"//---------------------------------------------------------------------------" << endl; | |
90293ba6 | 224 | fileout<<"// BP MANU CH. MEAN SIGMA"<<endl; |
3f637cfc | 225 | fileout<<"//---------------------------------------------------------------------------" << endl; |
226 | ||
227 | } | |
228 | ||
229 | // iterator over pedestal | |
90293ba6 | 230 | TIter next(gPedestalStore->CreateIterator()); |
2bf54af2 | 231 | AliMUONVCalibParam* ped; |
3f637cfc | 232 | |
2bf54af2 | 233 | while ( ( ped = dynamic_cast<AliMUONVCalibParam*>(next() ) ) ) |
3f637cfc | 234 | { |
2bf54af2 | 235 | busPatchId = ped->ID0(); |
236 | manuId = ped->ID1(); | |
3f637cfc | 237 | |
238 | for (channelId = 0; channelId < ped->Size() ; ++channelId) { | |
239 | ||
0aa03544 | 240 | pedMean = ped->ValueAsDouble(channelId, 0); |
3f637cfc | 241 | |
242 | if (pedMean > 0) { // connected channels | |
243 | ||
90293ba6 | 244 | ped->SetValueAsDouble(channelId, 0, pedMean/(Double_t)gNEvents); |
3f637cfc | 245 | |
0aa03544 | 246 | pedMean = ped->ValueAsDouble(channelId, 0); |
247 | pedSigma = ped->ValueAsDouble(channelId, 1); | |
3f637cfc | 248 | |
90293ba6 | 249 | ped->SetValueAsDouble(channelId, 1, TMath::Sqrt(TMath::Abs(pedSigma/(Double_t)gNEvents - pedMean*pedMean))); |
3f637cfc | 250 | |
0aa03544 | 251 | pedMean = ped->ValueAsDouble(channelId, 0) + 0.5 ; |
90293ba6 | 252 | // pedMean = ped->ValueAsDouble(channelId, 0) ; |
0aa03544 | 253 | pedSigma = ped->ValueAsDouble(channelId, 1); |
3f637cfc | 254 | |
255 | ||
256 | if (!flatOutputFile.IsNull()) { | |
257 | fileout << "\t" << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" | |
258 | << pedMean <<"\t"<< pedSigma << endl; | |
259 | } | |
260 | ||
53f515d3 | 261 | if (gCommand.CompareTo("ped") == 0) |
262 | { | |
90293ba6 | 263 | gPedMeanHisto->Fill(pedMean); |
264 | gPedSigmaHisto->Fill(pedSigma); | |
3f637cfc | 265 | |
266 | tree->Fill(); | |
267 | } | |
53f515d3 | 268 | } |
3f637cfc | 269 | } |
3f637cfc | 270 | } |
271 | ||
272 | // file outputs | |
273 | if (!flatOutputFile.IsNull()) | |
274 | fileout.close(); | |
275 | ||
53f515d3 | 276 | if (gCommand.CompareTo("ped") == 0) |
277 | { | |
278 | histoFile->Write(); | |
279 | histoFile->Close(); | |
280 | } | |
3f637cfc | 281 | |
282 | // delete tree; | |
283 | ||
3f637cfc | 284 | } |
285 | ||
286 | //________________ | |
90293ba6 | 287 | void MakePedStoreForGain(Int_t injCharge) |
3f637cfc | 288 | { |
289 | // store pedestal map in root file | |
290 | ||
90293ba6 | 291 | // Int_t injCharge = 200; |
3f637cfc | 292 | |
293 | TTree* tree = 0x0; | |
294 | ||
295 | // compute and store pedestals | |
53f515d3 | 296 | sprintf(flatFile,"%s/%s_%d_DAC_%d.ped",gOutFolder,filenam,gRunNumber,injCharge); |
297 | cout << "MUONTRKda : Flat file generated = " << flatFile << "\n"; | |
298 | MakePedStore(flatFile); | |
3f637cfc | 299 | |
300 | TString mode("UPDATE"); | |
301 | ||
90293ba6 | 302 | if (gCommand.Contains("cre")) { |
3f637cfc | 303 | mode = "RECREATE"; |
304 | } | |
90293ba6 | 305 | TFile* histoFile = new TFile(gHistoFileName_gain, mode.Data(), "MUON Tracking Gains"); |
3f637cfc | 306 | |
307 | // second argument should be the injected charge, taken from config crocus file | |
308 | // put also info about run number could be usefull | |
90293ba6 | 309 | AliMpIntPair* pair = new AliMpIntPair(gRunNumber, injCharge); |
3f637cfc | 310 | |
311 | if (mode.CompareTo("UPDATE") == 0) { | |
312 | tree = (TTree*)histoFile->Get("t"); | |
313 | tree->SetBranchAddress("run",&pair); | |
90293ba6 | 314 | tree->SetBranchAddress("ped",&gPedestalStore); |
3f637cfc | 315 | |
316 | } else { | |
317 | tree = new TTree("t","Pedestal tree"); | |
318 | tree->Branch("run", "AliMpIntPair",&pair); | |
90293ba6 | 319 | tree->Branch("ped", "AliMUON2DMap",&gPedestalStore); |
3f637cfc | 320 | tree->SetBranchAddress("run",&pair); |
90293ba6 | 321 | tree->SetBranchAddress("ped",&gPedestalStore); |
3f637cfc | 322 | |
323 | } | |
324 | ||
325 | tree->Fill(); | |
326 | tree->Write("t", TObject::kOverwrite); // overwrite the tree | |
327 | histoFile->Close(); | |
328 | ||
329 | delete pair; | |
330 | } | |
331 | ||
332 | //________________ | |
53f515d3 | 333 | // void MakeGainStore(TString flatOutputFile) |
334 | void MakeGainStore() | |
3f637cfc | 335 | { |
90293ba6 | 336 | Double_t goodA1Min = 0.5; |
337 | Double_t goodA1Max = 2.; | |
338 | Double_t goodA2Min = -0.5E-03; | |
339 | Double_t goodA2Max = 1.E-03; | |
340 | ||
3f637cfc | 341 | // open file mutrkgain.root |
342 | // read again the pedestal for the calibration runs (9 runs ?) | |
343 | // need the injection charge from config file (to be done) | |
344 | // For each channel make a TGraphErrors (mean, sigma) vs injected charge | |
345 | // Fit with a polynomial fct | |
346 | // store the result in a flat file. | |
347 | ||
53f515d3 | 348 | |
349 | TFile* histoFile = new TFile(gHistoFileName_gain); | |
350 | ||
351 | AliMUON2DMap* map[11]; | |
352 | AliMUONVCalibParam* ped[11]; | |
353 | AliMpIntPair* run[11]; | |
354 | ||
355 | //read back from root file | |
356 | TTree* tree = (TTree*)histoFile->Get("t"); | |
357 | Int_t nEntries = tree->GetEntries(); | |
358 | ||
359 | // read back info | |
360 | for (Int_t i = 0; i < nEntries; ++i) { | |
361 | map[i] = 0x0; | |
362 | run[i] = 0x0; | |
363 | tree->SetBranchAddress("ped",&map[i]); | |
364 | tree->SetBranchAddress("run",&run[i]); | |
365 | tree->GetEvent(i); | |
366 | // std::cout << map[i] << " " << run[i] << std::endl; | |
367 | } | |
368 | gRunNumber=(Double_t)run[0]->GetFirst(); | |
369 | ||
370 | // some print | |
371 | cout<<"\n ******** MUONTRKda for Gain computing (Run = " << gRunNumber << ")\n" << endl; | |
372 | cout<<" * Date : " << date.AsString("l") << "\n" << endl; | |
373 | cout << " Entries = " << nEntries << " DAC values \n" << endl; | |
374 | for (Int_t i = 0; i < nEntries; ++i) { | |
375 | cout<< " Run = " << (Double_t)run[i]->GetFirst() << " DAC = " << (Double_t)run[i]->GetSecond() << endl; | |
376 | } | |
377 | cout << "" << endl; | |
378 | ||
379 | ||
90293ba6 | 380 | Double_t pedMean[11]; |
381 | Double_t pedSigma[11]; | |
382 | Double_t injCharge[11]; | |
383 | Double_t injChargeErr[11]; | |
3f637cfc | 384 | |
90293ba6 | 385 | // full print out |
386 | ||
53f515d3 | 387 | sprintf(filename,"%s/%s_%d.log",gOutFolder,filenam,gRunNumber); |
388 | logOutputFile=filename; | |
389 | ||
390 | filcout.open(logOutputFile.Data()); | |
391 | filcout<<"//====================================================" << endl; | |
392 | filcout<<"// MUONTRKda for Gain computing (Run = " << gRunNumber << ")" << endl; | |
393 | filcout<<"//====================================================" << endl; | |
394 | filcout<<"// * Date : " << date.AsString("l") << "\n" << endl; | |
395 | ||
396 | // cout<<"\n ******** MUONTRKda for Gain computing (Run = " << gRunNumber << ")\n" << endl; | |
397 | // cout<<" * Date : " << date.AsString("l") << "\n" << endl; | |
398 | ||
399 | ||
400 | ||
90293ba6 | 401 | // why 2 files ? (Ch. F.) |
402 | FILE *pfilen = 0; | |
53f515d3 | 403 | FILE *pfilef = 0; |
404 | if(gPrintLevel>=2) | |
90293ba6 | 405 | { |
53f515d3 | 406 | sprintf(filename,"%s/%s_%d.param",gOutFolder,filenam,gRunNumber); |
407 | cout << " fit parameter file = " << filename << "\n"; | |
408 | pfilen = fopen (filename,"w"); | |
3f637cfc | 409 | |
90293ba6 | 410 | fprintf(pfilen,"//===================================================================\n"); |
411 | fprintf(pfilen,"// BP MANU CH. a0 a1 a2 xlim P(chi2) P(chi2)2 Q\n"); | |
412 | fprintf(pfilen,"//===================================================================\n"); | |
53f515d3 | 413 | fprintf(pfilen,"// * Run : %d \n",gRunNumber); |
414 | fprintf(pfilen,"//===================================================================\n"); | |
415 | ||
416 | sprintf(filename,"%s/%s_%d.bad",gOutFolder,filenam,gRunNumber); | |
417 | cout << " Bad channel file = " << filename << "\n"; | |
418 | pfilef = fopen (filename,"w"); | |
419 | ||
420 | fprintf(pfilef,"//=================================================\n"); | |
421 | fprintf(pfilef,"// Bad Channel file calculated by MUONTRKda \n"); | |
422 | fprintf(pfilef,"//=================================================\n"); | |
423 | fprintf(pfilef,"// * Run : %d \n",gRunNumber); | |
424 | fprintf(pfilef,"// * Date : %s \n",date.AsString("l")); | |
425 | fprintf(pfilef,"//=======================================\n"); | |
426 | fprintf(pfilef,"// BP MANU CH. a1 a2 thres. Q\n"); | |
427 | fprintf(pfilef,"//=======================================\n"); | |
3f637cfc | 428 | } |
429 | ||
90293ba6 | 430 | FILE *pfilew=0; |
53f515d3 | 431 | if(flatOutputFile.IsNull()) |
432 | { | |
433 | sprintf(filename,"%s_%d.par",filenam,gRunNumber); | |
434 | flatOutputFile=filename; | |
435 | } | |
436 | if(!flatOutputFile.IsNull()) | |
90293ba6 | 437 | { |
438 | pfilew = fopen (flatOutputFile.Data(),"w"); | |
439 | ||
440 | fprintf(pfilew,"//=================================================\n"); | |
441 | fprintf(pfilew,"// Calibration file calculated by MUONTRKda \n"); | |
442 | fprintf(pfilew,"//=================================================\n"); | |
443 | fprintf(pfilew,"// * Run : %d \n",gRunNumber); | |
444 | fprintf(pfilew,"// * Date : %s \n",date.AsString("l")); | |
445 | fprintf(pfilew,"// * Statictics : %d \n",gNEvents); | |
446 | fprintf(pfilew,"// * # of MANUS : %d \n",gNManu); | |
447 | fprintf(pfilew,"// * # of channels : %d \n",gNChannel); | |
448 | fprintf(pfilew,"//-------------------------------------------------\n"); | |
449 | fprintf(pfilew,"//=======================================\n"); | |
450 | fprintf(pfilew,"// BP MANU CH. a1 a2 thres. Q\n"); | |
451 | fprintf(pfilew,"//=======================================\n"); | |
452 | } | |
453 | ||
53f515d3 | 454 | FILE *pfilep = 0; |
455 | if(gPrintLevel==3) | |
456 | { | |
457 | sprintf(filename,"%s/%s_%d.peak",gOutFolder,filenam,gRunNumber); | |
458 | cout << " File containing Peak mean values = " << filename << "\n"; | |
459 | pfilep = fopen (filename,"w"); | |
460 | ||
461 | fprintf(pfilep,"//===============================================================================================================================\n"); | |
462 | fprintf(pfilep,"// * Run : %d \n",gRunNumber); | |
463 | fprintf(pfilep,"//===============================================================================================================================\n"); | |
464 | fprintf(pfilep,"// BP MANU CH. Ped. <0> <1> <2> <3> <4> <5> <6> <7> <8> <9> <10> \n"); | |
465 | // fprintf(pfilep,"// %9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f fC\n",level_fC[0],level_fC[1],level_fC[2],level_fC[3],level_fC[4],level_fC[5],level_fC[6],level_fC[7],level_fC[8],level_fC[9],level_fC[10]); | |
466 | // fprintf(pfilep,"// %9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f\n",level_err[0],level_err[1],level_err[2],level_err[3],level_err[4],level_err[5],level_err[6],level_err[7],level_err[8],level_err[9],level_err[10]); | |
467 | fprintf(pfilep,"//===============================================================================================================================\n"); | |
468 | } | |
90293ba6 | 469 | |
470 | ||
3f637cfc | 471 | |
90293ba6 | 472 | // plot out |
3f637cfc | 473 | |
90293ba6 | 474 | TFile* gainFile = 0x0; |
53f515d3 | 475 | sprintf(gRootFileName,"%s/%s_%d.root",gOutFolder,filenam,gRunNumber); |
90293ba6 | 476 | gainFile = new TFile(gRootFileName,"RECREATE"); |
477 | ||
478 | Double_t chi2 = 0.; | |
479 | Double_t chi2P2 = 0.; | |
480 | Double_t prChi2 = 0; | |
481 | Double_t prChi2P2 =0; | |
482 | Double_t a0,a1,a2; | |
483 | Int_t busPatchId ; | |
484 | Int_t manuId ; | |
485 | Int_t channelId ; | |
486 | Int_t threshold = 0; | |
487 | Int_t Q = 0; | |
53f515d3 | 488 | Int_t p1 ; |
489 | Int_t p2 ; | |
490 | Double_t gain; | |
491 | Double_t capa=0.2; // internal capacitor (pF) | |
90293ba6 | 492 | |
493 | TTree *tg = new TTree("tg","TTree avec class Manu_DiMu"); | |
494 | ||
495 | tg->Branch("bp",&busPatchId, "busPatchId/I"); | |
496 | tg->Branch("manu",&manuId, "manuId/I"); | |
497 | tg->Branch("channel",&channelId, "channelId/I"); | |
498 | ||
499 | tg->Branch("a0",&a0, "a0/D"); | |
500 | tg->Branch("a1",&a1, "a1/D"); | |
501 | tg->Branch("a2",&a2, "a2/D"); | |
502 | tg->Branch("Pchi2",&prChi2, "prChi2/D"); | |
503 | tg->Branch("Pchi2_2",&prChi2P2, "prChi2P2/D"); | |
504 | tg->Branch("Threshold",&threshold, "threshold/I"); | |
505 | tg->Branch("Q",&Q, "Q/I"); | |
53f515d3 | 506 | tg->Branch("p1",&p1, "p1/I"); |
507 | tg->Branch("p2",&p2, "p2/I"); | |
508 | tg->Branch("gain",&gain, "gain/D"); | |
509 | ||
510 | // bad BusPatch and manu | |
511 | Int_t num_tot_BP=800; | |
512 | Int_t num_tot_Manu=1500; | |
513 | // Int_t bad_channel[num_tot_BP][num_tot_Manu]; | |
514 | Int_t bad_channel[800][1500]; | |
515 | for ( Int_t i = 0; i < num_tot_BP ; i++ ) | |
516 | { for ( Int_t j = 0; j < num_tot_Manu ; j++ ) bad_channel[i][j]=0;} | |
3f637cfc | 517 | |
3f637cfc | 518 | |
90293ba6 | 519 | char graphName[256]; |
3f637cfc | 520 | |
521 | // iterates over the first pedestal run | |
0aa03544 | 522 | TIter next(map[0]->CreateIterator()); |
2bf54af2 | 523 | AliMUONVCalibParam* p; |
3f637cfc | 524 | |
90293ba6 | 525 | Int_t nmanu = 0; |
53f515d3 | 526 | Int_t nGoodChannel = 0; |
527 | Int_t nGoodChannel_a1 = 0; | |
90293ba6 | 528 | Int_t nBadChannel = 0; |
53f515d3 | 529 | Int_t nBadChannel_a1 = 0; |
530 | Int_t nBadChannel_a2 = 0; | |
531 | Int_t nplot=0; | |
90293ba6 | 532 | Double_t sumProbChi2 = 0.; |
533 | Double_t sumA1 = 0.; | |
534 | Double_t sumProbChi2P2 = 0.; | |
535 | Double_t sumA2 = 0.; | |
536 | ||
537 | Double_t x[11], xErr[11], y[11], yErr[11]; | |
53f515d3 | 538 | Double_t xp[11], xpErr[11], yp[11], ypErr[11]; |
90293ba6 | 539 | |
2bf54af2 | 540 | while ( ( p = dynamic_cast<AliMUONVCalibParam*>(next() ) ) ) |
3f637cfc | 541 | { |
2bf54af2 | 542 | ped[0] = p; |
3f637cfc | 543 | |
90293ba6 | 544 | busPatchId = p->ID0(); |
545 | manuId = p->ID1(); | |
3f637cfc | 546 | |
547 | // read back pedestal from the other runs for the given (bupatch, manu) | |
3f637cfc | 548 | for (Int_t i = 1; i < nEntries; ++i) { |
0aa03544 | 549 | ped[i] = static_cast<AliMUONVCalibParam*>(map[i]->FindObject(busPatchId, manuId)); |
3f637cfc | 550 | } |
551 | ||
552 | // compute for each channel the gain parameters | |
90293ba6 | 553 | for ( channelId = 0; channelId < ped[0]->Size() ; ++channelId ) { |
3f637cfc | 554 | |
53f515d3 | 555 | gain=0.4; |
556 | ||
3f637cfc | 557 | Int_t n = 0; |
558 | for (Int_t i = 0; i < nEntries; ++i) { | |
559 | ||
560 | if (!ped[i]) continue; //shouldn't happen. | |
0aa03544 | 561 | pedMean[i] = ped[i]->ValueAsDouble(channelId, 0); |
562 | pedSigma[i] = ped[i]->ValueAsDouble(channelId, 1); | |
563 | injCharge[i] = (Double_t)run[i]->GetSecond(); | |
90293ba6 | 564 | injChargeErr[i] = 0.01*injCharge[i]; |
565 | if(injChargeErr[i] <= 1.) injChargeErr[i]=1.; | |
566 | ||
567 | // if(n<2)cout << nEntries << " " << i << " " << injCharge[i] << endl; | |
568 | ||
569 | // cout << busPatchId << "\t" << manuId <<"\t"<< channelId << "\t" << n << " " << pedMean[i] << " " << pedSigma[i] << " " << injCharge[i] << " " << injChargeErr[i] << endl; | |
3f637cfc | 570 | |
fd99693f | 571 | if (pedMean[i] < 0) continue; // not connected |
3f637cfc | 572 | |
573 | if (pedSigma[i] <= 0) pedSigma[i] = 1.; // should not happen. | |
574 | n++; | |
575 | } | |
576 | ||
53f515d3 | 577 | |
578 | // print_peak_mean_values | |
579 | if(gPrintLevel==3) | |
580 | { | |
581 | ||
582 | fprintf(pfilep,"%4i%5i%5i%10.3f",busPatchId,manuId,channelId,0.); | |
583 | fprintf(pfilep,"%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f mV\n",pedMean[0],pedMean[1],pedMean[2],pedMean[3],pedMean[4],pedMean[5],pedMean[6],pedMean[7],pedMean[8],pedMean[9],pedMean[10]); | |
584 | fprintf(pfilep," %9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f%9.1f \n",pedSigma[0],pedSigma[1],pedSigma[2],pedSigma[3],pedSigma[4],pedSigma[5],pedSigma[6],pedSigma[7],pedSigma[8],pedSigma[9],pedSigma[10]); | |
585 | } | |
586 | ||
587 | // makegain | |
90293ba6 | 588 | |
589 | ||
590 | // Fit Method: Linear fit over 6 points + fit parabolic function over 3 points) | |
591 | ||
592 | // 1. - linear fit over 6 points | |
593 | ||
594 | Double_t par[4] = {0.,0.,0.,gkADCMax}; | |
595 | ||
596 | Int_t nInit = 1; | |
597 | Int_t nbs = nEntries - nInit; | |
598 | Int_t nbpf1 = 6; // linear fit over nbf1 points | |
599 | ||
600 | for (Int_t j = 0; j < nbs; ++j) | |
601 | { | |
602 | Int_t k = j + nInit; | |
603 | x[j] = pedMean[k]; | |
604 | xErr[j] = pedSigma[k]; | |
605 | y[j] = injCharge[k]; | |
606 | yErr[j] = injChargeErr[k]; | |
607 | ||
608 | } | |
609 | ||
53f515d3 | 610 | TF1 *f1 = new TF1("f1",funcLin,0.,gkADCMax,2); |
611 | TGraphErrors *graphErr = new TGraphErrors(nbpf1, x, y, xErr, yErr); | |
90293ba6 | 612 | |
613 | f1->SetParameters(0,0); | |
614 | ||
615 | graphErr->Fit("f1","RQ"); | |
616 | ||
617 | chi2 = f1->GetChisquare(); | |
618 | f1->GetParameters(par); | |
619 | ||
53f515d3 | 620 | delete graphErr; |
621 | graphErr=0; | |
622 | delete f1; | |
623 | ||
90293ba6 | 624 | prChi2 = TMath::Prob(chi2, nbpf1 - 2); |
625 | ||
626 | Double_t xLim = pedMean[nInit + nbpf1 - 1]; | |
627 | Double_t yLim = par[0]+par[1] * xLim; | |
628 | ||
629 | a0 = par[0]; | |
630 | a1 = par[1]; | |
631 | ||
90293ba6 | 632 | // 2. - Translation : new origin (xLim, yLim) + parabolic fit over nbf2 points |
633 | ||
53f515d3 | 634 | // par[0] = 0; |
90293ba6 | 635 | Int_t nbpf2 = nEntries - (nInit + nbpf1) + 1; |
636 | ||
637 | if(nbpf2 > 1) | |
638 | { | |
53f515d3 | 639 | for (Int_t j = 0; j < nbpf2; j++) |
90293ba6 | 640 | { |
641 | Int_t k = j + (nInit + nbpf1) - 1; | |
53f515d3 | 642 | xp[j] = pedMean[k] - xLim; |
643 | xpErr[j] = pedSigma[k]; | |
90293ba6 | 644 | |
53f515d3 | 645 | yp[j] = injCharge[k] - yLim - par[1]*xp[j]; |
646 | ypErr[j] = injChargeErr[k]; | |
90293ba6 | 647 | |
3f637cfc | 648 | } |
90293ba6 | 649 | |
53f515d3 | 650 | TF1 *f2 = new TF1("f2",funcParabolic,0.,gkADCMax,1); |
651 | TGraphErrors *graphErr = new TGraphErrors(nbpf2, xp, yp, xpErr, ypErr); | |
90293ba6 | 652 | |
653 | graphErr->Fit(f2,"RQ"); | |
654 | chi2P2 = f2->GetChisquare(); | |
655 | f2->GetParameters(par); | |
656 | ||
53f515d3 | 657 | delete graphErr; |
658 | graphErr=0; | |
659 | delete f2; | |
660 | ||
90293ba6 | 661 | prChi2P2 = TMath::Prob(chi2P2, nbpf2-1); |
662 | ||
53f515d3 | 663 | |
664 | // ------------- print out in log file | |
665 | // if (busPatchId == 6 && manuId == 116 && ( channelId >= 17 && channelId <= 20) ) | |
666 | // { | |
667 | // filcout << " \n ********! Print_out.: BP= " << busPatchId << " Manu_Id= " << manuId | |
668 | // << " Ch.= " << channelId << ":" << endl; | |
669 | ||
670 | // for (Int_t j = 0; j < nbpf1; ++j) | |
671 | // {filcout << j << " " << x[j] << " " << xErr[j] << " " << y[j] << " " << yErr[j] << endl;} | |
672 | // filcout << " a0,a1 = " << a0 << " , " << a1 << " pr_chi2 = " << prChi2 << endl ; | |
673 | ||
674 | // for (Int_t j = 0; j < nbpf2; ++j) | |
675 | // {filcout << j << " " << xp[j] << " " << xpErr[j] << " " << yp[j] << " " << ypErr[j] << endl;} | |
676 | // filcout << " a2 = " << par[0] << " pr_chi2_2 = " << prChi2P2 << endl; | |
677 | ||
678 | // } | |
679 | // ------------------------------------------ | |
680 | ||
681 | ||
682 | ||
90293ba6 | 683 | a2 = par[0]; |
684 | ||
685 | par[0] = a0; | |
686 | par[1] = a1; | |
687 | par[2] = a2; | |
688 | par[3] = xLim; | |
689 | ||
53f515d3 | 690 | // delete graphErr; |
90293ba6 | 691 | |
3f637cfc | 692 | } |
693 | ||
90293ba6 | 694 | // Prints |
695 | ||
53f515d3 | 696 | p1 = TMath::Nint(ceil(prChi2*14))+1; |
697 | p2 = TMath::Nint(ceil(prChi2P2*14))+1; | |
90293ba6 | 698 | |
699 | Double_t x0 = -par[0]/par[1]; // value of x corresponding to à 0 fC | |
700 | threshold = TMath::Nint(ceil(par[3]-x0)); // linear if x < threshold | |
701 | ||
53f515d3 | 702 | if(gPrintLevel>=2) |
90293ba6 | 703 | { |
704 | fprintf(pfilen,"%4i %4i %2i",busPatchId,manuId,channelId); | |
53f515d3 | 705 | fprintf(pfilen," %6.2f %6.4f %10.3e %4.2f %5.3f %x %5.3f %x\n", |
706 | par[0], par[1], par[2], par[3], prChi2, p1, prChi2P2, p2); | |
90293ba6 | 707 | } |
708 | ||
709 | // some tests | |
710 | ||
53f515d3 | 711 | if(par[1]< goodA1Min || par[1]> goodA1Max) |
90293ba6 | 712 | { |
53f515d3 | 713 | p1=0; |
714 | nBadChannel_a1++; | |
715 | if (gPrintLevel && nBadChannel_a1 < 1) | |
90293ba6 | 716 | { |
53f515d3 | 717 | cout << " !!!!! " << nBadChannel_a1 << " !!!!!!!! Bad Calib.: BP= " << busPatchId << " Manu_Id= " << manuId << |
90293ba6 | 718 | " Ch.= " << channelId << ":"; |
719 | cout << " a1 = " << par[1] << " out of limit : [" << goodA1Min << "," << goodA1Max << | |
720 | "]" << endl; | |
721 | } | |
90293ba6 | 722 | } |
53f515d3 | 723 | |
724 | if(par[2]< goodA2Min || par[2]> goodA2Max) | |
90293ba6 | 725 | { |
53f515d3 | 726 | p2=0; |
727 | nBadChannel_a2++; | |
728 | if (gPrintLevel && nBadChannel_a2 < 1) | |
90293ba6 | 729 | { |
53f515d3 | 730 | cout << " !!!!! " << nBadChannel_a2 << " !!!!!!!! Bad Calib.: BP= " << busPatchId << " Manu_Id= " << manuId |
90293ba6 | 731 | << " Ch.= " << channelId << ":"; |
732 | cout << " a2 = " << par[2] << " out of limit : [" << goodA2Min << "," << goodA2Max | |
733 | << "]" << endl; | |
53f515d3 | 734 | |
735 | for (Int_t j = 0; j < nbpf2; ++j) | |
736 | {cout << j << " " << x[j] << " " << xErr[j] << " " << y[j] << " " << yErr[j] << endl;} | |
90293ba6 | 737 | } |
90293ba6 | 738 | } |
90293ba6 | 739 | |
53f515d3 | 740 | Q = p1*16 + p2; // fit quality |
741 | if(p1==0)Q=0; // bad linear fit <=> bad calibration | |
742 | ||
743 | if(p1>0 && p2>0) | |
90293ba6 | 744 | { |
53f515d3 | 745 | nGoodChannel++; |
746 | sumProbChi2P2 += prChi2P2; | |
747 | sumA2 += par[2]; | |
748 | } | |
749 | else | |
750 | { | |
751 | nBadChannel++; | |
752 | if(busPatchId < num_tot_BP && manuId < num_tot_Manu) bad_channel[busPatchId][manuId]++; | |
753 | else{cout << " Warning : busPatch = " << busPatchId << " Manu = " << manuId << endl;} | |
754 | if(gPrintLevel>=2)fprintf(pfilef,"%4i %5i %2i %7.4f %10.3e %4i %2x\n",busPatchId,manuId,channelId,par[1],par[2],threshold,Q); | |
90293ba6 | 755 | } |
53f515d3 | 756 | |
757 | ||
758 | if(p1>0) | |
759 | { | |
760 | nGoodChannel_a1++; | |
761 | sumProbChi2 += prChi2; | |
762 | sumA1 += par[1]; | |
763 | gain=1./(par[1]*capa); | |
764 | } | |
765 | ||
90293ba6 | 766 | |
767 | tg->Fill(); | |
768 | ||
769 | if (!flatOutputFile.IsNull()) | |
770 | { | |
771 | fprintf(pfilew,"%4i %5i %2i %7.4f %10.3e %4i %2x\n",busPatchId,manuId,channelId,par[1],par[2],threshold,Q); | |
772 | } | |
773 | ||
774 | // Plots | |
775 | ||
776 | if(gPlotLevel){ | |
53f515d3 | 777 | if(Q==0 and nplot < 100) |
778 | // if(p1>1 && p2==0 and nplot < 100) | |
779 | // if(p1>1 && p2>1 and nplot < 100) | |
780 | { | |
781 | nplot++; | |
782 | TF1 *f2Calib = new TF1("f2Calib",funcCalib,0.,gkADCMax,NFITPARAMS); | |
783 | ||
784 | TGraphErrors *graphErr = new TGraphErrors(nEntries,pedMean,injCharge,pedSigma,injChargeErr); | |
785 | ||
786 | sprintf(graphName,"BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId); | |
787 | ||
788 | graphErr->SetTitle(graphName); | |
789 | graphErr->SetMarkerColor(3); | |
790 | graphErr->SetMarkerStyle(12); | |
791 | graphErr->Write(graphName); | |
792 | ||
793 | sprintf(graphName,"f2_BusPatch_%d_Manu_%d_Ch_%d",busPatchId, manuId,channelId); | |
794 | f2Calib->SetTitle(graphName); | |
795 | f2Calib->SetLineColor(4); | |
796 | f2Calib->SetParameters(par); | |
797 | f2Calib->Write(graphName); | |
798 | ||
799 | delete graphErr; | |
800 | delete f2Calib; | |
801 | } | |
90293ba6 | 802 | } |
3f637cfc | 803 | } |
90293ba6 | 804 | nmanu++; |
53f515d3 | 805 | if(fmod(nmanu,100)==0)std::cout << " Nb manu = " << nmanu << std::endl; |
3f637cfc | 806 | } |
807 | ||
808 | // file outputs for gain | |
53f515d3 | 809 | // if (!flatOutputFile.IsNull()) fileout.close(); |
810 | if (!flatOutputFile.IsNull()) fclose(pfilew); | |
811 | if(gPrintLevel==2){ fclose(pfilen); fclose(pfilef);} | |
812 | if(gPrintLevel==3) fclose(pfilep); | |
90293ba6 | 813 | |
814 | tg->Write(); | |
815 | histoFile->Close(); | |
3f637cfc | 816 | |
90293ba6 | 817 | //OutPut |
818 | if (gPrintLevel) | |
819 | { | |
53f515d3 | 820 | filcout << "\n List of problematic BusPatch and Manu " << endl; |
821 | filcout << " ========================================" << endl; | |
822 | filcout << " BP Manu Nb Channel " << endl ; | |
823 | filcout << " ========================================" << endl; | |
824 | for ( Int_t i = 0 ; i < num_tot_BP ; i++ ) | |
825 | { for ( Int_t j = 0 ; j < num_tot_Manu ; j++ ) | |
826 | if (bad_channel[i][j] != 0 ) filcout << "\t" << i << "\t " << j << "\t\t" << bad_channel[i][j] << endl;} | |
827 | filcout << " ========================================" << endl; | |
828 | ||
829 | ||
830 | filcout << "\n Nb of channels in raw data = " << nmanu*64 << " (" << nmanu << " Manu)" << endl; | |
831 | filcout << "\n Nb of fully calibrated channel = " << nGoodChannel << " (" << goodA1Min << "<a1<" << goodA1Max | |
90293ba6 | 832 | << " and " << goodA2Min << "<a2<" << goodA2Max << ") " << endl; |
53f515d3 | 833 | filcout << "\n Nb of Bad channel = " << nBadChannel << endl; |
834 | ||
835 | filcout << "\n Nb of Good a1 channels = " << nGoodChannel_a1 << " (" << goodA1Min << "<a1<" << goodA1Max << ") " << endl; | |
90293ba6 | 836 | |
53f515d3 | 837 | Double_t meanA1 = sumA1/(nGoodChannel_a1); |
838 | Double_t meanProbChi2 = sumProbChi2/(nGoodChannel_a1); | |
839 | Double_t meanA2 = sumA2/(nGoodChannel); | |
840 | Double_t meanProbChi2P2 = sumProbChi2P2/(nGoodChannel); | |
90293ba6 | 841 | |
842 | Double_t capaManu = 0.2; // pF | |
53f515d3 | 843 | filcout << "\n linear fit : <a1> = " << meanA1 << "\t <gain> = " << 1./(meanA1*capaManu) |
90293ba6 | 844 | << " mV/fC (capa= " << capaManu << " pF)" << endl; |
53f515d3 | 845 | filcout << " Prob(chi2)> = " << meanProbChi2 << endl; |
846 | filcout << "\n parabolic fit: <a2> = " << meanA2 << endl; | |
847 | filcout << " Prob(chi2)> = " << meanProbChi2P2 << "\n" << endl; | |
848 | ||
849 | } | |
850 | ||
851 | ||
852 | return ; | |
853 | ||
3f637cfc | 854 | } |
855 | ||
856 | //*************************************************************// | |
857 | ||
858 | // main routine | |
859 | int main(Int_t argc, Char_t **argv) | |
860 | { | |
861 | ||
fd99693f | 862 | // needed for streamer application |
863 | gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo", | |
864 | "*", | |
865 | "TStreamerInfo", | |
866 | "RIO", | |
867 | "TStreamerInfo()"); | |
868 | ||
a76f513d | 869 | TFitter *minuitFit = new TFitter(NFITPARAMS); |
870 | TVirtualFitter::SetFitter(minuitFit); | |
fd99693f | 871 | |
3f637cfc | 872 | Int_t skipEvents = 0; |
873 | Int_t maxEvents = 1000000; | |
90293ba6 | 874 | Int_t MaxDateEvents = 1000000; |
875 | Int_t injCharge = 0; | |
0aa03544 | 876 | Double_t nSigma = 3; |
3f637cfc | 877 | Int_t threshold = -1; |
878 | Char_t inputFile[256]; | |
90293ba6 | 879 | |
53f515d3 | 880 | // TString flatOutputFile; |
881 | TString fesOutputFile; | |
3f637cfc | 882 | TString crocusOutputFile; |
883 | TString crocusConfigFile; | |
884 | ||
3f637cfc | 885 | // option handler |
886 | ||
887 | // decode the input line | |
888 | for (Int_t i = 1; i < argc; i++) // argument 0 is the executable name | |
889 | { | |
890 | Char_t* arg; | |
891 | ||
892 | arg = argv[i]; | |
893 | if (arg[0] != '-') continue; | |
894 | switch (arg[1]) | |
895 | { | |
896 | case 'f' : | |
897 | i++; | |
898 | sprintf(inputFile,argv[i]); | |
899 | break; | |
900 | case 'a' : | |
901 | i++; | |
902 | flatOutputFile = argv[i]; | |
903 | break; | |
53f515d3 | 904 | case 'b' : |
905 | i++; | |
906 | sprintf(gOutFolder,argv[i]); | |
907 | break; | |
3f637cfc | 908 | case 'o' : |
909 | i++; | |
910 | crocusOutputFile = argv[i]; | |
911 | break; | |
912 | case 'c' : | |
913 | i++; | |
914 | crocusConfigFile = argv[i]; | |
915 | break; | |
53f515d3 | 916 | // case 'v' : |
917 | // i++; | |
918 | // logOutputFile = argv[i]; | |
919 | // break; | |
3f637cfc | 920 | case 'e' : |
921 | i++; | |
90293ba6 | 922 | gCommand = argv[i]; |
3f637cfc | 923 | break; |
924 | case 'd' : | |
925 | i++; | |
90293ba6 | 926 | gPrintLevel=atoi(argv[i]); |
927 | break; | |
928 | case 'g' : | |
929 | i++; | |
930 | gPlotLevel=atoi(argv[i]); | |
3f637cfc | 931 | break; |
932 | case 's' : | |
933 | i++; | |
934 | skipEvents=atoi(argv[i]); | |
935 | break; | |
90293ba6 | 936 | case 'l' : |
937 | i++; | |
938 | injCharge=atoi(argv[i]); | |
939 | break; | |
940 | case 'm' : | |
941 | i++; | |
942 | sscanf(argv[i],"%d",&MaxDateEvents); | |
943 | break; | |
3f637cfc | 944 | case 'n' : |
945 | i++; | |
946 | sscanf(argv[i],"%d",&maxEvents); | |
947 | break; | |
948 | case 'p' : | |
949 | i++; | |
0aa03544 | 950 | sscanf(argv[i],"%lf",&nSigma); |
3f637cfc | 951 | break; |
90293ba6 | 952 | case 'r' : |
953 | i++; | |
954 | sprintf(gHistoFileName_gain,argv[i]); | |
955 | break; | |
3f637cfc | 956 | case 't' : |
957 | i++; | |
958 | sscanf(argv[i],"%d",&threshold); | |
959 | break; | |
960 | case 'h' : | |
961 | i++; | |
962 | printf("\n******************* %s usage **********************",argv[0]); | |
963 | printf("\n%s -options, the available options are :",argv[0]); | |
964 | printf("\n-h help (this screen)"); | |
965 | printf("\n"); | |
966 | printf("\n Input"); | |
967 | printf("\n-f <raw data file> (default = %s)",inputFile); | |
968 | printf("\n-c <Crocus config. file> (default = %s)",crocusConfigFile.Data()); | |
969 | printf("\n"); | |
970 | printf("\n Output"); | |
971 | printf("\n-a <Flat ASCII file> (default = %s)",flatOutputFile.Data()); | |
53f515d3 | 972 | printf("\n-o <CROCUS cmd file> (default = %s)",crocusOutputFile.Data()); |
973 | // printf("\n-v <log output file> (default = %s)",logOutputFile.Data()); | |
3f637cfc | 974 | printf("\n"); |
975 | printf("\n Options"); | |
53f515d3 | 976 | printf("\n-b <output directory> (default = %s)",gOutFolder); |
90293ba6 | 977 | printf("\n-d <print level> (default = %d)",gPrintLevel); |
978 | printf("\n-g <plot level> (default = %d)",gPlotLevel); | |
979 | printf("\n-l <DAC level> (default = %d)",injCharge); | |
980 | printf("\n-m <max date events> (default = %d)",MaxDateEvents); | |
3f637cfc | 981 | printf("\n-s <skip events> (default = %d)",skipEvents); |
982 | printf("\n-n <max events> (default = %d)",maxEvents); | |
983 | printf("\n-p <n sigmas> (default = %f)",nSigma); | |
90293ba6 | 984 | printf("\n-r root file data for gain(default = %s)",gHistoFileName_gain); |
3f637cfc | 985 | printf("\n-t <threshold (-1 = no)> (default = %d)",threshold); |
90293ba6 | 986 | printf("\n-e <execute ped/gain> (default = %s)",gCommand.Data()); |
3f637cfc | 987 | printf("\n-e <gain create> make gain & create a new root file"); |
988 | printf("\n-e <gain> make gain & update root file"); | |
989 | printf("\n-e <gain compute> make gain & compute gains"); | |
990 | ||
991 | printf("\n\n"); | |
992 | exit(-1); | |
993 | default : | |
994 | printf("%s : bad argument %s (please check %s -h)\n",argv[0],argv[i],argv[0]); | |
995 | argc = 2; exit(-1); // exit if error | |
996 | } // end of switch | |
997 | } // end of for i | |
998 | ||
90293ba6 | 999 | // set gCommand to lower case |
1000 | gCommand.ToLower(); | |
3f637cfc | 1001 | |
53f515d3 | 1002 | |
3f637cfc | 1003 | // decoding the events |
1004 | ||
1005 | Int_t status; | |
3f637cfc | 1006 | void* event; |
1007 | ||
90293ba6 | 1008 | gPedMeanHisto = 0x0; |
1009 | gPedSigmaHisto = 0x0; | |
3f637cfc | 1010 | |
1011 | TStopwatch timers; | |
1012 | ||
1013 | timers.Start(kTRUE); | |
1014 | ||
1015 | // once we have a configuration file in db | |
1016 | // copy locally a file from daq detector config db | |
1017 | // The current detector is identified by detector code in variable | |
1018 | // DATE_DETECTOR_CODE. It must be defined. | |
1019 | // If environment variable DAQDA_TEST_DIR is defined, files are copied from DAQDA_TEST_DIR | |
1020 | // instead of the database. The usual environment variables are not needed. | |
1021 | if (!crocusConfigFile.IsNull()) { | |
1022 | status = daqDA_DB_getFile("myconfig", crocusConfigFile.Data()); | |
1023 | if (status) { | |
1024 | printf("Failed to get config file : %d\n",status); | |
1025 | return -1; | |
1026 | } | |
1027 | } | |
1028 | ||
1029 | ||
1030 | status = monitorSetDataSource(inputFile); | |
1031 | if (status) { | |
1032 | cerr << "ERROR : monitorSetDataSource status (hex) = " << hex << status | |
1033 | << " " << monitorDecodeError(status) << endl; | |
1034 | return -1; | |
1035 | } | |
1036 | status = monitorDeclareMp("MUON Tracking monitoring"); | |
1037 | if (status) { | |
1038 | cerr << "ERROR : monitorDeclareMp status (hex) = " << hex << status | |
1039 | << " " << monitorDecodeError(status) << endl; | |
1040 | return -1; | |
1041 | } | |
1042 | ||
0aa03544 | 1043 | Int_t busPatchId; |
1044 | UShort_t manuId; | |
1045 | UChar_t channelId; | |
1046 | UShort_t charge; | |
1047 | ||
53f515d3 | 1048 | if (gCommand.CompareTo("comp") != 0) |
1049 | { | |
1050 | cout << "\nMUONTRKda : Reading data from file " << inputFile <<endl; | |
3f637cfc | 1051 | |
53f515d3 | 1052 | while(1) |
1053 | { | |
1054 | if (gNDateEvents >= MaxDateEvents) break; | |
1055 | if (gNEvents >= maxEvents) break; | |
1056 | if (gNEvents && gNEvents % 100 == 0) | |
1057 | cout<<"Cumulated events " << gNEvents << endl; | |
1058 | ||
1059 | // check shutdown condition | |
1060 | if (daqDA_checkShutdown()) | |
1061 | break; | |
1062 | ||
1063 | // Skip Events if needed | |
1064 | while (skipEvents) { | |
1065 | status = monitorGetEventDynamic(&event); | |
1066 | skipEvents--; | |
1067 | } | |
3f637cfc | 1068 | |
53f515d3 | 1069 | // starts reading |
1070 | status = monitorGetEventDynamic(&event); | |
1071 | if (status < 0) { | |
1072 | cout<<"EOF found"<<endl; | |
1073 | break; | |
1074 | } | |
1075 | ||
1076 | // gNDateEvents++; | |
3f637cfc | 1077 | |
53f515d3 | 1078 | // decoding rawdata headers |
1079 | AliRawReader *rawReader = new AliRawReaderDate(event); | |
3f637cfc | 1080 | |
53f515d3 | 1081 | Int_t eventType = rawReader->GetType(); |
1082 | gRunNumber = rawReader->GetRunNumber(); | |
1083 | ||
1084 | // Output log file initialisations | |
1085 | ||
1086 | if(gNDateEvents==0) | |
1087 | { | |
1088 | if (gCommand.CompareTo("ped") == 0){ | |
1089 | sprintf(flatFile,"%s/MUONTRKda_ped_%d.log",gOutFolder,gRunNumber); | |
1090 | logOutputFile=flatFile; | |
1091 | ||
1092 | filcout.open(logOutputFile.Data()); | |
1093 | filcout<<"//=================================================" << endl; | |
1094 | filcout<<"// MUONTRKda for Pedestal run = " << gRunNumber << endl; | |
1095 | cout<<"\n ******** MUONTRKda for Pedestal run = " << gRunNumber << "\n" << endl; | |
1096 | } | |
1097 | ||
1098 | if (gCommand.Contains("gain")){ | |
1099 | sprintf(flatFile,"%s/%s_%d_DAC_%d.log",gOutFolder,filenam,gRunNumber,injCharge); | |
1100 | logOutputFile=flatFile; | |
1101 | ||
1102 | filcout.open(logOutputFile.Data()); | |
1103 | filcout<<"//=================================================" << endl; | |
1104 | filcout<<"// MUONTRKda for Gain run = " << gRunNumber << " (DAC=" << injCharge << ")" << endl; | |
1105 | cout<<"\n ******** MUONTRKda for Gain run = " << gRunNumber << " (DAC=" << injCharge << ")\n" << endl; | |
1106 | } | |
1107 | ||
1108 | filcout<<"//=================================================" << endl; | |
1109 | filcout<<"// * Date : " << date.AsString("l") << "\n" << endl; | |
1110 | cout<<" * Date : " << date.AsString("l") << "\n" << endl; | |
1111 | ||
1112 | } | |
1113 | ||
1114 | gNDateEvents++; | |
1115 | ||
3f637cfc | 1116 | |
3f637cfc | 1117 | |
53f515d3 | 1118 | if (eventType != PHYSICS_EVENT) |
1119 | continue; // for the moment | |
3f637cfc | 1120 | |
53f515d3 | 1121 | // decoding MUON payload |
1122 | AliMUONRawStreamTracker* rawStream = new AliMUONRawStreamTracker(rawReader); | |
3f637cfc | 1123 | |
53f515d3 | 1124 | // loops over DDL |
1125 | rawStream->First(); // if GlitchError ? what we are doing ? | |
1126 | while( (status = rawStream->Next(busPatchId, manuId, channelId, charge)) ) { | |
0aa03544 | 1127 | |
53f515d3 | 1128 | if (gNEvents == 0) |
1129 | gNChannel++; | |
0aa03544 | 1130 | |
53f515d3 | 1131 | // if (gPrintLevel) printf("manuId: %d, channelId: %d charge: %d\n", manuId, |
1132 | // channelId, charge); | |
0aa03544 | 1133 | |
53f515d3 | 1134 | MakePed(busPatchId, (Int_t)manuId, (Int_t)channelId, (Int_t)charge); |
3f637cfc | 1135 | |
53f515d3 | 1136 | } // Next digit |
1137 | ||
1138 | // if parity error do not count event // JLC (modified: Christian 11/10/07) | |
1139 | if (!rawStream->GetPayLoad()->GetParityErrBus().GetSize() && !rawStream->GetPayLoad()->GetGlitchErrors() | |
1140 | && !rawStream->GetPayLoad()->GetPaddingErrors()) { | |
1141 | gNEvents++; | |
1142 | } | |
1143 | if (rawStream->GetPayLoad()->GetGlitchErrors())gGlitchErrors++; | |
1144 | if (rawStream->GetPayLoad()->GetParityErrBus().GetSize())gParityErrors++; | |
1145 | if (rawStream->GetPayLoad()->GetPaddingErrors())gPaddingErrors++; | |
1146 | ||
1147 | delete rawReader; | |
1148 | delete rawStream; | |
3f637cfc | 1149 | |
53f515d3 | 1150 | } // while (1) |
1151 | } | |
3f637cfc | 1152 | |
3f637cfc | 1153 | |
90293ba6 | 1154 | |
1155 | if (gCommand.CompareTo("ped") == 0) | |
1156 | { | |
53f515d3 | 1157 | sprintf(flatFile,"MUONTRKda_ped_%d.ped",gRunNumber); |
90293ba6 | 1158 | if(flatOutputFile.IsNull())flatOutputFile=flatFile; |
1159 | MakePedStore(flatOutputFile); | |
1160 | } | |
3f637cfc | 1161 | |
1162 | // option gain -> update root file with pedestal results | |
1163 | // gain + create -> recreate root file | |
1164 | // gain + comp -> update root file and compute gain parameters | |
1165 | ||
53f515d3 | 1166 | if (gCommand.Contains("gain")) |
1167 | { | |
1168 | MakePedStoreForGain(injCharge); | |
1169 | } | |
3f637cfc | 1170 | |
53f515d3 | 1171 | if (gCommand.Contains("comp")) |
1172 | { | |
1173 | ||
1174 | // if(flatOutputFile.IsNull())flatOutputFile="MUONTRKda_gain.par"; | |
1175 | // MakeGainStore(flatOutputFile); | |
1176 | MakeGainStore(); | |
1177 | } | |
3f637cfc | 1178 | |
1179 | ||
90293ba6 | 1180 | delete gPedestalStore; |
3f637cfc | 1181 | |
a76f513d | 1182 | delete minuitFit; |
1183 | TVirtualFitter::SetFitter(0); | |
1184 | ||
3f637cfc | 1185 | timers.Stop(); |
1186 | ||
90293ba6 | 1187 | if (gCommand.CompareTo("comp") != 0) |
1188 | { | |
1189 | cout << "MUONTRKda : Nb of DATE events = " << gNDateEvents << endl; | |
53f515d3 | 1190 | cout << "MUONTRKda : Nb of Glitch errors = " << gGlitchErrors << endl; |
1191 | cout << "MUONTRKda : Nb of Parity errors = " << gParityErrors << endl; | |
1192 | cout << "MUONTRKda : Nb of Padding errors = " << gPaddingErrors << endl; | |
90293ba6 | 1193 | cout << "MUONTRKda : Nb of events used = " << gNEvents << endl; |
53f515d3 | 1194 | |
1195 | filcout << "MUONTRKda : Nb of DATE events = " << gNDateEvents << endl; | |
1196 | filcout << "MUONTRKda : Nb of Glitch errors = " << gGlitchErrors << endl; | |
1197 | filcout << "MUONTRKda : Nb of Parity errors = " << gParityErrors << endl; | |
1198 | filcout << "MUONTRKda : Nb of Padding errors = " << gPaddingErrors << endl; | |
1199 | filcout << "MUONTRKda : Nb of events used = " << gNEvents << endl; | |
1200 | ||
90293ba6 | 1201 | } |
53f515d3 | 1202 | |
1203 | ||
1204 | cout << "\nMUONTRKda : Output file generated : " << logOutputFile << endl; | |
1205 | ||
90293ba6 | 1206 | if (gCommand.CompareTo("ped") == 0) |
1207 | { | |
1208 | if (!(crocusConfigFile.IsNull())) | |
1209 | cout << "MUONTRKda : CROCUS command file generated : " << crocusOutputFile.Data() << endl; | |
1210 | else | |
1211 | cout << "MUONTRKda : WARNING no CROCUS command file generated" << endl; | |
1212 | cout << "\nMUONTRKda : Histo file generated for pedestal : " << gHistoFileName << endl; | |
1213 | } | |
3f637cfc | 1214 | else |
90293ba6 | 1215 | { |
53f515d3 | 1216 | cout << "\nMUONTRKda : Root data file for gain calculation : " << gHistoFileName_gain << endl; |
1217 | } | |
1218 | ||
1219 | if (gCommand.CompareTo("comp") == 0) | |
1220 | { | |
1221 | cout << "MUONTRKda : Root Histo. file generated : " << gRootFileName << endl; | |
1222 | cout << "MUONTRKda : Flat ASCII file generated (to SHUTTLE) : " << flatOutputFile << endl; | |
1223 | } | |
1224 | ||
1225 | ||
1226 | ||
1227 | // Store IN FES | |
1228 | ||
1229 | if (gCommand.CompareTo("comp") == 0 || gCommand.CompareTo("ped") == 0) | |
1230 | { | |
1231 | printf("\n ***** STORE FILE in FES ****** \n"); | |
1232 | ||
1233 | // to be sure that env variable is set | |
1234 | gSystem->Setenv("DAQDALIB_PATH", "$DATE_SITE/infoLogger"); | |
1235 | ||
1236 | if (!flatOutputFile.IsNull()) | |
1237 | { | |
1238 | ||
1239 | // flatOutputFile.Prepend("./"); | |
1240 | ||
1241 | status = daqDA_FES_storeFile(flatOutputFile.Data(),flatOutputFile.Data()); | |
1242 | if (status) | |
1243 | { | |
1244 | printf(" Failed to export file : %d\n",status); | |
1245 | } | |
1246 | else if(gPrintLevel) printf("Export file: %s\n",flatOutputFile.Data()); | |
1247 | } | |
90293ba6 | 1248 | } |
3f637cfc | 1249 | |
53f515d3 | 1250 | filcout.close(); |
90293ba6 | 1251 | printf("\nExecution time : R:%7.2fs C:%7.2fs\n", timers.RealTime(), timers.CpuTime()); |
3f637cfc | 1252 | |
1253 | return status; | |
1254 | } |