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