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 //____________________________________________________________________
20 // This is a class that constructs AliFMDMult (reconstructed
21 // multiplicity) from of Digits
23 // This class reads either digits from a TClonesArray or raw data from
24 // a DDL file (or similar), and stores the read ADC counts in an
25 // internal cache (fAdcs).
27 // From the cached values it then calculates the number of particles
28 // that hit a region of the FMDs, as specified by the user.
30 // The reconstruction can be done in two ways: Either via counting the
31 // number of empty strips (Poisson method), or by converting the ADC
32 // signal to an energy deposition, and then dividing by the typical
33 // energy loss of a particle.
35 // +---------------------+ +---------------------+
36 // | AliFMDReconstructor |<>-----| AliFMDMultAlgorithm |
37 // +---------------------+ +---------------------+
40 // +-----------+---------+
42 // +-------------------+ +------------------+
43 // | AliFMDMultPoisson | | AliFMDMultNaiive |
44 // +-------------------+ +------------------+
46 // AliFMDReconstructor acts as a manager class. It contains a list of
47 // AliFMDMultAlgorithm objects. The call graph looks something like
50 // +----------------------+ +----------------------+
51 // | :AliFMDReconstructor | | :AliFMDMultAlgorithm |
52 // +----------------------+ +----------------------+
55 // ------------>| | PreRun +-+
56 // | |------------------------------->| |
58 // | |-----+ (for each event) |
59 // | | | *ProcessEvent |
61 // || |<---+ PreEvent +-+
62 // || |------------------------------>| |
65 // || | | ProcessDigits |
68 // ||| | *ProcessDigit(digit) +-+
69 // ||| |----------------------------->| |
73 // || |------------------------------>| |
77 // | |------------------------------->| |
84 //-- Authors: Evgeny Karpechev(INR) and Alla Maevsksia
85 // Latest changes by Christian Holm Christensen <cholm@nbi.dk>
88 //____________________________________________________________________
90 #include <AliLog.h> // ALILOG_H
91 #include <AliRun.h> // ALIRUN_H
92 #include <AliRunLoader.h> // ALIRUNLOADER_H
93 #include <AliLoader.h> // ALILOADER_H
94 #include <AliHeader.h> // ALIHEADER_H
95 #include <AliRawReader.h> // ALIRAWREADER_H
96 #include <AliGenEventHeader.h> // ALIGENEVENTHEADER_H
97 #include "AliFMD.h" // ALIFMD_H
98 #include "AliFMDDigit.h" // ALIFMDDIGIT_H
99 #include "AliFMDReconstructor.h" // ALIFMDRECONSTRUCTOR_H
100 #include "AliFMDRawStream.h" // ALIFMDRAWSTREAM_H
101 #include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
102 #include "AliFMDMultAlgorithm.h" // ALIFMDMULTALGORITHM_H
103 #include "AliFMDMultPoisson.h" // ALIFMDMULTPOISSON_H
104 #include "AliFMDMultNaiive.h" // ALIFMDMULTNAIIVE_H
106 //____________________________________________________________________
107 ClassImp(AliFMDReconstructor);
109 //____________________________________________________________________
110 AliFMDReconstructor::AliFMDReconstructor()
111 : AliReconstructor(),
116 // Make a new FMD reconstructor object - default CTOR.
122 fAlgorithms.Add(new AliFMDMultNaiive);
123 fAlgorithms.Add(new AliFMDMultPoisson);
127 //____________________________________________________________________
128 AliFMDReconstructor::AliFMDReconstructor(const AliFMDReconstructor& other)
129 : AliReconstructor(),
135 SetPedestal(other.fPedestal, other.fPedestalWidth, other.fPedestalFactor);
137 fFMDLoader = other.fFMDLoader;
138 fRunLoader = other.fRunLoader;
141 fAlgorithms.Delete();
142 TIter next(&(other.fAlgorithms));
143 AliFMDMultAlgorithm* algorithm = 0;
144 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
145 fAlgorithms.Add(algorithm);
146 fAlgorithms.SetOwner(kFALSE);
150 //____________________________________________________________________
152 AliFMDReconstructor::operator=(const AliFMDReconstructor& other)
154 // Assignment operator
155 SetPedestal(other.fPedestal, other.fPedestalWidth, other.fPedestalFactor);
157 fFMDLoader = other.fFMDLoader;
158 fRunLoader = other.fRunLoader;
161 fAlgorithms.Delete();
162 TIter next(&(other.fAlgorithms));
163 AliFMDMultAlgorithm* algorithm = 0;
164 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
165 fAlgorithms.Add(algorithm);
166 fAlgorithms.SetOwner(kFALSE);
171 //____________________________________________________________________
172 AliFMDReconstructor::~AliFMDReconstructor()
175 fAlgorithms.Delete();
178 //____________________________________________________________________
180 AliFMDReconstructor::SetPedestal(Float_t mean, Float_t width, Float_t factor)
182 // Set the pedestal, and pedestal width
184 fPedestalWidth = width;
185 fPedestalFactor = factor;
188 //____________________________________________________________________
190 AliFMDReconstructor::Reconstruct(AliRunLoader* runLoader,
191 AliRawReader* rawReader) const
193 // Collects all digits in the same active volume into number of
196 // Reconstruct number of particles in given group of pads for given
197 // FMDvolume determined by numberOfVolume,
198 // numberOfMinSector, numberOfMaxSector, numberOfMinRing,
201 // The reconstruction method is choosen based on the number of empty
204 Error("Exec","Run Loader loader is NULL - Session not opened");
207 fRunLoader = runLoader;
208 fFMDLoader = runLoader->GetLoader("FMDLoader");
210 Fatal("AliFMDReconstructor","Can not find FMD (loader) "
211 "in specified event");
213 // Get the AliRun object
214 if (!fRunLoader->GetAliRun()) fRunLoader->LoadgAlice();
216 // Get the AliFMD object
217 fFMD = static_cast<AliFMD*>(fRunLoader->GetAliRun()->GetDetector("FMD"));
219 AliError("Can not get FMD from gAlice");
222 fFMDLoader->LoadRecPoints("RECREATE");
224 if (!fRunLoader->TreeE()) fRunLoader->LoadHeader();
226 TIter next(&fAlgorithms);
227 AliFMDMultAlgorithm* algorithm = 0;
228 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
229 algorithm->PreRun(fFMD);
233 while (rawReader->NextEvent()) {
234 ProcessEvent(event, rawReader);
239 Int_t nEvents= Int_t(fRunLoader->TreeE()->GetEntries());
240 for(Int_t event = 0; event < nEvents; event++)
241 ProcessEvent(event, 0);
246 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
247 algorithm->PostRun();
249 fFMDLoader->UnloadRecPoints();
255 //____________________________________________________________________
257 AliFMDReconstructor::Reconstruct(AliRunLoader* runLoader) const
259 // Collects all digits in the same active volume into number of
262 // Reconstruct number of particles in given group of pads for given
263 // FMDvolume determined by numberOfVolume,
264 // numberOfMinSector, numberOfMaxSector, numberOfMinRing,
267 // The reconstruction method is choosen based on the number of empty
269 Reconstruct(runLoader, 0);
273 //____________________________________________________________________
275 AliFMDReconstructor::ProcessEvent(Int_t event,
276 AliRawReader* reader) const
278 // Process one event read from either a clones array or from a a raw
280 fRunLoader->GetEvent(event) ;
281 //event z-vertex for correction eta-rad dependence
282 AliHeader *header = fRunLoader->GetHeader();
283 if (!header) Warning("ProcessEvent", "no AliHeader found!");
284 AliGenEventHeader* genHeader = (header ? header->GenEventHeader() : 0);
286 // Get the Z--coordinate from the event header
288 if (genHeader) genHeader->PrimaryVertex(o);
289 else Warning("ProcessEvent", "No GenEventHeader Found");
290 fCurrentVertex = o.At(2);
292 // If the recontruction tree isn't loaded, load it
293 if(fFMDLoader->TreeR()==0) fFMDLoader->MakeTree("R");
295 // Load or recreate the digits
296 if (fFMDLoader->LoadDigits((reader ? "UPDATE" : "READ"))) {
298 Error("Exec","Error occured while loading digits. Exiting.");
302 // Get the digits tree
303 TTree* digitTree = fFMDLoader->TreeD();
306 Error("Exec","Can not get Tree with Digits. "
307 "Nothing to reconstruct - Exiting");
310 fFMDLoader->MakeTree("D");
311 digitTree = fFMDLoader->TreeD();
314 // Get the FMD branch holding the digits.
315 TBranch *digitBranch = digitTree->GetBranch("FMD");
316 TClonesArray* digits = fFMD->Digits();
319 Error("Exec", "No digit branch for the FMD found");
322 fFMD->MakeBranchInTree(digitTree, fFMD->GetName(), &(digits), 4000, 0);
324 if (!reader) digitBranch->SetAddress(&digits);
327 AliFMDRawReader rawRead(fFMD, reader);
328 AliDebug(10, Form("Making raw reader with sample rate: %d",
329 fFMD->GetSampleRate()));
330 rawRead.SetSampleRate(fFMD->GetSampleRate());
334 // Read the ADC values from a clones array.
335 AliDebug(10, "Reading ADCs from Digits array");
336 // read Digits, and reconstruct the particles
337 if (!fFMDLoader->TreeD()->GetEvent(0)) return;
340 TIter next(&fAlgorithms);
341 AliFMDMultAlgorithm* algorithm = 0;
342 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
343 algorithm->PreEvent(fFMDLoader->TreeR(), fCurrentVertex);
345 ProcessDigits(digits);
349 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
350 algorithm->PostEvent();
354 fFMDLoader->WriteDigits("OVERWRITE");
356 fFMDLoader->UnloadDigits();
357 fFMDLoader->TreeR()->Reset();
358 fFMDLoader->TreeR()->Fill();
359 fFMDLoader->WriteRecPoints("OVERWRITE");
362 //____________________________________________________________________
364 AliFMDReconstructor::SubtractPedestal(AliFMDDigit* digit) const
366 // Member function to subtract the pedestal from a digit
367 // This implementation does nothing, but a derived class could over
368 // load this to subtract a pedestal that was given in a database or
369 // something like that.
372 Float_t ped = fPedestal + fPedestalFactor * fPedestalWidth;
373 if (digit->Count3() > 0) counts = digit->Count3();
374 else if (digit->Count2() > 0) counts = digit->Count2();
375 else counts = digit->Count1();
376 counts = TMath::Max(Int_t(counts - ped), 0);
377 return UShort_t(counts);
380 //____________________________________________________________________
382 AliFMDReconstructor::ProcessDigits(TClonesArray* digits) const
384 Int_t nDigits = digits->GetEntries();
385 for (Int_t i = 0; i < nDigits; i++) {
386 AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
387 AliFMDSubDetector* subDetector = 0;
388 switch (digit->Detector()) {
389 case 1: subDetector = fFMD->GetFMD1(); break;
390 case 2: subDetector = fFMD->GetFMD2(); break;
391 case 3: subDetector = fFMD->GetFMD3(); break;
394 Warning("ProcessDigits", "Unknown detector: FMD%d" , digit->Detector());
398 AliFMDRing* ring = 0;
400 switch(digit->Ring()) {
403 ring = subDetector->GetInner();
404 ringZ = subDetector->GetInnerZ();
408 ring = subDetector->GetOuter();
409 ringZ = subDetector->GetOuterZ();
413 Warning("ProcessDigits", "Unknown ring: FMD%d%c", digit->Detector(),
418 Float_t realZ = fCurrentVertex + ringZ;
419 Float_t stripR = ((ring->GetHighR() - ring->GetLowR())
420 / ring->GetNStrips() * (digit->Strip() + .5)
422 Float_t theta = TMath::ATan2(stripR, realZ);
423 Float_t phi = (2 * TMath::Pi() / ring->GetNSectors()
424 * (digit->Sector() + .5));
425 Float_t eta = -TMath::Log(TMath::Tan(theta / 2));
426 UShort_t counts = SubtractPedestal(digit);
428 TIter next(&fAlgorithms);
429 AliFMDMultAlgorithm* algorithm = 0;
430 while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
431 algorithm->ProcessDigit(digit, eta, phi, counts);
436 //____________________________________________________________________
438 AliFMDReconstructor::FillESD(AliRunLoader* /*fRunLoader*/,
439 AliESD* /*esd*/) const
441 // nothing to be done
445 //____________________________________________________________________