]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/MUONTRKda.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Version number incremented
[u/mrichter/AliRoot.git] / MUON / MUONTRKda.cxx
CommitLineData
830e3f26 1/*
2Contact: Jean-Luc Charvet <jean-luc.charvet@cea.fr>
3Link: http://aliceinfo.cern.ch/static/Offline/dimuon/muon_html/README_Mchda
4Run Type: PEDESTAL, CALIBRATION
5DA Type: LDC
6Number of events needed: 400 events for pedestal and each calibration run
7Input Files: Rawdata file (DATE format)
8Output Files: local dir (not persistent) -> MUONTRKda_ped_<run#>.ped , MUONTRKda_gain_<run#>.par
9 FXS -> run<#>_MCH_<ldc>_PEDESTALS, run<#>_MCH_<ldc>_GAINS
10Trigger 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 35Version for MUONTRKda MUON tracking
36(A. Baldisseri, J.-L. Charvet & Ch. Finck)
3f637cfc 37
38
39Rem: AliMUON2DMap stores all channels, even those which are not connected
40 if pedMean == -1, channel not connected to a pad
41
3f637cfc 42
43*/
44extern "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 88const Int_t gkNChannels = AliMpConstants::ManuNofChannels();
89const Int_t gkADCMax = 4095;
90
91AliMUONVStore* gPedestalStore = new AliMUON2DMap(kFALSE);
92
93Int_t gNManu = 0;
94Int_t gNChannel = 0;
95UInt_t gRunNumber = 0;
96Int_t gNEvents = 0;
97Int_t gNDateEvents = 0;
666b5407 98Int_t gPrintLevel = 1; // global printout variable (others: 2 and 3)
90293ba6 99Int_t gPlotLevel = 0; // global plot variable
830e3f26 100Int_t gFES = 1; // by default FES is used
90293ba6 101
102TH1F* gPedMeanHisto = 0x0;
103TH1F* gPedSigmaHisto = 0x0;
104Char_t gHistoFileName[256];
105
81028269 106// used for computing gain parameters
107Int_t nbpf1 = 6; // linear fit over nbf1 points
108
109//Char_t gHistoFileName_gain[256]="MUONTRKda_gain_data.root";
110Char_t gHistoFileName_gain[256]="MUONTRKda_gain.data";
53f515d3 111Char_t gRootFileName[256];
112Char_t gOutFolder[256]=".";
113Char_t filename[256];
114Char_t filenam[256]="MUONTRKda_gain";
81028269 115Char_t flatFile[256]="";
53f515d3 116
81028269 117//ofstream filcout;
53f515d3 118
119TString flatOutputFile;
120TString logOutputFile;
81028269 121TString logOutputFile_comp;
90293ba6 122TString gCommand("ped");
53f515d3 123TTimeStamp date;
3f637cfc 124
125// funtions
126
90293ba6 127
128//________________
129Double_t funcLin(Double_t *x, Double_t *par)
130{
131 return par[0] + par[1]*x[0];
132}
133
3f637cfc 134//________________
90293ba6 135Double_t funcParabolic(Double_t *x, Double_t *par)
3f637cfc 136{
90293ba6 137 return par[0]*x[0]*x[0];
138}
3f637cfc 139
90293ba6 140//________________
141Double_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//__________
153void 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//________________
179void 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 298void 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)
373void 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
920int 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}