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