ALIROOT-5433 Transition to CDHv3 in HLT
[u/mrichter/AliRoot.git] / HLT / MUON / utils / AliHLTMUONRawDataHistoComponent.cxx
CommitLineData
a63da6d6 1/**************************************************************************
2 * This file is property of and copyright by the ALICE HLT Project *
3 * All rights reserved. *
4 * *
5 * Primary Authors: *
6 * Artur Szostak <artursz@iafrica.com> *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
1d8ae082 17// $Id: $
a63da6d6 18
19///
20/// @file AliHLTMUONRawDataHistoComponent.cxx
21/// @author Artur Szostak <artursz@iafrica.com>
22/// @date 30 April 2008
23/// @brief Implementation of the raw data histogramming component for dHLT.
24///
1d8ae082 25/// The class implements a component for checking basic statistics and errors in
26/// raw data from the muon spectrometer. It is useful for basic monitoring.
a63da6d6 27
28#include "AliHLTMUONRawDataHistoComponent.h"
29#include "AliHLTMUONConstants.h"
30#include "AliHLTMUONUtils.h"
31#include "AliHLTDataTypes.h"
32#include "AliCDBEntry.h"
33#include "AliCDBManager.h"
34#include "AliRawDataHeader.h"
16e6f752 35#include "AliHLTCDHWrapper.h"
a63da6d6 36#include "TTimeStamp.h"
37#include <cstdlib>
38#include <cstring>
39#include <cerrno>
40#include <cmath>
41#include <new>
42
43
44// Helper type for memory allocation.
45typedef const AliHLTMUONMansoTrackStruct* AliHLTMUONMansoTrackStructP;
46
47
48ClassImp(AliHLTMUONRawDataHistoComponent);
49
50
51AliHLTMUONRawDataHistoComponent::AliHLTMUONRawDataHistoComponent() :
52 AliHLTMUONProcessor(),
7989fd8e 53 fTrackerDecoder(),
54 fTriggerDecoder(),
a63da6d6 55 fLastPublishTime(-1),
56 fCurrentEventTime(-1),
a35ef960 57 fPublishDelay(0),
7989fd8e 58 fSuppressEmptyHists(false),
a35ef960 59 fProcessDataEventsOnly(false),
60 fClearAfterPublish(false)
a63da6d6 61{
62 /// Default constructor initialises all histogram object pointers to NULL.
63
64 for (int i = 0; i < 22; i++)
65 {
66 fErrorHist[i] = NULL;
67 }
68 for (int i = 0; i < 20; i++)
69 {
70 fManuHist[i] = NULL;
71 fSignalHist[i] = NULL;
72 }
7989fd8e 73
74 fTrackerDecoder.ExitOnError(false);
75 fTrackerDecoder.TryRecover(false);
76 fTrackerDecoder.SendDataOnParityError(true);
77 fTrackerDecoder.AutoDetectTrailer(true);
78 fTrackerDecoder.CheckForTrailer(true);
79
80 fTriggerDecoder.ExitOnError(false);
81 fTriggerDecoder.TryRecover(false);
82 fTriggerDecoder.AutoDetectScalars(false);
a63da6d6 83}
84
85
86AliHLTMUONRawDataHistoComponent::~AliHLTMUONRawDataHistoComponent()
87{
88 /// Default destructor deletes all histogram objects if they are still allocated.
89
90}
91
92
93const char* AliHLTMUONRawDataHistoComponent::GetComponentID()
94{
a63da6d6 95 /// Inherited from AliHLTComponent. Returns the component ID.
a63da6d6 96
97 return AliHLTMUONConstants::RawDataHistogrammerId();
98}
99
100
101void AliHLTMUONRawDataHistoComponent::GetInputDataTypes(AliHLTComponentDataTypeList& list)
102{
a63da6d6 103 /// Inherited from AliHLTProcessor. Returns the list of expected input data types.
a63da6d6 104
105 assert( list.empty() );
106 list.push_back( AliHLTMUONConstants::DDLRawDataType() );
107}
108
109
110AliHLTComponentDataType AliHLTMUONRawDataHistoComponent::GetOutputDataType()
111{
112 /// Inherited from AliHLTComponent. Returns kAliHLTHistogramDataTypeID.
113
114 return AliHLTMUONConstants::HistogramDataType();
115}
116
117
118void AliHLTMUONRawDataHistoComponent::GetOutputDataSize(
119 unsigned long& constBase, double& inputMultiplier
120 )
121{
a63da6d6 122 /// Inherited from AliHLTComponent. Returns an estimate of the expected output data size.
a63da6d6 123
95c6f20f 124 constBase = (sizeof(TH1D)+50*sizeof(double))*22 + (sizeof(TH1D)+1024*4*sizeof(double))*20*2;
a63da6d6 125 inputMultiplier = 0;
126}
127
128
129AliHLTComponent* AliHLTMUONRawDataHistoComponent::Spawn()
130{
a63da6d6 131 /// Inherited from AliHLTComponent. Creates a new object instance.
a63da6d6 132
133 return new AliHLTMUONRawDataHistoComponent;
134}
135
136
137bool AliHLTMUONRawDataHistoComponent::IgnoreArgument(const char* arg) const
138{
139 /// Return true if the argument is one of -cdbpath -run or -delaysetup
140 /// to prevent the parent class from parsing these arguments in DoInit.
141
142 if (strcmp(arg, "-cdbpath") == 0 or strcmp(arg, "-run") == 0 or
143 strcmp(arg, "-delaysetup") == 0)
144 {
145 return true;
146 }
147 else
148 {
149 return false;
150 }
151}
152
153
154int AliHLTMUONRawDataHistoComponent::DoInit(int argc, const char** argv)
155{
156 /// Inherited from AliHLTComponent.
157 /// Parses the command line parameters and initialises the component.
158
159 HLTInfo("Initialising dHLT raw data histogrammer component.");
160
161 // Inherit the parents functionality.
162 int result = AliHLTMUONProcessor::DoInit(argc, argv);
163 if (result != 0) return result;
164
165 fLastPublishTime = fCurrentEventTime = -1;
a35ef960 166 fPublishDelay = 0;
a63da6d6 167 bool pubDelaySet = false;
7989fd8e 168 fSuppressEmptyHists = false;
169 fProcessDataEventsOnly = false;
a35ef960 170 fClearAfterPublish = false;
7989fd8e 171 fTrackerDecoder.TryRecover(false);
172 fTriggerDecoder.TryRecover(false);
a63da6d6 173
174 for (int i = 0; i < argc; i++)
175 {
176 if (ArgumentAlreadyHandled(i, argv[i])) continue;
177
178 if (strcmp(argv[i], "-pubdelay") == 0)
179 {
180 if (pubDelaySet)
181 {
182 HLTWarning("The publishing delay value was already specified."
183 " Will replace previous value given by -pubdelay."
184 );
185 }
186
187 if (argc <= i+1)
188 {
189 HLTError("The value for the publishing delay was not specified.");
190 return -EINVAL;
191 }
192
193 char* cpErr = NULL;
194 double num = strtod(argv[i+1], &cpErr);
195 if (cpErr == NULL or *cpErr != '\0' or num < 0)
196 {
197 HLTError("Cannot convert '%s' to a positive floating point number.",
198 argv[i+1]
199 );
200 return -EINVAL;
201 }
202 fPublishDelay = num;
203 pubDelaySet = true;
204
205 i++;
206 continue;
207 }
7989fd8e 208
209 if (strcmp(argv[i], "-noemptyhists") == 0)
210 {
211 fSuppressEmptyHists = true;
212 continue;
213 }
214
215 if (strcmp(argv[i], "-onlydataevents") == 0)
216 {
217 fProcessDataEventsOnly = true;
218 continue;
219 }
220
a35ef960 221 if (strcmp(argv[i], "-clearafterpub") == 0)
222 {
223 fClearAfterPublish = true;
224 continue;
225 }
226
7989fd8e 227 if (strcmp(argv[i], "-tryrecover") == 0)
228 {
229 fTrackerDecoder.TryRecover(true);
230 fTriggerDecoder.TryRecover(true);
231 continue;
232 }
a63da6d6 233
234 HLTError("Unknown option '%s'.", argv[i]);
235 return -EINVAL;
236 }
237
238 try
239 {
240 char name[256];
241 char title[1024];
242
243 // Do not add to current directory to prevent memory leak warning.
244 // We will not be leaking any memory if we dont add to the directory.
245 TH1::AddDirectory(kFALSE);
246
247 for (int i = 0; i < 22; i++)
248 {
249 AliHLTInt32_t equipId = AliHLTMUONUtils::DDLNumberToEquipId(i);
250 sprintf(name, "rawDataErrors_%d", equipId);
251 sprintf(title, "Distribution of errors found in raw data from DDL %d.", equipId);
a63da6d6 252 fErrorHist[i] = new TH1D(name, title, 40, 0.5, 40.5);
253 fErrorHist[i]->SetXTitle("Error code");
254 fErrorHist[i]->SetYTitle("Number of errors");
255 }
256 for (int i = 0; i < 20; i++)
257 {
258 AliHLTInt32_t equipId = AliHLTMUONUtils::DDLNumberToEquipId(i);
259 sprintf(name, "manuDistrib_%d", equipId);
260 sprintf(title, "Distribution of MANUs containing raw data in DDL %d.", equipId);
261 fManuHist[i] = new TH1D(name, title, 2048, -0.5, 2047.5);
262 fManuHist[i]->SetXTitle("MANU number (as seen in raw data)");
263 fManuHist[i]->SetYTitle("Number of signals read.");
264 sprintf(name, "signalDistrib_%d", equipId);
265 sprintf(title, "Distribution of signals in raw data from DDL %d.", equipId);
266 fSignalHist[i] = new TH1D(name, title, 4096, -0.5, 4095.5);
267 fSignalHist[i]->SetXTitle("Channels");
268 fSignalHist[i]->SetYTitle("dN/dChannel");
269 }
270 }
271 catch (const std::bad_alloc&)
272 {
273 HLTError("Could not allocate more memory for histogram objects.");
274 FreeObjects();
275 return -ENOMEM;
276 }
277
278 return 0;
279}
280
281
282int AliHLTMUONRawDataHistoComponent::DoDeinit()
283{
284 /// Inherited from AliHLTComponent. Performs a cleanup of the component.
285 /// Will delete all histogram objects.
286
287 HLTInfo("Deinitialising dHLT raw data histogrammer component.");
288 fCurrentEventTime = -1;
289 FreeObjects();
290 return 0;
291}
292
293
294int AliHLTMUONRawDataHistoComponent::DoEvent(
1d8ae082 295 const AliHLTComponentEventData& evtData,
296 AliHLTComponentTriggerData& trigData
a63da6d6 297 )
298{
299 /// Inherited from AliHLTProcessor.
300 /// Processes the new event data and generates summary histograms.
301
7989fd8e 302 if (fProcessDataEventsOnly and not IsDataEvent()) return 0; // Only process data events.
303
a63da6d6 304 fCurrentEventTime = TTimeStamp().AsDouble();
305
306 const AliHLTComponentBlockData* block = GetFirstInputBlock(AliHLTMUONConstants::DDLRawDataType());
307 for ( ; block != NULL; block = GetNextInputBlock())
308 {
309 HLTDebug("Handling block with fDataType = '%s', fPtr = %p,"
310 " fSize = %u bytes and fSpecification = 0x%8.8X.",
311 DataType2Text(block->fDataType).c_str(), block->fPtr,
312 block->fSize, block->fSpecification
313 );
314
315 if (AliHLTMUONUtils::IsTrackerDDL(block->fSpecification))
316 {
1d8ae082 317 bool decodeOk = ProcessTrackerDDL(block);
318 if (not decodeOk and DumpDataOnError())
319 {
320 DumpEvent(evtData, trigData);
321 return -EIO;
322 }
a63da6d6 323 }
324 else if (AliHLTMUONUtils::IsTriggerDDL(block->fSpecification))
325 {
1d8ae082 326 bool decodeOk = ProcessTriggerDDL(block);
327 if (not decodeOk and DumpDataOnError())
328 {
329 DumpEvent(evtData, trigData);
330 return -EIO;
331 }
a63da6d6 332 }
333 else
334 {
335 HLTError("Received a raw data block with an invalid specification of"
336 " 0x%8.8X. Expected raw data only from one DDL and not multiple"
337 " DDLs as indicated by the specification.",
338 block->fSpecification
339 );
340 }
341 }
342
7989fd8e 343 // See if 'fPublishDelay' number of seconds has elapsed or this is the first event,
344 // in that case publish the histograms. Do not publish histograms that are empty
345 // if the fSuppressEmptyHists flag is set.
346 if (fLastPublishTime == -1 or fCurrentEventTime - fLastPublishTime >= fPublishDelay)
a63da6d6 347 {
348 for (int i = 0; i < 22; i++)
349 {
7989fd8e 350 if (fSuppressEmptyHists and fErrorHist[i]->GetEntries() == 0) continue;
a63da6d6 351 PushBack(fErrorHist[i],
352 AliHLTMUONConstants::HistogramDataType(),
353 AliHLTMUONUtils::DDLNumberToSpec(i)
354 );
a35ef960 355 // If requested, clear histogram when published.
356 if (fClearAfterPublish) fErrorHist[i]->Reset("M");
a63da6d6 357 }
358 for (int i = 0; i < 20; i++)
359 {
360 AliHLTUInt32_t spec = AliHLTMUONUtils::DDLNumberToSpec(i);
7989fd8e 361 if (not (fSuppressEmptyHists and fManuHist[i]->GetEntries() == 0))
362 {
363 PushBack(fManuHist[i], AliHLTMUONConstants::HistogramDataType(), spec);
a35ef960 364 if (fClearAfterPublish) fManuHist[i]->Reset("M");
7989fd8e 365 }
366 if (not (fSuppressEmptyHists and fSignalHist[i]->GetEntries() == 0))
367 {
368 PushBack(fSignalHist[i], AliHLTMUONConstants::HistogramDataType(), spec);
a35ef960 369 if (fClearAfterPublish) fSignalHist[i]->Reset("M");
7989fd8e 370 }
a63da6d6 371 }
372 fLastPublishTime = fCurrentEventTime;
373 }
374
375 return 0;
376}
377
378
1d8ae082 379bool AliHLTMUONRawDataHistoComponent::ProcessTrackerDDL(const AliHLTComponentBlockData* block)
a63da6d6 380{
381 /// Processes a raw data block from the tracker stations.
382
383 AliHLTInt32_t ddl = AliHLTMUONUtils::SpecToDDLNumber(block->fSpecification);
384 assert(0 <= ddl and ddl < 20);
385
7989fd8e 386 fTrackerDecoder.GetHandler().ErrorHist(fErrorHist[ddl]);
387 fTrackerDecoder.GetHandler().ManuHist(fManuHist[ddl]);
388 fTrackerDecoder.GetHandler().SignalHist(fSignalHist[ddl]);
389
16e6f752 390 AliHLTCDHWrapper cdh(block->fPtr);
391 if (block->fSize >= sizeof(AliRawDataHeader) &&
392 block->fSize >= cdh.GetHeaderSize()) // in case if cdh v3
a63da6d6 393 {
7989fd8e 394 AliHLTUInt8_t* payload = reinterpret_cast<AliHLTUInt8_t*>(block->fPtr)
16e6f752 395 + cdh.GetHeaderSize();
396 UInt_t payloadSize = UInt_t(block->fSize) - cdh.GetHeaderSize();
1d8ae082 397 return fTrackerDecoder.Decode(payload, payloadSize);
a63da6d6 398 }
399 else
400 {
401 HLTError("Received a raw data block that is too short to be valid."
402 " Its size is only %d bytes",
403 block->fSize
404 );
405 fErrorHist[ddl]->Fill(40);
406 }
1d8ae082 407 return false;
a63da6d6 408}
409
410
1d8ae082 411bool AliHLTMUONRawDataHistoComponent::ProcessTriggerDDL(const AliHLTComponentBlockData* block)
a63da6d6 412{
413 /// Processes a raw data block from the trigger stations.
414
a63da6d6 415 AliHLTInt32_t ddl = AliHLTMUONUtils::SpecToDDLNumber(block->fSpecification);
416 assert(21 <= ddl and ddl < 22);
417
7989fd8e 418 fTriggerDecoder.GetHandler().ErrorHist(fErrorHist[ddl]);
419
16e6f752 420 AliHLTCDHWrapper cdh(block->fPtr);
421 if (block->fSize >= sizeof(AliRawDataHeader) &&
422 block->fSize >= cdh.GetHeaderSize()) // in case if cdh v3)
a63da6d6 423 {
16e6f752 424 AliHLTUInt8_t* payload = reinterpret_cast<AliHLTUInt8_t*>(block->fPtr);
425 payload += cdh.GetHeaderSize();
426 UInt_t payloadSize = UInt_t(block->fSize) - cdh.GetHeaderSize();
427 bool scalarEvent = ((cdh.GetL1TriggerMessage() & 0x1) == 0x1);
1d8ae082 428 return fTriggerDecoder.Decode(payload, payloadSize, scalarEvent);
a63da6d6 429 }
430 else
431 {
432 HLTError("Received a raw data block that is too short to be valid."
433 " Its size is only %d bytes",
434 block->fSize
435 );
436 fErrorHist[ddl]->Fill(40);
437 }
1d8ae082 438 return false;
a63da6d6 439}
440
441
442void AliHLTMUONRawDataHistoComponent::FreeObjects()
443{
444 /// Deletes all the histogram objects that were allocated.
445
446 for (int i = 0; i < 22; i++)
447 {
448 if (fErrorHist[i] != NULL)
449 {
450 delete fErrorHist[i];
451 fErrorHist[i] = NULL;
452 }
453 }
454 for (int i = 0; i < 20; i++)
455 {
456 if (fManuHist[i] != NULL)
457 {
458 delete fManuHist[i];
459 fManuHist[i] = NULL;
460 }
461 if (fSignalHist[i] != NULL)
462 {
463 delete fSignalHist[i];
464 fSignalHist[i] = NULL;
465 }
466 }
467}
7989fd8e 468