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 #include "AliMUONPedestalEventGenerator.h"
20 #include "AliCodeTimer.h"
22 #include "AliHeader.h"
24 #include "AliMUONCalibrationData.h"
25 #include "AliMUONRawWriter.h"
26 #include "AliMUONVCalibParam.h"
27 #include "AliMUONVDigit.h"
28 #include "AliMUONVDigitStore.h"
29 #include "AliMUONVStore.h"
30 #include "AliMpCathodType.h"
31 #include "AliMpConstants.h"
32 #include "AliMpDEStore.h"
33 #include "AliMpDetElement.h"
34 #include "AliMpPlaneType.h"
35 #include "AliRawDataHeaderSim.h"
36 #include "AliLoader.h"
37 #include "AliRunLoader.h"
38 #include <TClonesArray.h>
42 #include <TStopwatch.h>
48 //-----------------------------------------------------------------------------
49 /// \class AliMUONPedestalEventGenerator
51 /// Generate simulated pedestal events for MUON TRK, to be able to e.g. test
52 /// online calibration routines.
54 /// The pedestals themselves are taken from the CDB. What we get from the CDB
55 /// is, per channel, the mean and the sigma of the pedestal. We then use
56 /// those informations to randomly get the pedestals for each channel, for
57 /// each event (picking in a gaus(mean,sigma)).
59 /// Output can be just digits, or digits + raw (ddl), or digits + raw (ddl)
60 /// + raw (date files, one per LDC), depending of ctor and MakeDDL() method.
62 /// \author L. Aphecetche
63 //-----------------------------------------------------------------------------
66 ClassImp(AliMUONPedestalEventGenerator)
69 Int_t AliMUONPedestalEventGenerator::fgCounter(0);
71 //std::streambuf* RedirectTo(std::ostream& what, std::ostream& to)
73 // std::streambuf* old = what.rdbuf();
75 // std::streambuf* psbuf = to.rdbuf();
81 //_____________________________________________________________________________
82 AliMUONPedestalEventGenerator::AliMUONPedestalEventGenerator(Int_t runNumber,
85 : TTask("AliMUONPedestalEventGenerator","Generate fake pedestal events"),
86 fCalibrationData(new AliMUONCalibrationData(runNumber)),
87 fDateFileName(filename),
88 fGAliceFileName("galice.root"),
91 fPedestals(fCalibrationData->Pedestals()),
95 /// Will generate pedestals according to (mean,sigma)s found in CDB
96 /// for run runNumber.
97 /// Will generate nevents events
98 /// If filename is != "", it will be the basename of the output LDC files
100 if (!gSystem->IsAbsoluteFileName(fGAliceFileName))
102 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
104 fGAliceFileName = absFileName;
105 delete[] absFileName;
108 AliRunLoader* runLoader = LoadRun("recreate");
111 AliError("Could not create RunLoader");
115 runLoader->SetNumberOfEventsPerFile(nevents);
117 // Initialize event headers.
118 runLoader->MakeTree("E");
120 for ( Int_t iEvent = 0; iEvent < nevents; ++iEvent )
122 runLoader->SetEventNumber(iEvent);
123 runLoader->GetHeader()->Reset(runNumber,iEvent);
124 runLoader->TreeE()->Fill();
126 runLoader->WriteHeader("OVERWRITE");
127 runLoader->CdGAFile();
128 runLoader->Write(0, TObject::kOverwrite);
135 //_____________________________________________________________________________
136 AliMUONPedestalEventGenerator::~AliMUONPedestalEventGenerator()
139 delete fCalibrationData;
140 AliInfo(Form("make a digit counter %d",fgCounter));
145 //_____________________________________________________________________________
147 AliMUONPedestalEventGenerator::ConvertRawFilesToDate()
149 /// convert raw data DDL files to DATE files with the program "dateStream".
150 /// we make one file per LDC
152 AliCodeTimerAuto("",0)
154 AliInfo("Converting raw to date");
156 const Int_t kIDet = AliDAQ::DetectorID("MUONTRK");
158 const Int_t kNLDCs = 5;//TMath::CeilNint(AliDAQ::NumberOfLdcs(kIDet));
160 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
163 AliError("the program dateStream was not found");
169 AliRunLoader* runLoader = LoadRun("read");
170 if (!runLoader) return kFALSE;
172 AliInfo(Form("converting raw data DDL files to DATE files %s", fDateFileName.Data()));
173 FILE** pipe = new FILE*[kNLDCs];
175 for ( Int_t iFile = 0; iFile < kNLDCs; ++iFile)
178 // Note the option -s. It is used in order to avoid
179 // the generation of SOR/EOR events.
180 snprintf(command, 256, "dateStream -c -D -o %s.LDC%d -# %d -C",
181 fDateFileName.Data(), iFile, runLoader->GetNumberOfEvents());
182 pipe[iFile] = gSystem->OpenPipe(command, "w");
185 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); ++iEvent)
190 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(kIDet); ++iDDL)
192 Int_t ddlID = AliDAQ::DdlID(kIDet,iDDL);
193 Int_t ldcID = Int_t(ldc + 0.0001);
194 ldc += AliDAQ::NumberOfLdcs(kIDet) / AliDAQ::NumberOfDdls(kIDet);
196 char rawFileName[256];
197 snprintf(rawFileName, 256, "raw%d/%s",
198 iEvent, AliDAQ::DdlFileName(kIDet,iDDL));
200 // check existence and size of raw data file
201 FILE* file = fopen(rawFileName, "rb");
203 fseek(file, 0, SEEK_END);
204 unsigned long size = ftell(file);
208 if (ldcID != prevLDC) {
209 fprintf(pipe[ldcID], " LDC Id %d\n", ldcID);
212 fprintf(pipe[ldcID], " Equipment Id %d Payload %s\n", ddlID, rawFileName);
218 for ( Int_t iFile = 0; iFile < kNLDCs; ++iFile)
220 result += gSystem->ClosePipe(pipe[iFile]);
223 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); ++iEvent)
226 snprintf(command, 256, "rm -r raw%d", iEvent);
227 gSystem->Exec(command);
233 return (result == 0);
236 //_____________________________________________________________________________
238 AliMUONPedestalEventGenerator::DigitStore()
240 /// Return digt container; create it if it does not exist
242 if (!fDigitStore) fDigitStore = AliMUONVDigitStore::Create("AliMUONDigitStoreV2R");
246 //_____________________________________________________________________________
248 AliMUONPedestalEventGenerator::Exec(Option_t*)
250 /// Main steering method
252 AliCodeTimerAuto("",0)
256 AliError("No pedestal store. Cannot proceed.");
260 AliRunLoader* runLoader = LoadRun("update");
262 Int_t nevents = runLoader->GetNumberOfEvents();
264 for ( Int_t i = 0; i < nevents ; ++i )
266 runLoader->GetEvent(i);
268 fLoader->MakeDigitsContainer();
269 TTree* treeD = fLoader->TreeD();
272 AliError(Form("Could not get TreeD for event %d",i));
276 DigitStore()->Connect(*treeD);
278 GenerateDigits(*(DigitStore()));
280 // Fill the output treeD
283 // Write to the output tree(D).
284 // Please note that as GlobalTrigger, LocalTrigger and Digits are in the same
285 // tree (=TreeD) in different branches, this WriteDigits in fact writes all of
288 AliCodeTimerStart("WriteDigits")
289 fLoader->WriteDigits("OVERWRITE");
290 AliCodeTimerStop("WriteDigits")
292 fLoader->UnloadDigits();
296 AliCodeTimerAuto("Digits2Raw",1);
301 runLoader->WriteRunLoader("OVERWRITE");
305 // Finally, if instructed to do so, convert DDL files to DATE file(s)
306 if ( fMakeDDL && fDateFileName.Length() > 0 )
308 AliCodeTimerAuto("ConvertRawFilesToDate",1)
309 Bool_t dateOutput = ConvertRawFilesToDate();
312 AliError("DATE output failed. Exiting.");
318 //_____________________________________________________________________________
320 AliMUONPedestalEventGenerator::Digits2Raw(Int_t event)
322 /// Converts digits (from MUON.Digits.root file) to Raw DDL ascii files.
324 AliCodeTimerAuto("",0)
328 AliRawDataHeaderSim header;
329 fRawWriter = new AliMUONRawWriter;
330 fRawWriter->SetHeader(header);
333 // Generate RAW data from the digits
334 // Be carefull to create&change to the correct directory first...
336 TString baseDir = gSystem->WorkingDirectory();
339 snprintf(dirName, 256, "raw%d", event);
340 gSystem->MakeDirectory(dirName);
341 if (!gSystem->ChangeDirectory(dirName))
343 AliError(Form("couldn't change to directory %s", dirName));
347 fRawWriter->Digits2Raw(DigitStore(),0);
349 gSystem->ChangeDirectory(baseDir);
352 //_____________________________________________________________________________
354 AliMUONPedestalEventGenerator::GenerateDigits(AliMUONVDigitStore& digitStore)
356 /// Generate digits (where ADC is set to pedestal value) for all MUON TRK
359 AliCodeTimerAuto("",0)
365 TIter next(fPedestals->CreateIterator());
366 AliMUONVCalibParam* pedestals;
368 while ( ( pedestals = static_cast<AliMUONVCalibParam*>(next())) )
370 Int_t detElemId = pedestals->ID0();
371 Int_t manuId = pedestals->ID1();
373 AliMpDetElement* de = AliMpDEStore::Instance()->GetDetElement(detElemId);
374 AliMp::PlaneType planeType = AliMp::kBendingPlane;
375 if ( manuId & AliMpConstants::ManuMask(AliMp::kNonBendingPlane) )
377 planeType = AliMp::kNonBendingPlane;
379 AliMp::CathodType cathode = de->GetCathodType(planeType);
383 for ( Int_t manuChannel = 0; manuChannel < pedestals->Size(); ++manuChannel )
385 Float_t mean = pedestals->ValueAsFloat(manuChannel,0);
386 if (mean == AliMUONVCalibParam::InvalidFloatValue())
388 // This is a poor's man way of knowing if that channel really exists.
389 // Better and safer way (but much slower too) would be to check pad existence
390 // using AliMpVSegmentation::PadByLocation(manuId,manuChannel)
393 else if ( mean < 1 || mean > 4095 )
395 AliFatal(Form("Got an invalid mean pedestal value for DE %d Manu %d"
396 " channel %d : mean = %e",detElemId,manuId,manuChannel,
401 Float_t sigma = pedestals->ValueAsFloat(manuChannel,1);
405 AliWarning(Form("Got a negative sigma pedestal value for DE %d Manu %d"
406 " channel %d : sigma = %e, will use Abs()=%e",
407 detElemId,manuId,manuChannel,
412 AliMUONVDigit* d = digitStore.Add(detElemId,manuId,manuChannel,
414 AliMUONVDigitStore::kIgnore);
419 ped = gRandom->Gaus(mean,sigma);
421 Int_t pedADC = TMath::FloorNint(ped);
425 // we do not set the remaining parts of the digit, as in principle
426 // this is all we need : manuId, manuChannel and ADC, as far as
427 // real data is concerned.
433 AliDebug(1,Form("ngenerated=%d nmanus=%d",ngenerated,nmanus));
436 //_____________________________________________________________________________
438 AliMUONPedestalEventGenerator::LoadRun(const char* mode)
440 /// Get access to AliRunLoader object
441 while (AliRunLoader::Instance())
443 AliDebug(1,Form("Deleting AliRunLoader %p",AliRunLoader::Instance()));
444 delete AliRunLoader::Instance();
447 AliRunLoader* runLoader =
448 AliRunLoader::Open(fGAliceFileName,AliConfig::GetDefaultEventFolderName(),
451 AliDebug(1,Form("AliRunLoader(%s)=%p",mode,runLoader));
455 AliError("No run loader found in file galice.root");
462 if (smode.Contains("RECREATE"))
464 AliInfo("Creating folder structure");
465 AliConfig::Instance()
466 ->CreateDetectorFolders(runLoader->GetEventFolder(),
468 fLoader = new AliLoader("MUON",runLoader->GetEventFolder());
469 runLoader->AddLoader(fLoader);
472 fLoader = static_cast<AliLoader*>(runLoader->GetDetectorLoader("MUON"));