2 //**************************************************************************
3 //* This file is property of and copyright by the ALICE HLT Project *
4 //* ALICE Experiment at CERN, All rights reserved. *
6 //* Primary Authors: Kalliopi Kanaki <Kalliopi.Kanaki@ift.uib.no> *
7 //* for The ALICE HLT Project. *
9 //* Permission to use, copy, modify and distribute this software and its *
10 //* documentation strictly for non-commercial purposes is hereby granted *
11 //* without fee, provided that the above copyright notice appears in all *
12 //* copies and that both the copyright notice and this permission notice *
13 //* appear in the supporting documentation. The authors make no claims *
14 //* about the suitability of this software for any purpose. It is *
15 //* provided "as is" without express or implied warranty. *
16 //**************************************************************************
18 /** @file AliHLTTPCHistogramHandlerComponent.cxx
19 @author Kalliopi Kanaki
21 @brief The Histogram Handler component
27 #include "AliHLTTPCHistogramHandlerComponent.h"
28 #include "AliHLTTPCDefinitions.h"
29 #include "AliCDBEntry.h"
30 #include "AliCDBManager.h"
31 #include "AliHLTTPCTransform.h"
37 #include "TObjArray.h"
38 #include "TObjString.h"
45 ClassImp(AliHLTTPCHistogramHandlerComponent) //ROOT macro for the implementation of ROOT specific class methods
47 AliHLTTPCHistogramHandlerComponent::AliHLTTPCHistogramHandlerComponent()
51 fKryptonHistograms(0),
53 fIgnoreSpecification(kFALSE),
57 fTotalClusterChargeIROCAll(NULL),
58 fTotalClusterChargeOROCAll(NULL),
59 fQMaxPartitionAll(NULL),
60 fPlotQmaxROCAll(NULL),
61 fNumberOfClusters(NULL),
64 fHistTPCSideAmax(NULL),
65 fHistTPCSideCmax(NULL),
66 fHistTPCSideAtot(NULL),
67 fHistTPCSideCtot(NULL),
68 fHistTPCSideArms(NULL),
69 fHistTPCSideCrms(NULL),
73 // see header file for class documentation
75 // refer to README to build package
77 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
80 AliHLTTPCHistogramHandlerComponent::~AliHLTTPCHistogramHandlerComponent() {
81 // see header file for class documentation
85 // Public functions to implement AliHLTComponent's interface.
86 // These functions are required for the registration process
88 const char* AliHLTTPCHistogramHandlerComponent::GetComponentID() {
89 // see header file for class documentation
91 return "TPCHistogramHandler";
94 void AliHLTTPCHistogramHandlerComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list) {
95 // see header file for class documentation
98 list.push_back( kAliHLTDataTypeHistogram );
101 AliHLTComponentDataType AliHLTTPCHistogramHandlerComponent::GetOutputDataType() {
102 // see header file for class documentation
104 return kAliHLTDataTypeHistogram;
107 int AliHLTTPCHistogramHandlerComponent::GetOutputDataTypes(AliHLTComponentDataTypeList& tgtList) {
108 // see header file for class documentation
111 tgtList.push_back(kAliHLTDataTypeHistogram);
112 return tgtList.size();
115 void AliHLTTPCHistogramHandlerComponent::GetOutputDataSize( unsigned long& constBase, double& inputMultiplier ) {
116 // see header file for class documentation
122 AliHLTComponent* AliHLTTPCHistogramHandlerComponent::Spawn() {
123 // see header file for class documentation
125 return new AliHLTTPCHistogramHandlerComponent();
128 int AliHLTTPCHistogramHandlerComponent::DoInit( int argc, const char** argv ) {
129 // see header file for class documentation
136 TString configuration="";
138 for (int j=0; j<argc && iResult>=0; j++) {
141 if (!configuration.IsNull()) configuration+=" ";
142 configuration+=argument;
145 if (!configuration.IsNull()) {
146 iResult=Configure(configuration.Data());
148 iResult=Reconfigure(NULL, NULL);
151 if(fUseGeneral == kFALSE){
152 fHistTPCSideAmax = new TH2F("fHistTPCSideAmax","TPC side A (max signal)",250,-250,250,250,-250,250);
153 fHistTPCSideAmax->SetXTitle("global X (cm)"); fHistTPCSideAmax->SetYTitle("global Y (cm)");
154 fHistTPCSideCmax = new TH2F("fHistTPCSideCmax","TPC side C (max signal)",250,-250,250,250,-250,250);
155 fHistTPCSideCmax->SetXTitle("global X (cm)"); fHistTPCSideCmax->SetYTitle("global Y (cm)");
157 fHistTPCSideAtot = new TH2F("fHistTPCSideAtot","TPC side A (total signal)",250,-250,250,250,-250,250);
158 fHistTPCSideAtot->SetXTitle("global X (cm)"); fHistTPCSideAtot->SetYTitle("global Y (cm)");
159 fHistTPCSideCtot = new TH2F("fHistTPCSideCtot","TPC side C (total signal)",250,-250,250,250,-250,250);
160 fHistTPCSideCtot->SetXTitle("global X (cm)"); fHistTPCSideCtot->SetYTitle("global Y (cm)");
162 fHistTPCSideArms = new TH2F("fHistTPCSideArms","TPC side A (baseline RMS)",250,-250,250,250,-250,250);
163 fHistTPCSideArms->SetXTitle("global X (cm)"); fHistTPCSideArms->SetYTitle("global Y (cm)");
164 fHistTPCSideCrms = new TH2F("fHistTPCSideCrms","TPC side C (baseline RMS)",250,-250,250,250,-250,250);
165 fHistTPCSideCrms->SetXTitle("global X (cm)"); fHistTPCSideCrms->SetYTitle("global Y (cm)");
170 int AliHLTTPCHistogramHandlerComponent::DoDeinit() {
171 // see header file for class documentation
173 if(fHistTPCSideAmax){ delete fHistTPCSideAmax; fHistTPCSideAmax=NULL; }
174 if(fHistTPCSideCmax){ delete fHistTPCSideCmax; fHistTPCSideCmax=NULL; }
176 if(fHistTPCSideAtot){ delete fHistTPCSideAtot; fHistTPCSideAtot=NULL; }
177 if(fHistTPCSideCtot){ delete fHistTPCSideCtot; fHistTPCSideCtot=NULL; }
179 if(fHistTPCSideArms){ delete fHistTPCSideArms; fHistTPCSideArms=NULL; }
180 if(fHistTPCSideCrms){ delete fHistTPCSideCrms; fHistTPCSideCrms=NULL; }
185 int AliHLTTPCHistogramHandlerComponent::DoEvent(const AliHLTComponentEventData&/* evtData*/, AliHLTComponentTriggerData& /*trigData*/){
186 // see header file for class documentation
188 //HLTInfo("--- Entering DoEvent() in TPCHistogramHandler ---");
190 if(GetFirstInputBlock( kAliHLTDataTypeSOR ) || GetFirstInputBlock( kAliHLTDataTypeEOR )) return 0;
192 if(fUseGeneral == kFALSE){
193 fTotalClusterChargeIROCAll = new TH1F("fTotalClusterChargeIROCAll","Total Charge of clusters in all IROC",4000,0,4000);
194 fTotalClusterChargeOROCAll = new TH1F("fTotalClusterChargeOROCAll","Total Charge of clusters in all OROC",4000,0,4000);
195 fQMaxPartitionAll = new TH1F("fQMaxPartitionAll", "QMax for All Partitions", 216,0,216);
196 fPlotQmaxROCAll = new TH1F("fQMaxROCAll", "QMax for All ROC", 72,0,72);
197 fNumberOfClusters = new TH1F("fNumberOfClusters", "Total Number of Clusters", 1,0,1);
200 const TObject *iter = NULL;
201 for(iter = GetFirstInputObject(kAliHLTDataTypeHistogram|kAliHLTDataOriginTPC); iter != NULL; iter = GetNextInputObject()){
203 // HLTInfo("Event 0x%08LX (%Lu) received datatype: %s - required datatype: %s",
204 // evtData.fEventID, evtData.fEventID,
205 // DataType2Text(GetDataType(iter)).c_str(),
206 // DataType2Text(kAliHLTDataTypeHistogram | kAliHLTDataOriginTPC).c_str());
208 // if (GetDataType(iter) == (kAliHLTDataTypeHistogram | kAliHLTDataOriginTPC) && GetEventCount()<2){
209 // HLTWarning("data type %s is depricated, use %s (kAliHLTDataTypeHistogram)!",
210 // DataType2Text(kAliHLTDataTypeHistogram).c_str(),
211 // DataType2Text(kAliHLTDataTypeHistogram | kAliHLTDataOriginTPC).c_str());
214 if (GetDataType(iter) != (kAliHLTDataTypeHistogram | kAliHLTDataOriginTPC)) continue;
216 // Summing the output histograms of the AliHLTTPCNoiseMapComponent (from partition to TPC sides)
218 if(fUseGeneral == kFALSE){
219 if(fNoiseHistograms == kTRUE){
221 fHistTH2Tmp = (TH2F*)iter;
222 TString histName = fHistTH2Tmp->GetName();
224 UInt_t minSlice = AliHLTTPCDefinitions::GetMinSliceNr(GetSpecification(iter));
225 UInt_t maxSlice = AliHLTTPCDefinitions::GetMaxSliceNr(GetSpecification(iter));
226 UInt_t minPartition = AliHLTTPCDefinitions::GetMinPatchNr(GetSpecification(iter));
227 UInt_t maxPartition = AliHLTTPCDefinitions::GetMaxPatchNr(GetSpecification(iter));
229 if((minSlice!=maxSlice) || (minPartition!=maxPartition)){
230 HLTWarning("TPCHistogramHandler::The Noise Map component is not running on partition level!");
233 // minSlice=maxSlice, when the Noise Map component runs on partition level (as it should)
235 if (histName.Contains("fHistSideAMaxSignal")) fHistTPCSideAmax->Add(fHistTH2Tmp,1);
236 else if(histName.Contains("fHistSideATotSignal")) fHistTPCSideAtot->Add(fHistTH2Tmp,1);
237 else if(histName.Contains("fHistSideAPadRMS")) fHistTPCSideArms->Add(fHistTH2Tmp,1);
238 else if(histName.Contains("fHistSideCMaxSignal")) fHistTPCSideCmax->Add(fHistTH2Tmp,1);
239 else if(histName.Contains("fHistSideCTotSignal")) fHistTPCSideCtot->Add(fHistTH2Tmp,1);
240 else if(histName.Contains("fHistSideCPadRMS")) fHistTPCSideCrms->Add(fHistTH2Tmp,1);
243 } // endif fNoiseHistograms==kTRUE
245 // Summing the output of AliHLTTPCClusterHistoComponent
246 if(fKryptonHistograms){
247 Int_t thisrow=-1,thissector=-1,row=-1;
249 AliHLTUInt8_t slice = AliHLTTPCDefinitions::GetMinSliceNr(GetSpecification(iter));
250 AliHLTUInt8_t patch = AliHLTTPCDefinitions::GetMinPatchNr(GetSpecification(iter));
251 row = AliHLTTPCTransform::GetFirstRow(patch);
252 AliHLTTPCTransform::Slice2Sector(slice,row,thissector,thisrow);
254 fHistTH1Tmp = (TH1F*)iter;
255 TString name = fHistTH1Tmp->GetName();
257 if(name=="fTotalClusterChargeIROCAll"){
258 fTotalClusterChargeIROCAll->Add(fTotalClusterChargeIROCAll,fHistTH1Tmp,1,1);
260 else if(name=="fTotalClusterChargeOROCAll"){
261 fTotalClusterChargeOROCAll->Add(fTotalClusterChargeOROCAll,fHistTH1Tmp,1,1);
263 else if(name=="fQMaxPartitionAll"){
264 for(Int_t t=0;t<216;t++){
265 if(fHistTH1Tmp->GetBinContent(t)>fQMaxPartitionAll->GetBinContent(t)){
266 fQMaxPartitionAll->SetBinContent(t,fHistTH1Tmp->GetBinContent(t));
270 else if(name=="fQMaxROCAll"){
271 for(Int_t t=0;t<72;t++){
272 if(fHistTH1Tmp->GetBinContent(t)>fPlotQmaxROCAll->GetBinContent(t)){
273 fPlotQmaxROCAll->SetBinContent(t,fHistTH1Tmp->GetBinContent(t));
277 else if(name=="fNumberOfClusters"){
278 fNumberOfClusters->Add(fNumberOfClusters,fHistTH1Tmp,1,1);
281 HLTWarning("No histogram names match. %s",name.Data());
284 } //endif fKryptonHistograms==kTRUE
287 else { // means fUseGeneral ==kTRUE
289 TH1 * tmp = (TH1*)iter;
290 TString histName = tmp->GetName();
291 UInt_t minSlice = AliHLTTPCDefinitions::GetMinSliceNr(GetSpecification(iter));
292 UInt_t maxSlice = AliHLTTPCDefinitions::GetMaxSliceNr(GetSpecification(iter));
293 UInt_t minPartition = AliHLTTPCDefinitions::GetMinPatchNr(GetSpecification(iter));
294 UInt_t maxPartition = AliHLTTPCDefinitions::GetMaxPatchNr(GetSpecification(iter));
296 Bool_t histogramNotAdded = kTRUE;
298 for(UInt_t i=0;i<fHistogramData.size();i++){
299 if(fIgnoreSpecification == kTRUE){
302 if(histName.Contains(Form("_Slice_%.2d%.2d_Partition_%.2d%.2d", minSlice, maxSlice, minPartition, maxPartition))){
303 cout << "HistogramContains the given string." << endl;
306 histName.ReplaceAll(Form("_Slice_%.2d%.2d_Partition_%.2d%.2d", minSlice, maxSlice, minPartition, maxPartition),"");
309 if(histName.CompareTo(fHistogramData.at(i).fHistogram->GetName())==0){
310 if((minSlice==fHistogramData.at(i).fMinSlice && maxSlice == fHistogramData.at(i).fMaxSlice) || fIgnoreSpecification == kTRUE){
311 if((minPartition==fHistogramData.at(i).fMinPartition && maxPartition == fHistogramData.at(i).fMaxPartition) || fIgnoreSpecification == kTRUE){
312 fHistogramData.at(i).fHistogram->Add(tmp);
313 histogramNotAdded = kFALSE;
317 HLTWarning("Histogram with same name does not have the same partition specification");
321 HLTWarning("Histogram with same name does not have the same slice specification");
326 if(fHistogramData.size()==0){
327 if(fIgnoreSpecification == kTRUE){
328 if(histName.Contains(Form("_Slice_%.2d%.2d_Partition_%.2d%.2d", minSlice, maxSlice, minPartition, maxPartition))){
329 cout << "HistogramContains the given string." << endl;
331 histName.ReplaceAll(Form("_Slice_%.2d%.2d_Partition_%.2d%.2d", minSlice, maxSlice, minPartition, maxPartition),"");
335 if(histogramNotAdded == kTRUE){
336 tmp->SetName(histName);
337 AliHLTHistogramData histogramData;
338 histogramData.fHistogram = tmp;
339 histogramData.fMinSlice = minSlice;
340 histogramData.fMaxSlice = maxSlice;
341 histogramData.fMinPartition = minPartition;
342 histogramData.fMaxPartition = maxPartition;
343 fHistogramData.push_back(histogramData);
346 } // end for loop over histogram blocks
352 void AliHLTTPCHistogramHandlerComponent::MakeHistosPublic() {
353 // see header file for class documentation
355 if(fUseGeneral == kFALSE){
356 if(fNoiseHistograms){
357 PushBack((TObject*)fHistTPCSideAmax,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification( 0,17,0,5));
358 PushBack((TObject*)fHistTPCSideCmax,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(18,35,0,5));
360 PushBack((TObject*)fHistTPCSideAtot,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification( 0,17,0,5));
361 PushBack((TObject*)fHistTPCSideCtot,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(18,35,0,5));
363 PushBack((TObject*)fHistTPCSideArms,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification( 0,17,0,5));
364 PushBack((TObject*)fHistTPCSideCrms,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(18,35,0,5));
366 //if(fHistTH2Tmp) { delete fHistTH2Tmp; fHistTH2Tmp=NULL; }
367 // if(fHistTPCSideA){ delete fHistTPCSideA; fHistTPCSideA=NULL; }
368 // if(fHistTPCSideC){ delete fHistTPCSideC; fHistTPCSideC=NULL; }
371 if(fKryptonHistograms){
372 PushBack((TObject*)fTotalClusterChargeIROCAll,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(0,17,0,1));
373 PushBack((TObject*)fTotalClusterChargeOROCAll,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(0,17,2,5));
374 PushBack((TObject*)fQMaxPartitionAll, kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(0,35,0,5));
375 PushBack((TObject*)fPlotQmaxROCAll, kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(0,35,0,5));
376 PushBack((TObject*)fNumberOfClusters, kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(0,35,0,5));
378 if(fTotalClusterChargeIROCAll){ delete fTotalClusterChargeIROCAll; fTotalClusterChargeIROCAll=NULL; }
379 if(fTotalClusterChargeOROCAll){ delete fTotalClusterChargeOROCAll; fTotalClusterChargeOROCAll=NULL; }
380 if(fQMaxPartitionAll) { delete fQMaxPartitionAll; fQMaxPartitionAll=NULL; }
381 if(fPlotQmaxROCAll) { delete fPlotQmaxROCAll; fPlotQmaxROCAll=NULL; }
382 if(fNumberOfClusters) { delete fNumberOfClusters; fNumberOfClusters=NULL; }
383 if(fHistTH1Tmp) { delete fHistTH1Tmp; fHistTH1Tmp=NULL; }
386 else{ // means fUseGeneral == kTRUE
387 for(UInt_t i=0;i<fHistogramData.size();i++){
388 PushBack((TObject*)fHistogramData.at(i).fHistogram,kAliHLTDataTypeHistogram,AliHLTTPCDefinitions::EncodeDataSpecification(
389 fHistogramData.at(i).fMinSlice,fHistogramData.at(i).fMaxSlice,fHistogramData.at(i).fMinPartition,fHistogramData.at(i).fMaxPartition));
391 fHistogramData.clear();
395 // if(fPlotSideA) histos.Add(fHistSideA);
396 // if(fPlotSideC) histos.Add(fHistSideC);
397 // if(fApplyNoiseMap) histos.Add(fHistCDBMap);
399 // TIter iterator(&histos);
400 // while(TObject *pObj=iterator.Next()){
402 // PushBack(pObj, kAliHLTDataTypeHistogram|kAliHLTDataOriginTPC, fSpecification);
406 // //PushBack( (TObject*) &histos, kAliHLTDataTypeHistogram, fSpecification);
409 int AliHLTTPCHistogramHandlerComponent::Configure(const char* arguments) {
410 // see header file for class documentation
413 if (!arguments) return iResult;
414 HLTInfo("parsing configuration string \'%s\'", arguments);
416 TString allArgs=arguments;
420 TObjArray* pTokens=allArgs.Tokenize(" ");
422 for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
423 argument=((TObjString*)pTokens->At(i))->GetString();
424 if (argument.IsNull()) continue;
426 if (argument.CompareTo("-sum-noise-histograms")==0) {
427 fNoiseHistograms = kTRUE;
428 HLTInfo("got \'-sum-noise-histograms\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
430 } else if (argument.CompareTo("-sum-krypton-histograms")==0) {
431 fKryptonHistograms = kTRUE;
432 HLTInfo("got \'-sum-krypton-histograms\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
434 } else if (argument.CompareTo("-use-general")==0) {
436 HLTInfo("got \'-use-general\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
438 } else if (argument.CompareTo("-ignore-specification")==0) {
439 fIgnoreSpecification = kTRUE;
440 HLTInfo("got \'-ignore-specification\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
443 HLTError("unknown argument %s", argument.Data());
449 /*for (int i=0; i<pTokens->GetEntries() && iResult>=0; i++) {
450 argument=((TObjString*)pTokens->At(i))->GetString();
451 if (argument.IsNull()) continue;
453 if (argument.CompareTo("-sum-noise-histograms")==0) {
454 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
455 HLTInfo("got \'-sum-noise-histograms\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
457 } else if (argument.CompareTo("-sum-krypton-histograms")==0) {
458 if ((bMissingParam=(++i>=pTokens->GetEntries()))) break;
459 HLTInfo("got \'-sum-krypton-histograms\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
461 } else if (argument.CompareTo("-use-general")==0) {
463 HLTInfo("got \'-use-general\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
465 } else if (argument.CompareTo("-ignore-specification")==0) {
466 fIgnoreSpecification = kTRUE;
467 HLTInfo("got \'-ignore-specification\': %s", ((TObjString*)pTokens->At(i))->GetString().Data());
470 HLTError("unknown argument %s", argument.Data());
481 HLTError("missing parameter for argument %s", argument.Data());
487 int AliHLTTPCHistogramHandlerComponent::Reconfigure(const char* cdbEntry, const char* chainId) {
488 // see header file for class documentation
491 const char* path="HLT/ConfigTPC/TPCHistogramHandlerComponent";
492 const char* defaultNotify="";
495 defaultNotify=" (default)";
499 HLTInfo("reconfigure from entry %s%s, chain id %s", path, defaultNotify,(chainId!=NULL && chainId[0]!=0)?chainId:"<none>");
500 AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path/*,GetRunNo()*/);
502 TObjString* pString=dynamic_cast<TObjString*>(pEntry->GetObject());
504 HLTInfo("received configuration object string: \'%s\'", pString->GetString().Data());
505 iResult=Configure(pString->GetString().Data());
507 HLTError("configuration object \"%s\" has wrong type, required TObjString", path);
510 HLTError("cannot fetch object \"%s\" from CDB", path);