activating individual HLT simulations from digits and raw data
[u/mrichter/AliRoot.git] / HLT / sim / AliHLTSimulation.cxx
CommitLineData
4cbaf07b 1// $Id$
2
c5123824 3//**************************************************************************
4//* This file is property of and copyright by the ALICE HLT Project *
5//* ALICE Experiment at CERN, All rights reserved. *
6//* *
7//* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> *
8//* for The ALICE HLT Project. *
9//* *
10//* Permission to use, copy, modify and distribute this software and its *
11//* documentation strictly for non-commercial purposes is hereby granted *
12//* without fee, provided that the above copyright notice appears in all *
13//* copies and that both the copyright notice and this permission notice *
14//* appear in the supporting documentation. The authors make no claims *
15//* about the suitability of this software for any purpose. It is *
16//* provided "as is" without express or implied warranty. *
17//**************************************************************************
4cbaf07b 18
19/** @file AliHLTSimulation.cxx
20 @author Matthias Richter
21 @date
22 @brief Binding class for HLT simulation in AliRoot. */
23
24#include <cassert>
25#include <cerrno>
7777fa45 26#include "TObjArray.h"
27#include "TObjString.h"
4cbaf07b 28#include "AliHLTSimulation.h"
5a29e5cb 29#include "AliSimulation.h"
4cbaf07b 30#include "AliLog.h"
90c37647 31#include "AliRun.h"
4cbaf07b 32#include "AliRunLoader.h"
90c37647 33#include "AliHeader.h"
90c37647 34#include "AliCDBManager.h"
35#include "AliCDBEntry.h"
36#include "AliCDBPath.h"
37#include "AliCDBId.h"
38#include "AliCDBMetaData.h"
8dd79235 39#include "AliCDBStorage.h"
40#include "AliGRPObject.h"
9f199018 41#include "AliGRPManager.h"
4cbaf07b 42#include "AliHLTSystem.h"
5879742f 43#include "AliHLTConfigurationHandler.h"
7bf6c76d 44#include "AliHLTPluginBase.h"
7777fa45 45#include "AliRawReaderFile.h"
46#include "AliRawReaderDate.h"
47#include "AliRawReaderRoot.h"
a8420176 48#include "AliESDEvent.h"
2c0e5942 49#include "AliHLTOUTComponent.h"
9f199018 50#include "AliTracker.h"
8dd79235 51#include "TGeoGlobalMagField.h"
52#include "TSystem.h"
4e1f4995 53#include "TMath.h"
9f199018 54#include "TGeoGlobalMagField.h"
4cbaf07b 55
56#if ALIHLTSIMULATION_LIBRARY_VERSION != LIBHLTSIM_VERSION
57#error library version in header file and lib*.pkg do not match
58#endif
59
60/** ROOT macro for the implementation of ROOT specific class methods */
61ClassImp(AliHLTSimulation);
62
63AliHLTSimulation::AliHLTSimulation()
5879742f 64 : fOptions()
65 , fpPluginBase(new AliHLTPluginBase)
66 , fpRawReader(NULL)
4cbaf07b 67{
68 // see header file for class documentation
69 // or
70 // refer to README to build package
71 // or
72 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
73}
74
75AliHLTSimulation::~AliHLTSimulation()
76{
77 // see header file for function documentation
7bf6c76d 78 if (fpPluginBase) delete fpPluginBase;
79 fpPluginBase=NULL;
80
7777fa45 81 if (fpRawReader) {
82 delete fpRawReader;
83 }
84 fpRawReader=NULL;
4cbaf07b 85}
86
87AliHLTSimulation* AliHLTSimulation::CreateInstance()
88{
89 // see header file for function documentation
90 return new AliHLTSimulation;
91}
92
93int AliHLTSimulation::DeleteInstance(AliHLTSimulation* pSim)
94{
95 // see header file for function documentation
96 assert(pSim!=NULL);
97 delete pSim;
98 return 0;
99}
100
101int AliHLTSimulation::Init(AliRunLoader* pRunLoader, const char* options)
102{
103 // init the simulation
104 fOptions=options;
7777fa45 105 TString sysOp;
4cbaf07b 106
7bf6c76d 107 if(!fpPluginBase) {
108 AliError("internal initialization failed");
109 return -EINVAL;
110 }
111
112 AliHLTSystem* pSystem=fpPluginBase->GetInstance();
113 if (!pSystem) {
114 AliError("can not get AliHLTSystem instance");
4cbaf07b 115 return -ENOMEM;
116 }
7bf6c76d 117 if (pSystem->CheckStatus(AliHLTSystem::kError)) {
4cbaf07b 118 AliError("HLT system in error state");
119 return -EFAULT;
120 }
121
7777fa45 122 // scan options for specific entries
123 TObjArray* pTokens=fOptions.Tokenize(" ");
124 if (pTokens) {
125 int iEntries=pTokens->GetEntries();
126 for (int i=0; i<iEntries; i++) {
127 TString token=(((TObjString*)pTokens->At(i))->GetString());
128 if (token.Contains("rawfile=")) {
129 TString param=token.ReplaceAll("rawfile=", "");
130 if (param.EndsWith("/")) {
131 AliInfo(Form("creating AliRawReaderFile (%s)", param.Data()));
132 fpRawReader = new AliRawReaderFile(param);
133 } else if (param.EndsWith(".root")) {
134 AliInfo(Form("creating AliRawReaderRoot (%s)", param.Data()));
135 fpRawReader = new AliRawReaderRoot(param);
136 } else if (!param.IsNull()) {
137 AliInfo(Form("creating AliRawReaderDate (%s)", param.Data()));
138 fpRawReader = new AliRawReaderDate(param);
7777fa45 139 }
140 if (fpRawReader) {
141 fpRawReader->RewindEvents();
142 int count=0;
32f88fc5 143 for ( ; fpRawReader->NextEvent(); count++) {/* empty body */};
7777fa45 144 if (count!=pRunLoader->GetNumberOfEvents()) {
c5123824 145 AliError(Form("mismatch in event count: runloader %d, rawreader %d; ignoring rawreader",
7777fa45 146 pRunLoader->GetNumberOfEvents(), count));
147 count=0;
148 }
149 if (count>0) {
150 fpRawReader->RewindEvents();
151 fpRawReader->NextEvent();
152 } else {
153 delete fpRawReader;
154 fpRawReader=NULL;
155 }
156 }
2c0e5942 157 } else if (token.Contains("writerawfiles=")) {
158 if (!token.ReplaceAll("writerawfiles=", "").Contains("HLT")) {
5879742f 159 if (TestBit(kOneChain) && AliHLTOUTComponent::TestGlobalOption(AliHLTOUTComponent::kWriteRawFiles)) {
160 AliWarning("empty argument 'writerawfiles=' disables HLTOUTComponent mode 'raw' which was set by argument 'hltout-mode'");
161 }
2c0e5942 162 AliHLTOUTComponent::ClearGlobalOption(AliHLTOUTComponent::kWriteRawFiles);
163 }
5879742f 164 } else if (token.BeginsWith("hltout-mode=")) {
165 // this is a legacy mode to emulate the behavior before Dec 2010 where only
166 // one chain was executed on either digits or simulated raw data and the output
167 // was controlled via global flags
168 // add to the arguments for AliHLTSystem as also there the information is needed
169 if (sysOp.Length()>0) sysOp+=" ";
170 sysOp+=token;
171 TString param=token.ReplaceAll("hltout-mode=", "");
172 SetBit(kOneChain);
173 if (param.CompareTo("raw")==0) {
174 // please note that this option
175 AliHLTOUTComponent::SetGlobalOption(AliHLTOUTComponent::kWriteRawFiles);
176 AliHLTOUTComponent::ClearGlobalOption(AliHLTOUTComponent::kWriteDigits);
177 } else if (param.CompareTo("digits")==0) {
178 // please note that this option
179 AliHLTOUTComponent::ClearGlobalOption(AliHLTOUTComponent::kWriteRawFiles);
180 AliHLTOUTComponent::SetGlobalOption(AliHLTOUTComponent::kWriteDigits);
181 } else if (param.CompareTo("legacy")==0) {
182 AliHLTOUTComponent::SetGlobalOption(AliHLTOUTComponent::kWriteRawFiles);
183 AliHLTOUTComponent::SetGlobalOption(AliHLTOUTComponent::kWriteDigits);
184 } else {
185 AliError(Form("invalid parameter for argument 'hltout-mode=' %s, allowed: raw, digits, legacy ... ignoring argument and using the standard simulation", param.Data()));
186 ResetBit(kOneChain);
187 }
7777fa45 188 } else {
189 if (sysOp.Length()>0) sysOp+=" ";
190 sysOp+=token;
191 }
192 }
193 delete pTokens;
194 }
5879742f 195 // only store the options for AliHLTSystem
196 fOptions=sysOp;
197
198 // if no specific hltout-mode has been chosen set the split mode for
199 // running separate chains for digits and raw data
200 if (!fOptions.Contains("hltout-mode=")) fOptions+=" hltout-mode=split";
7777fa45 201
90c37647 202 AliCDBManager* man = AliCDBManager::Instance();
203 if (man && man->IsDefaultStorageSet())
204 {
8dd79235 205 // init solenoid field
9f199018 206 // 2009-11-07 magnetic field handling fo HLT components has been switched to the
207 // global AliMagF instance, the HLT/ConfigHLT/SolenoidBz entry is obsolete
208 // The global instance is either established by the AliRoot environment or the
209 // component external interface.
210 if (TGeoGlobalMagField::Instance()->GetField()) {
211 AliDebug(0, Form("magnetic field: %f", AliTracker::GetBz()));
8dd79235 212 } else {
213 // workaround for bug #51285
9f199018 214 AliGRPManager grpman;
215 if (grpman.ReadGRPEntry() &&
216 grpman.SetMagField()) {
217 // nothing to do any more
8dd79235 218 }
9f199018 219 AliError(Form("can not get the AliMagF instance, falling back to GRP entry (%f)", AliTracker::GetBz()));
90c37647 220 }
8dd79235 221 } else if (man) {
222 AliError("OCDB default storage not yet set, can not prepare OCDB entries");
223 } else {
224 AliError("unable to get instance of AliCDBMetaData, can not prepare OCDB entries");
90c37647 225 }
226
5879742f 227 // configure the main HLTSystem instance for digit simulation (pRawReader NULL)
228 return ConfigureHLTSystem(pSystem, fOptions.Data(), pRunLoader, TestBit(kOneChain)?fpRawReader:NULL);
229}
230
231int AliHLTSimulation::ConfigureHLTSystem(AliHLTSystem* pSystem, const char* options, AliRunLoader* pRunLoader, AliRawReader* pRawReader) const
232{
233 // scan options and configure AliHLTSystem
234 if (pSystem->ScanOptions(options)<0) {
4cbaf07b 235 AliError("error setting options for HLT system");
236 return -EINVAL;
237 }
238
7bf6c76d 239 if (!pSystem->CheckStatus(AliHLTSystem::kReady)) {
5879742f 240 if ((pSystem->Configure(pRawReader, pRunLoader))<0) {
4cbaf07b 241 AliError("error during HLT system configuration");
242 return -EFAULT;
243 }
244 }
245
246 return 0;
247}
248
4cbaf07b 249int AliHLTSimulation::Run(AliRunLoader* pRunLoader)
250{
251 // HLT reconstruction for simulated data
7bf6c76d 252 if(!fpPluginBase) {
253 AliError("internal initialization failed");
254 return -EINVAL;
255 }
256
4cbaf07b 257 if(!pRunLoader) {
258 AliError("Missing RunLoader! 0x0");
259 return -EINVAL;
260 }
261
4cbaf07b 262 int iResult=0;
263
7bf6c76d 264 AliHLTSystem* pSystem=fpPluginBase->GetInstance();
265 if (!pSystem) {
266 AliError("can not get AliHLTSystem instance");
267 return -ENOMEM;
268 }
269
270 if (pSystem->CheckStatus(AliHLTSystem::kError)) {
4cbaf07b 271 AliError("HLT system in error state");
272 return -EFAULT;
273 }
7777fa45 274
5879742f 275 // run the main HLTSystem instance for digit simulation (pRawReader NULL)
276 // in legacy mode only one chain is run and the output is controlled via
277 // global flags
278 if (!TestBit(kOneChain)) AliInfo("running HLT simulation for digits");
279 iResult=RunHLTSystem(pSystem, pRunLoader, TestBit(kOneChain)?fpRawReader:NULL);
280
281 // now run once again with the raw data as input, a completely new HLT system
282 // with new configurations is used
283 if (fpRawReader && !TestBit(kOneChain)) {
284 AliInfo("running HLT simulation for raw data");
285 int iLocalResult=0;
286 AliHLTConfigurationHandler* confHandler=new AliHLTConfigurationHandler;
287 // note that the configuration handler is owned by the
288 // AliHLTSystem instance from now on
289 AliHLTSystem rawSimulation(kHLTLogDefault, "", NULL, confHandler);
290 if ((iLocalResult=ConfigureHLTSystem(&rawSimulation, fOptions.Data(), pRunLoader, fpRawReader))>=0) {
291 iLocalResult=RunHLTSystem(&rawSimulation, pRunLoader, fpRawReader);
292 }
293 if (iResult>=0) iResult=iLocalResult;
294 }
295
296 return iResult;
297}
298
299int AliHLTSimulation::RunHLTSystem(AliHLTSystem* pSystem, AliRunLoader* pRunLoader, AliRawReader* pRawReader) const
300{
301 // run reconstruction cycle for AliHLTSystem
302 int nEvents = pRunLoader->GetNumberOfEvents();
303 int iResult=0;
304
7777fa45 305 // Note: the rawreader is already placed at the first event
5879742f 306 if ((iResult=pSystem->Reconstruct(1, pRunLoader, pRawReader))>=0) {
7bf6c76d 307 pSystem->FillESD(0, pRunLoader, NULL);
7777fa45 308 for (int i=1; i<nEvents; i++) {
5879742f 309 if (pRawReader && !pRawReader->NextEvent()) {
c5123824 310 AliError("mismatch in event count, rawreader corrupted");
7777fa45 311 break;
312 }
5879742f 313 pSystem->Reconstruct(1, pRunLoader, pRawReader);
7bf6c76d 314 pSystem->FillESD(i, pRunLoader, NULL);
7777fa45 315 }
dee38f1b 316 // send specific 'event' to execute the stop sequence
7bf6c76d 317 pSystem->Reconstruct(0, NULL, NULL);
4cbaf07b 318 }
319 return iResult;
320}
321
4cbaf07b 322AliHLTSimulation* AliHLTSimulationCreateInstance()
323{
324 // see header file for function documentation
325 return AliHLTSimulation::CreateInstance();
326}
327
328int AliHLTSimulationDeleteInstance(AliHLTSimulation* pSim)
329{
330 // see header file for function documentation
331 return AliHLTSimulation::DeleteInstance(pSim);
332}
333
334int AliHLTSimulationInit(AliHLTSimulation* pSim, AliRunLoader* pRunLoader, const char* options)
335{
336 assert(pSim!=NULL);
337 if (pSim) {
338 return pSim->Init(pRunLoader, options);
339 }
340 return -ENODEV;
341}
342
343int AliHLTSimulationRun(AliHLTSimulation* pSim, AliRunLoader* pRunLoader)
344{
345 assert(pSim!=NULL);
346 if (pSim) {
347 return pSim->Run(pRunLoader);
348 }
349 return -ENODEV;
350}
351
352int AliHLTSimulationGetLibraryVersion()
353{
354 // see header file for function documentation
355 return LIBHLTSIM_VERSION;
356}
b50357c9 357
5a29e5cb 358int AliHLTSimulationSetup(AliHLTSimulation* /*pHLTSim*/, AliSimulation* pSim, const char* specificObjects)
8dd79235 359{
360 // see header file for function documentation
361
362 // this is an attempt to solve issue #48360
363 // since there are many jobs running in parallel during the production,
364 // all the jobs want to put entries into the OCDB. The solution is to
365 // make them temporary, since they are only used to propagate information
366 // from the simulation to the reconstruction.
5a29e5cb 367
368 if (!pSim) return -EINVAL;
369 const char* entries[]={
5a29e5cb 370 NULL
371 };
372
373 TString specificStorage;
374 specificStorage.Form("local://%s",gSystem->pwd());
375 for (const char** pEntry=entries; *pEntry!=NULL; pEntry++) {
376 const char* pObject=specificObjects?strstr(specificObjects, *pEntry):NULL;
377 if (pObject) {
378 // skip this entry if it is found in the list and either
379 // last one or separated by a blank
380 pObject+=strlen(*pEntry);
381 if (*pObject==0 || *pObject==' ') continue;
382 }
383 pSim->SetSpecificStorage(*pEntry, specificStorage.Data());
8dd79235 384 }
1daa8d18 385
386 return 0;
8dd79235 387}
388
5a29e5cb 389#ifndef HAVE_COMPILEINFO
390extern "C" void CompileInfo(const char*& date, const char*& time)
b50357c9 391{
392 // the fall back compile info of the HLTsim library
393 // this is not up-to-date if other files have been changed and recompiled
394 date=__DATE__; time=__TIME__;
395 return;
396}
5a29e5cb 397#endif