1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
19 // EMCAL Preprocessor class. It runs by Shuttle at the end of the run,
20 // calculates stuff to be posted in OCDB
22 // Author: Boris Polichtchouk, 4 October 2006
23 // Adapted for EMCAL by Gustavo Conesa Balbastre, October 2006
24 // Updated by David Silvermyr May 2008, based on TPC code
25 ///////////////////////////////////////////////////////////////////////////////
31 #include "TParameter.h"
33 #include <TTimeStamp.h>
36 #include "AliShuttleInterface.h"
37 #include "AliEMCALPreprocessor.h"
39 #include "AliDCSValue.h"
40 #include "AliCDBMetaData.h"
41 #include "AliEMCALTriggerTRUDCSConfig.h"
42 #include "AliEMCALTriggerDCSConfig.h"
43 #include "AliCaloCalibPedestal.h"
44 #include "AliCaloCalibSignal.h"
45 #include "AliEMCALSensorTempArray.h"
47 const Int_t kValCutTemp = 100; // discard temperatures > 100 degrees
48 const Int_t kDiffCutTemp = 5; // discard temperature differences > 5 degrees
49 const TString kPedestalRunType = "PEDESTAL"; // pedestal run identifier
50 const TString kPhysicsRunType = "PHYSICS"; // physics run identifier
51 const TString kStandAloneRunType = "STANDALONE_BC"; // standalone run identifier
52 const TString kAmandaTemp = "PT_%02d.Temperature"; // Amanda string for temperature entries
53 //const Double_t kFitFraction = 0.7; // Fraction of DCS sensor fits required
54 const Double_t kFitFraction = -1.0; // Don't require minimum number of fits during commissioning
56 const TString kMetaResponsible = "David Silvermyr";
57 //legacy comments and return codes from TPC
58 const TString kMetaComment = "Preprocessor AliEMCAL data base entries.";
59 const int kReturnCodeNoInfo = 9;
60 const int kReturnCodeNoObject = 2;
61 const int kReturnCodeNoEntries = 1;
63 const int kNTRU = 30; // From 2011; 10 SuperModules (SM) * 3 TRU per SM
65 ClassImp(AliEMCALPreprocessor)
67 //_______________________________________________________________________________________
68 AliEMCALPreprocessor::AliEMCALPreprocessor() :
69 AliPreprocessor("EMC",0),
77 //_______________________________________________________________________________________
78 AliEMCALPreprocessor::AliEMCALPreprocessor(AliShuttleInterface* shuttle):
79 AliPreprocessor("EMC",shuttle),
84 // Constructor AddRunType(kPedestalRunType);
86 // define run types to be processed
87 AddRunType(kPedestalRunType);
88 AddRunType(kPhysicsRunType);
91 //______________________________________________________________________________________________
92 AliEMCALPreprocessor::AliEMCALPreprocessor(const AliEMCALPreprocessor& ) :
93 AliPreprocessor("EMCAL",0),
94 fConfEnv(0), fTemp(0), fConfigOK(kTRUE)
96 Fatal("AliEMCALPreprocessor", "copy constructor not implemented");
99 // assignment operator; use copy ctor to make life easy.
100 //______________________________________________________________________________________________
101 AliEMCALPreprocessor& AliEMCALPreprocessor::operator = (const AliEMCALPreprocessor &source )
103 // assignment operator; use copy ctor
104 if (&source == this) return *this;
106 new (this) AliEMCALPreprocessor(source);
110 //____________________________________________________________________________
111 AliEMCALPreprocessor::~AliEMCALPreprocessor()
114 if (fTemp) delete fTemp;
117 //______________________________________________________________________________________________
118 void AliEMCALPreprocessor::Initialize(Int_t run, UInt_t startTime,
121 // Creates AliTestDataDCS object -- start maps half an hour beforre actual run start
122 UInt_t startTimeLocal = startTime-1800;
123 AliPreprocessor::Initialize(run, startTimeLocal, endTime);
125 AliInfo(Form("\n\tRun %d \n\tStartTime %s \n\tEndTime %s", run,
126 TTimeStamp((time_t)startTime,0).AsString(),
127 TTimeStamp((time_t)endTime,0).AsString()));
129 // Preprocessor configuration
130 AliCDBEntry* entry = GetFromOCDB("Config", "Preprocessor");
131 if (entry) fConfEnv = (TEnv*) entry->GetObject();
133 Log("AliEMCALPreprocessor: Preprocessor Config OCDB entry missing.\n");
138 // Temperature sensors
141 TString tempConf = fConfEnv->GetValue("Temperature","ON");
143 if (tempConf != "OFF" ) {
144 entry = GetFromOCDB("Config", "Temperature");
145 if (entry) confTree = (TTree*) entry->GetObject();
147 Log("AliEMCALPreprocessor: Temperature Config OCDB entry missing.\n");
151 fTemp = new AliEMCALSensorTempArray(startTimeLocal, fEndTime, confTree, kAmandaTemp);
152 fTemp->SetValCut(kValCutTemp);
153 fTemp->SetDiffCut(kDiffCutTemp);
159 //______________________________________________________________________________________________
160 UInt_t AliEMCALPreprocessor::Process(TMap* dcsAliasMap)
162 // Fills data into EMCAL calibrations objects
163 // Amanda servers provide information directly through dcsAliasMap
165 if (!fConfigOK) return kReturnCodeNoInfo;
167 TObjArray *resultArray = new TObjArray();
168 TString errorHandling = fConfEnv->GetValue("ErrorHandling","ON");
169 errorHandling.ToUpper();
173 if (errorHandling == "OFF" ) {
174 if (!dcsAliasMap) dcsResult = kReturnCodeNoEntries;
175 else if (dcsAliasMap->GetEntries() == 0 ) dcsResult = kReturnCodeNoEntries;
176 status = new TParameter<int>("dcsResult",dcsResult);
177 resultArray->Add(status);
180 if (!dcsAliasMap) return kReturnCodeNoInfo;
181 else if (dcsAliasMap->GetEntries() == 0 ) return kReturnCodeNoInfo;
185 TString runType = GetRunType();
187 // Temperature sensors are processed by AliEMCALCalTemp
188 TString tempConf = fConfEnv->GetValue("Temperature","ON");
190 if (tempConf != "OFF" && dcsAliasMap ) {
191 UInt_t tempResult = MapTemperature(dcsAliasMap);
193 status = new TParameter<int>("tempResult",tempResult);
194 resultArray->Add(status);
196 // Trigger configuration processing: only for Physics runs
197 TString triggerConf = fConfEnv->GetValue("Trigger","ON");
198 triggerConf.ToUpper();
199 if( runType == kPhysicsRunType ) {
200 if (triggerConf != "OFF" && dcsAliasMap ) {
201 UInt_t triggerResult = MapTriggerConfig(dcsAliasMap);
202 result+=triggerResult;
203 status = new TParameter<int>("triggerResult",triggerResult);
204 resultArray->Add(status);
208 // Other calibration information will be retrieved through FXS files
210 // TList* fileSourcesDAQ = GetFile(AliShuttleInterface::kDAQ, "pedestals");
211 // const char* fileNamePed = GetFile(AliShuttleInterface::kDAQ, "pedestals", "LDC1");
213 // TList* fileSourcesHLT = GetFile(AliShuttleInterface::kHLT, "calib");
214 // const char* fileNameHLT = GetFile(AliShuttleInterface::kHLT, "calib", "LDC1");
218 if ( runType == kPedestalRunType ) {
219 Int_t numSources = 1;
220 Int_t pedestalSource[2] = {AliShuttleInterface::kDAQ, AliShuttleInterface::kHLT} ;
221 TString source = fConfEnv->GetValue("Pedestal","DAQ");
223 if (source != "OFF" ) {
224 if ( source == "HLT") pedestalSource[0] = AliShuttleInterface::kHLT;
225 if (!GetHLTStatus()) pedestalSource[0] = AliShuttleInterface::kDAQ;
226 if (source == "HLTDAQ" ) {
228 pedestalSource[0] = AliShuttleInterface::kHLT;
229 pedestalSource[1] = AliShuttleInterface::kDAQ;
231 if (source == "DAQHLT" ) numSources=2;
232 UInt_t pedestalResult=0;
233 for (Int_t i=0; i<numSources; i++ ) {
234 pedestalResult = ExtractPedestals(pedestalSource[i]);
235 if ( pedestalResult == 0 ) break;
237 result += pedestalResult;
238 status = new TParameter<int>("pedestalResult",pedestalResult);
239 resultArray->Add(status);
243 // SIGNAL/LED ENTRIES:
244 if( runType == kPhysicsRunType ) {
245 Int_t numSources = 1;
246 Int_t signalSource[2] = {AliShuttleInterface::kDAQ,AliShuttleInterface::kHLT} ;
247 TString source = fConfEnv->GetValue("Signal","DAQ");
249 if ( source != "OFF") {
250 if ( source == "HLT") signalSource[0] = AliShuttleInterface::kHLT;
251 if (!GetHLTStatus()) signalSource[0] = AliShuttleInterface::kDAQ;
252 if (source == "HLTDAQ" ) {
254 signalSource[0] = AliShuttleInterface::kHLT;
255 signalSource[1] = AliShuttleInterface::kDAQ;
257 if (source == "DAQHLT" ) numSources=2;
258 UInt_t signalResult=0;
259 for (Int_t i=0; i<numSources; i++ ) {
260 signalResult = ExtractSignal(signalSource[i]);
261 if ( signalResult == 0 ) break;
263 result += signalResult;
264 status = new TParameter<int>("signalResult",signalResult);
265 resultArray->Add(status);
270 // overall status at the end
271 if (errorHandling == "OFF" ) {
272 AliCDBMetaData metaData;
273 metaData.SetBeamPeriod(0);
274 metaData.SetResponsible(kMetaResponsible);
275 metaData.SetComment("Preprocessor AliEMCAL status.");
276 Bool_t storeOK = Store("Calib", "PreprocStatus", resultArray, &metaData, 0, kFALSE);
277 resultArray->Delete();
279 if ( !storeOK ) result=1;
287 //______________________________________________________________________________________________
288 UInt_t AliEMCALPreprocessor::MapTemperature(TMap* dcsAliasMap)
289 { // extract DCS temperature maps. Perform fits to save space
292 TMap *map = fTemp->ExtractDCS(dcsAliasMap);
294 fTemp->MakeSplineFit(map);
295 Double_t fitFraction = 1.0*fTemp->NumFits()/fTemp->NumSensors();
296 if (fitFraction > kFitFraction ) {
297 AliInfo(Form("Temperature values extracted, fits performed.\n"));
300 Log ("Too few temperature maps fitted. \n");
301 result = kReturnCodeNoInfo;
305 Log("No temperature map extracted. \n");
306 result = kReturnCodeNoInfo;
309 // Now store the final CDB file
311 if ( result == 0 ) { // some info was found
312 AliCDBMetaData metaData;
313 metaData.SetBeamPeriod(0);
314 metaData.SetResponsible(kMetaResponsible);
315 metaData.SetComment(kMetaComment);
317 Bool_t storeOK = Store("Calib", "Temperature", fTemp, &metaData, 0, kFALSE);
318 if ( !storeOK ) result=1;
319 AliInfo(Form("Temperature info stored. result %d\n", result));
325 //______________________________________________________________________________________________
326 UInt_t AliEMCALPreprocessor::MapTriggerConfig(TMap* dcsAliasMap)
327 { // extract DCS trigger info
328 AliInfo(Form("Get TRU info from DCS DPs.\n"));
330 const Int_t bufsize = 1000;
334 TObjArray *arrL0ALGSEL, *arrPEAKFINDER, *arrGLOBALTHRESH, *arrCOSMTHRESH;
335 TObjArray *arrMASK[6];
337 // overall object to hold STU and DCS config info
338 // DS comment: for now only holds TRU info, i.e. only partially filled
339 // (STU info only in raw data header; unfortunately not also picked up via DCS DPs)
340 AliEMCALTriggerDCSConfig *trigConfig = new AliEMCALTriggerDCSConfig();
341 // allocate space for TRU objects
342 TClonesArray *truArr = new TClonesArray("AliEMCALTriggerTRUDCSConfig", kNTRU);
343 for( iTRU = 0; iTRU < kNTRU; iTRU++){
344 new((*truArr)[iTRU]) AliEMCALTriggerTRUDCSConfig();
346 trigConfig->SetTRUArr(truArr);
348 // loop through all TRUs
349 bool debug = false; // debug flag for AliInfo printouts for each TRU
350 for( iTRU = 0; iTRU < kNTRU; iTRU++){
351 if (debug) AliInfo( Form("iTRU %d \n", iTRU) );
352 // get the shuttled values
353 snprintf( buf, bufsize, "EMC_TRU%02d_L0ALGSEL", iTRU );
354 arrL0ALGSEL = (TObjArray*) dcsAliasMap->GetValue( buf );
355 snprintf( buf, bufsize, "EMC_TRU%02d_PEAKFINDER", iTRU );
356 arrPEAKFINDER = (TObjArray*) dcsAliasMap->GetValue( buf );
357 snprintf( buf, bufsize, "EMC_TRU%02d_GLOBALTHRESH", iTRU );
358 arrGLOBALTHRESH = (TObjArray*) dcsAliasMap->GetValue( buf );
359 snprintf( buf, bufsize, "EMC_TRU%02d_COSMTHRESH", iTRU );
360 arrCOSMTHRESH = (TObjArray*) dcsAliasMap->GetValue( buf );
362 for( i = 0; i < 6; i++ ){
363 snprintf( buf, bufsize, "EMC_TRU%02d_MASK%d", iTRU, i );
364 arrMASK[i] = (TObjArray*) dcsAliasMap->GetValue( buf );
368 AliEMCALTriggerTRUDCSConfig* truConfig = trigConfig->GetTRUDCSConfig(iTRU);
370 AliWarning( Form("EMC DCS TRU%02d config not retrieved!\n", iTRU ));
374 // get last entries. fill the TRU object
376 AliWarning( Form("EMC DCS TRU%02d L0ALGSEL alias not found!\n", iTRU ));
379 if (debug) AliInfo( Form("arrL0ALGSEL has %d entries \n", arrL0ALGSEL->GetEntries()) );
380 if ( arrL0ALGSEL->GetEntries() > 0 ) {
381 dcsVal = (AliDCSValue *) arrL0ALGSEL->At( arrL0ALGSEL->GetEntries() - 1 );
382 if (dcsVal) truConfig->SetL0SEL( dcsVal->GetUInt() );
386 if( ! arrPEAKFINDER ){
387 AliWarning( Form("EMC DCS TRU%02d PEAKFINDER alias not found!\n", iTRU ));
390 if (debug) AliInfo( Form("arrPEAKFINDER has %d entries \n", arrPEAKFINDER->GetEntries()) );
391 if ( arrPEAKFINDER->GetEntries() > 0 ) {
392 dcsVal = (AliDCSValue *) arrPEAKFINDER->At( arrPEAKFINDER->GetEntries() - 1 );
393 if (dcsVal) truConfig->SetSELPF( dcsVal->GetUInt() );
397 if( ! arrGLOBALTHRESH ){
398 AliWarning( Form("EMC DCS TRU%02d GLOBALTHRESH alias not found!\n", iTRU ));
401 if (debug) AliInfo( Form("arrGLOBALTHRESH has %d entries \n", arrGLOBALTHRESH->GetEntries()) );
402 if ( arrGLOBALTHRESH->GetEntries() > 0 ) {
403 dcsVal = (AliDCSValue *) arrGLOBALTHRESH->At( arrGLOBALTHRESH->GetEntries() - 1 );
404 if (dcsVal) truConfig->SetGTHRL0( dcsVal->GetUInt() );
408 if( ! arrCOSMTHRESH ){
409 AliWarning( Form("EMC DCS TRU%02d COSMTHRESH alias not found!\n", iTRU ));
412 if (debug) AliInfo( Form("arrCOSMTHRESH has %d entries \n", arrCOSMTHRESH->GetEntries()) );
413 if ( arrCOSMTHRESH->GetEntries() > 0 ) {
414 dcsVal = (AliDCSValue *) arrCOSMTHRESH->At( arrCOSMTHRESH->GetEntries() - 1 );
415 if (dcsVal) truConfig->SetL0COSM( dcsVal->GetUInt() );
419 for( i = 0; i < 6; i++ ){
421 AliWarning( Form("EMC DCS TRU%02d MASK%d alias not found!\n", iTRU, i ));
424 if (debug) AliInfo( Form("arrMASK[%d] has %d entries \n", i, arrMASK[i]->GetEntries()) );
425 if ( arrMASK[i]->GetEntries() > 0 ) {
426 dcsVal = (AliDCSValue *) arrMASK[i]->At( arrMASK[i]->GetEntries() - 1 );
427 if (dcsVal) truConfig->SetMaskReg( dcsVal->GetUInt(), i );
433 AliInfo(Form("TRU info retrieved.\n"));
435 AliCDBMetaData metaData;
436 metaData.SetBeamPeriod(0);
437 metaData.SetResponsible(kMetaResponsible);
438 metaData.SetComment(kMetaComment);
441 Bool_t storeOK = Store("Calib", "Trigger", trigConfig, &metaData, 0, kFALSE);
442 if ( !storeOK ) result=1;
443 AliInfo(Form("TRU info stored. result %d\n", result));
448 //______________________________________________________________________________________________
449 UInt_t AliEMCALPreprocessor::ExtractPedestals(Int_t sourceFXS)
451 // Read pedestal file from file exchange server
452 // Only store if new pedestal info is available
456 AliCaloCalibPedestal *calibPed = new AliCaloCalibPedestal(AliCaloCalibPedestal::kEmCal);
459 TList* list = GetFileSources(sourceFXS,"pedestals");
460 if (list && list->GetEntries()>0) {
462 // loop through all files from LDCs
466 while (list->At(index)!=NULL) {
467 TObjString* fileNameEntry = (TObjString*) list->At(index);
468 if (fileNameEntry!=NULL) {
469 TString fileName = GetFile(sourceFXS, "pedestals",
470 fileNameEntry->GetString().Data());
471 TFile *f = TFile::Open(fileName);
473 Log ("Error opening pedestal file.");
474 result = kReturnCodeNoObject;
477 AliCaloCalibPedestal *calPed;
478 f->GetObject("emcCalibPedestal",calPed);
480 Log ("No pedestal calibration object in file.");
481 result = kReturnCodeNoObject;
484 if ( calPed->GetNEvents()>0 && calPed->GetNChanFills()>0 ) {
485 // add info for the modules available in the present file
486 Bool_t status = calibPed->AddInfo(calPed);
487 if (status) { changes++; }
497 // Store updated pedestal entry to OCDB
500 AliCDBMetaData metaData;
501 metaData.SetBeamPeriod(0);
502 metaData.SetResponsible(kMetaResponsible);
503 metaData.SetComment(kMetaComment);
505 Bool_t storeOK = StoreReferenceData("Calib", "Pedestals", calibPed, &metaData);
506 if ( !storeOK ) result++;
510 Log ("Error: no entries in input file list!");
511 result = kReturnCodeNoEntries;
517 //______________________________________________________________________________________________
518 UInt_t AliEMCALPreprocessor::ExtractSignal(Int_t sourceFXS)
519 { // Read signal file from file exchange server
520 // Only store if new signal info is available
523 AliCaloCalibSignal *calibSig = new AliCaloCalibSignal(AliCaloCalibSignal::kEmCal);
525 TList* list = GetFileSources(sourceFXS,"signal");
526 if (list && list->GetEntries()>0) {
528 // loop through all files from LDCs
532 while (list->At(index)!=NULL) {
533 TObjString* fileNameEntry = (TObjString*) list->At(index);
534 if (fileNameEntry!=NULL) {
535 TString fileName = GetFile(sourceFXS, "signal",
536 fileNameEntry->GetString().Data());
537 TFile *f = TFile::Open(fileName);
539 Log ("Error opening signal file.");
540 result = kReturnCodeNoObject;
543 AliCaloCalibSignal *calSig;
544 f->GetObject("emcCalibSignal",calSig);
546 Log ("No signal calibration object in file.");
547 result = kReturnCodeNoObject;
550 if ( calSig->GetNEvents()>0 ) {
551 // add info for the modules available in the present file
552 Bool_t status = calibSig->AddInfo(calSig);
553 if (status) { changes++; }
563 // Store updated signal entry to OCDB
566 AliCDBMetaData metaData;
567 metaData.SetBeamPeriod(0);
568 metaData.SetResponsible(kMetaResponsible);
569 metaData.SetComment(kMetaComment);
571 Bool_t storeOK = Store("Calib", "LED", calibSig, &metaData, 0, kFALSE);
572 if ( !storeOK ) result++;
576 Log ("Error: no entries in input file list!");
577 result = kReturnCodeNoEntries;