Removing unused runloader
[u/mrichter/AliRoot.git] / FMD / AliFMDReconstructor.cxx
CommitLineData
37c55dc0 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
0d0e6995 15/* $Id$ */
c2fc1258 16/** @file AliFMDReconstructor.cxx
17 @author Christian Holm Christensen <cholm@nbi.dk>
18 @date Mon Mar 27 12:47:09 2006
19 @brief FMD reconstruction
20*/
4347b38f 21//____________________________________________________________________
42403906 22//
02a27b50 23// This is a class that constructs AliFMDRecPoint objects from of Digits
6169f936 24// This class reads either digits from a TClonesArray or raw data from
4347b38f 25// a DDL file (or similar), and stores the read ADC counts in an
6169f936 26// internal cache (fAdcs). The rec-points are made via the naiive
27// method.
4347b38f 28//
37c55dc0 29//-- Authors: Evgeny Karpechev(INR) and Alla Maevsksia
4347b38f 30// Latest changes by Christian Holm Christensen <cholm@nbi.dk>
31//
32//
33//____________________________________________________________________
37c55dc0 34
f95a63c4 35// #include <AliLog.h> // ALILOG_H
6169f936 36// #include <AliRun.h> // ALIRUN_H
f95a63c4 37#include "AliFMDDebug.h"
56b1929b 38#include <AliRunLoader.h> // ALIRUNLOADER_H
56b1929b 39#include <AliHeader.h> // ALIHEADER_H
56b1929b 40#include <AliGenEventHeader.h> // ALIGENEVENTHEADER_H
1a1fdef7 41#include "AliFMDGeometry.h" // ALIFMDGEOMETRY_H
42#include "AliFMDParameters.h" // ALIFMDPARAMETERS_H
56b1929b 43#include "AliFMDDigit.h" // ALIFMDDIGIT_H
44#include "AliFMDReconstructor.h" // ALIFMDRECONSTRUCTOR_H
56b1929b 45#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
bf000c32 46#include "AliFMDRecPoint.h" // ALIFMDMULTNAIIVE_H
af885e0f 47#include "AliESDEvent.h" // ALIESDEVENT_H
aad72f45 48#include "AliESDVertex.h" // ALIESDVERTEX_H
8f6ee336 49#include <AliESDFMD.h> // ALIESDFMD_H
aad72f45 50#include <TMath.h>
9684be2f 51#include <TH1.h>
52#include <TH2.h>
53#include <TFile.h>
02a27b50 54class AliRawReader;
4347b38f 55
56//____________________________________________________________________
925e6570 57ClassImp(AliFMDReconstructor)
1a1fdef7 58#if 0
59 ; // This is here to keep Emacs for indenting the next line
60#endif
37c55dc0 61
4347b38f 62//____________________________________________________________________
63AliFMDReconstructor::AliFMDReconstructor()
64 : AliReconstructor(),
e9fd1e20 65 fMult(0x0),
66 fNMult(0),
67 fTreeR(0x0),
68 fCurrentVertex(0),
69 fESDObj(0x0),
9684be2f 70 fNoiseFactor(0),
71 fAngleCorrect(kTRUE),
72 fVertexType(kNoVertex),
73 fRunLoader(0x0),
74 fESD(0x0),
75 fDiagnostics(kFALSE),
76 fDiagStep1(0),
77 fDiagStep2(0),
78 fDiagStep3(0),
79 fDiagStep4(0),
80 fDiagAll(0)
4347b38f 81{
8f6ee336 82 // Make a new FMD reconstructor object - default CTOR.
a9579262 83 SetNoiseFactor();
84 SetAngleCorrect();
4347b38f 85}
86
42403906 87
88//____________________________________________________________________
89AliFMDReconstructor::AliFMDReconstructor(const AliFMDReconstructor& other)
8f6ee336 90 : AliReconstructor(),
91 fMult(other.fMult),
e9fd1e20 92 fNMult(other.fNMult),
93 fTreeR(other.fTreeR),
94 fCurrentVertex(other.fCurrentVertex),
95 fESDObj(other.fESDObj),
a9579262 96 fNoiseFactor(other.fNoiseFactor),
9684be2f 97 fAngleCorrect(other.fAngleCorrect),
98 fVertexType(other.fVertexType),
99 fRunLoader(other.fRunLoader),
100 fESD(other.fESD),
101 fDiagnostics(other.fDiagnostics),
102 fDiagStep1(other.fDiagStep1),
103 fDiagStep2(other.fDiagStep2),
104 fDiagStep3(other.fDiagStep3),
105 fDiagStep4(other.fDiagStep4),
106 fDiagAll(other.fDiagAll)
42403906 107{
56b1929b 108 // Copy constructor
42403906 109}
110
111
112//____________________________________________________________________
113AliFMDReconstructor&
114AliFMDReconstructor::operator=(const AliFMDReconstructor& other)
115{
56b1929b 116 // Assignment operator
a9579262 117 fMult = other.fMult;
118 fNMult = other.fNMult;
119 fTreeR = other.fTreeR;
e9fd1e20 120 fCurrentVertex = other.fCurrentVertex;
a9579262 121 fESDObj = other.fESDObj;
a9579262 122 fNoiseFactor = other.fNoiseFactor;
123 fAngleCorrect = other.fAngleCorrect;
9684be2f 124 fVertexType = other.fVertexType;
125 fRunLoader = other.fRunLoader;
126 fESD = other.fESD;
127 fDiagnostics = other.fDiagnostics;
128 fDiagStep1 = other.fDiagStep1;
129 fDiagStep2 = other.fDiagStep2;
130 fDiagStep3 = other.fDiagStep3;
131 fDiagStep4 = other.fDiagStep4;
132 fDiagAll = other.fDiagAll;
42403906 133 return *this;
134}
56b1929b 135
136//____________________________________________________________________
137AliFMDReconstructor::~AliFMDReconstructor()
138{
139 // Destructor
8f6ee336 140 if (fMult) fMult->Delete();
141 if (fMult) delete fMult;
142 if (fESDObj) delete fESDObj;
56b1929b 143}
4347b38f 144
145//____________________________________________________________________
146void
1a1fdef7 147AliFMDReconstructor::Init(AliRunLoader* runLoader)
148{
149 // Initialize the reconstructor
f95a63c4 150 AliFMDDebug(2, ("Init called with runloader 0x%x", runLoader));
bf000c32 151 // Initialize the geometry
9b48326f 152 AliFMDGeometry* geom = AliFMDGeometry::Instance();
153 geom->Init();
154 geom->InitTransformations();
155
156 // Initialize the parameters
157 AliFMDParameters* param = AliFMDParameters::Instance();
158 param->Init();
bf000c32 159
8f6ee336 160 // Current vertex position
1a1fdef7 161 fCurrentVertex = 0;
8f6ee336 162 // Create array of reconstructed strip multiplicities
bf000c32 163 fMult = new TClonesArray("AliFMDRecPoint", 51200);
8f6ee336 164 // Create ESD output object
165 fESDObj = new AliESDFMD;
166
167 // Check that we have a run loader
9684be2f 168 fRunLoader = runLoader;
8f6ee336 169
9684be2f 170 // Check if we need diagnostics histograms
171 if (!fDiagnostics) return;
172 fDiagStep1 = new TH2I("diagStep1", "Read ADC vs. Noise surpressed ADC",
173 1024, -.5, 1023.5, 1024, -.5, 1023.5);
174 fDiagStep1->SetDirectory(0);
175 fDiagStep1->GetXaxis()->SetTitle("ADC (read)");
176 fDiagStep1->GetYaxis()->SetTitle(Form("ADC (noise surpressed %4.f)",
177 fNoiseFactor));
178 fDiagStep2 = new TH2F("diagStep2", "ADC vs Edep deduced",
179 1024, -.5, 1023.5, 100, 0, 2);
180 fDiagStep2->SetDirectory(0);
181 fDiagStep2->GetXaxis()->SetTitle("ADC (noise surpressed)");
182 fDiagStep2->GetYaxis()->SetTitle("#Delta E [GeV]");
183 fDiagStep3 = new TH2F("diagStep3", "Edep vs Edep path corrected",
184 100, 0., 2., 100, 0., 2.);
185 fDiagStep3->SetDirectory(0);
186 fDiagStep3->GetXaxis()->SetTitle("#Delta E [GeV]");
187 fDiagStep3->GetYaxis()->SetTitle("#Delta E/#Delta x #times #delta x [GeV]");
188 fDiagStep4 = new TH2F("diagStep4", "Edep vs Multiplicity deduced",
189 100, 0., 2., 100, -.1, 19.9);
190 fDiagStep4->SetDirectory(0);
191 fDiagStep4->GetXaxis()->SetTitle("#Delta E/#Delta x #times #delta x [GeV]");
192 fDiagStep4->GetYaxis()->SetTitle("Multiplicity");
193 fDiagAll = new TH2F("diagAll", "Read ADC vs Multiplicity deduced",
194 1024, -.5, 1023.5, 100, -.1, 19.9);
195 fDiagAll->SetDirectory(0);
196 fDiagAll->GetXaxis()->SetTitle("ADC (read)");
197 fDiagAll->GetYaxis()->SetTitle("Multiplicity");
4347b38f 198}
dc8af42e 199
4347b38f 200//____________________________________________________________________
201void
1a1fdef7 202AliFMDReconstructor::ConvertDigits(AliRawReader* reader,
203 TTree* digitsTree) const
204{
205 // Convert Raw digits to AliFMDDigit's in a tree
f95a63c4 206 AliFMDDebug(2, ("Reading raw data into digits tree"));
1a1fdef7 207 AliFMDRawReader rawRead(reader, digitsTree);
208 // rawRead.SetSampleRate(fFMD->GetSampleRate());
209 rawRead.Exec();
4347b38f 210}
211
4347b38f 212//____________________________________________________________________
213void
9684be2f 214AliFMDReconstructor::GetVertex() const
4347b38f 215{
9684be2f 216 fVertexType = kNoVertex;
217 fCurrentVertex = 0;
a3537838 218 if (fESD) {
8f6ee336 219 const AliESDVertex* vertex = fESD->GetVertex();
220 if (vertex) {
f95a63c4 221 AliFMDDebug(2, ("Got vertex from ESD: %f", vertex->GetZv()));
8f6ee336 222 fCurrentVertex = vertex->GetZv();
9684be2f 223 fVertexType = kESDVertex;
224 return;
8f6ee336 225 }
a3537838 226 }
9684be2f 227 // Check if we can get the header tree
228 AliGenEventHeader* genHeader = ((!fRunLoader ||
229 !fRunLoader->GetHeader() ||
230 !fRunLoader->GetHeader()->GenEventHeader())
231 ? 0
232 : fRunLoader->GetHeader()->GenEventHeader());
233 if (genHeader) {
234 TArrayF vtx;
235 genHeader->PrimaryVertex(vtx);
236 fCurrentVertex = vtx[2];
237 fVertexType = kGenVertex;
f95a63c4 238 AliFMDDebug(2, ("Got vertex from generator: %f", fCurrentVertex));
9684be2f 239 AliWarning("Got vertex from generator event header");
240 return;
241 }
242 AliWarning("Didn't get any vertex from ESD or generator");
243}
244
245
246//____________________________________________________________________
247void
248AliFMDReconstructor::Reconstruct(TTree* digitsTree,
249 TTree* clusterTree) const
250{
251 // Reconstruct event from digits in tree
252 // Get the FMD branch holding the digits.
253 // FIXME: The vertex may not be known yet, so we may have to move
254 // some of this to FillESD.
f95a63c4 255 AliFMDDebug(2, ("Reconstructing from digits in a tree"));
9684be2f 256 GetVertex();
257
1a1fdef7 258 TBranch *digitBranch = digitsTree->GetBranch("FMD");
4347b38f 259 if (!digitBranch) {
1a1fdef7 260 Error("Exec", "No digit branch for the FMD found");
261 return;
4347b38f 262 }
0abe7182 263 TClonesArray* digits = new TClonesArray("AliFMDDigit");
1a1fdef7 264 digitBranch->SetAddress(&digits);
4347b38f 265
8f6ee336 266 if (fMult) fMult->Clear();
267 if (fESDObj) fESDObj->Clear();
268
269 fNMult = 0;
270 fTreeR = clusterTree;
271 fTreeR->Branch("FMD", &fMult);
272
f95a63c4 273 AliFMDDebug(5, ("Getting entry 0 from digit branch"));
1a1fdef7 274 digitBranch->GetEntry(0);
275
f95a63c4 276 AliFMDDebug(5, ("Processing digits"));
e802be3e 277 ProcessDigits(digits);
278
8f6ee336 279 Int_t written = clusterTree->Fill();
f95a63c4 280 AliFMDDebug(10, ("Filled %d bytes into cluster tree", written));
0abe7182 281 digits->Delete();
282 delete digits;
4347b38f 283}
1a1fdef7 284
8f6ee336 285
4347b38f 286//____________________________________________________________________
e802be3e 287void
288AliFMDReconstructor::ProcessDigits(TClonesArray* digits) const
4347b38f 289{
69b696b9 290 // For each digit, find the pseudo rapdity, azimuthal angle, and
291 // number of corrected ADC counts, and pass it on to the algorithms
292 // used.
e802be3e 293 Int_t nDigits = digits->GetEntries();
f95a63c4 294 AliFMDDebug(1, ("Got %d digits", nDigits));
a9579262 295 fESDObj->SetNoiseFactor(fNoiseFactor);
296 fESDObj->SetAngleCorrected(fAngleCorrect);
e802be3e 297 for (Int_t i = 0; i < nDigits; i++) {
298 AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
8f6ee336 299 AliFMDParameters* param = AliFMDParameters::Instance();
300 // Check that the strip is not marked as dead
301 if (param->IsDead(digit->Detector(), digit->Ring(),
15b17c89 302 digit->Sector(), digit->Strip())) {
f95a63c4 303 AliFMDDebug(10, ("FMD%d%c[%2d,%3d] is dead", digit->Detector(),
15b17c89 304 digit->Ring(), digit->Sector(), digit->Strip()));
305 continue;
306 }
8f6ee336 307
308 // digit->Print();
309 // Get eta and phi
310 Float_t eta, phi;
311 PhysicalCoordinates(digit, eta, phi);
4347b38f 312
8f6ee336 313 // Substract pedestal.
e802be3e 314 UShort_t counts = SubtractPedestal(digit);
4347b38f 315
8f6ee336 316 // Gain match digits.
317 Double_t edep = Adc2Energy(digit, eta, counts);
318
319 // Make rough multiplicity
320 Double_t mult = Energy2Multiplicity(digit, edep);
321
f95a63c4 322 AliFMDDebug(10, ("FMD%d%c[%2d,%3d]: "
8f6ee336 323 "ADC: %d, Counts: %d, Energy: %f, Mult: %f",
324 digit->Detector(), digit->Ring(), digit->Sector(),
325 digit->Strip(), digit->Counts(), counts, edep, mult));
326
327 // Create a `RecPoint' on the output branch.
bf000c32 328 AliFMDRecPoint* m =
329 new ((*fMult)[fNMult]) AliFMDRecPoint(digit->Detector(),
330 digit->Ring(),
331 digit->Sector(),
332 digit->Strip(),
333 eta, phi,
334 edep, mult);
8f6ee336 335 (void)m; // Suppress warnings about unused variables.
336 fNMult++;
337
338 fESDObj->SetMultiplicity(digit->Detector(), digit->Ring(),
339 digit->Sector(), digit->Strip(), mult);
340 fESDObj->SetEta(digit->Detector(), digit->Ring(),
341 digit->Sector(), digit->Strip(), eta);
9684be2f 342
343 if (fDiagAll) fDiagAll->Fill(digit->Counts(), mult);
4347b38f 344 }
4347b38f 345}
8f6ee336 346
1a1fdef7 347//____________________________________________________________________
348UShort_t
349AliFMDReconstructor::SubtractPedestal(AliFMDDigit* digit) const
350{
351 // Member function to subtract the pedestal from a digit
352 // This implementation does nothing, but a derived class could over
353 // load this to subtract a pedestal that was given in a database or
354 // something like that.
355
8f6ee336 356 Int_t counts = 0;
9684be2f 357 Int_t adc = 0;
8f6ee336 358 AliFMDParameters* param = AliFMDParameters::Instance();
a9579262 359 Float_t ped = param->GetPedestal(digit->Detector(),
8f6ee336 360 digit->Ring(),
361 digit->Sector(),
362 digit->Strip());
a9579262 363 Float_t noise = param->GetPedestalWidth(digit->Detector(),
364 digit->Ring(),
365 digit->Sector(),
366 digit->Strip());
f95a63c4 367 AliFMDDebug(15, ("Subtracting pedestal %f from signal %d",
a9579262 368 ped, digit->Counts()));
9684be2f 369 if (digit->Count3() > 0) adc = digit->Count3();
370 else if (digit->Count2() > 0) adc = digit->Count2();
371 else adc = digit->Count1();
372 counts = TMath::Max(Int_t(adc - ped), 0);
a9579262 373 if (counts < noise * fNoiseFactor) counts = 0;
f95a63c4 374 if (counts > 0) AliFMDDebug(15, ("Got a hit strip"));
9684be2f 375 if (fDiagStep1) fDiagStep1->Fill(adc, counts);
8f6ee336 376
1a1fdef7 377 return UShort_t(counts);
378}
379
4347b38f 380//____________________________________________________________________
8f6ee336 381Float_t
382AliFMDReconstructor::Adc2Energy(AliFMDDigit* digit,
a9579262 383 Float_t eta,
8f6ee336 384 UShort_t count) const
385{
386 // Converts number of ADC counts to energy deposited.
387 // Note, that this member function can be overloaded by derived
388 // classes to do strip-specific look-ups in databases or the like,
389 // to find the proper gain for a strip.
390 //
391 // In this simple version, we calculate the energy deposited as
392 //
393 // EnergyDeposited = cos(theta) * gain * count
394 //
395 // where
396 //
397 // Pre_amp_MIP_Range
398 // gain = ----------------- * Energy_deposited_per_MIP
399 // ADC_channel_size
400 //
401 // is constant and the same for all strips.
a9579262 402 if (count <= 0) return 0;
8f6ee336 403 AliFMDParameters* param = AliFMDParameters::Instance();
404 Float_t gain = param->GetPulseGain(digit->Detector(),
405 digit->Ring(),
406 digit->Sector(),
407 digit->Strip());
f95a63c4 408 AliFMDDebug(15, ("Converting counts %d to energy via factor %f",
8f6ee336 409 count, gain));
a9579262 410
411 Double_t edep = count * gain;
9684be2f 412 if (fDiagStep2) fDiagStep2->Fill(count, edep);
a9579262 413 if (fAngleCorrect) {
9684be2f 414 Double_t theta = 2 * TMath::ATan(TMath::Exp(-eta));
415 Double_t corr = TMath::Abs(TMath::Cos(theta));
416 Double_t cedep = corr * edep;
f95a63c4 417 AliFMDDebug(10, ("correcting for path %f * %f = %f (eta=%f, theta=%f)",
9684be2f 418 edep, corr, cedep, eta, theta));
419 if (fDiagStep3) fDiagStep3->Fill(edep, cedep);
420 edep = cedep;
a9579262 421 }
8f6ee336 422 return edep;
423}
424
425//____________________________________________________________________
426Float_t
427AliFMDReconstructor::Energy2Multiplicity(AliFMDDigit* /* digit */,
428 Float_t edep) const
429{
430 // Converts an energy signal to number of particles.
431 // Note, that this member function can be overloaded by derived
432 // classes to do strip-specific look-ups in databases or the like,
433 // to find the proper gain for a strip.
434 //
435 // In this simple version, we calculate the multiplicity as
436 //
437 // multiplicity = Energy_deposited / Energy_deposited_per_MIP
438 //
439 // where
440 //
441 // Energy_deposited_per_MIP = 1.664 * SI_density * SI_thickness
442 //
443 // is constant and the same for all strips
444 AliFMDParameters* param = AliFMDParameters::Instance();
445 Double_t edepMIP = param->GetEdepMip();
446 Float_t mult = edep / edepMIP;
447 if (edep > 0)
f95a63c4 448 AliFMDDebug(15, ("Translating energy %f to multiplicity via "
8f6ee336 449 "divider %f->%f", edep, edepMIP, mult));
9684be2f 450 if (fDiagStep4) fDiagStep4->Fill(edep, mult);
8f6ee336 451 return mult;
452}
453
454//____________________________________________________________________
455void
456AliFMDReconstructor::PhysicalCoordinates(AliFMDDigit* digit,
457 Float_t& eta,
458 Float_t& phi) const
459{
460 // Get the eta and phi of a digit
461 //
462 // Get geometry.
9b48326f 463 AliFMDGeometry* geom = AliFMDGeometry::Instance();
bf000c32 464 Double_t x, y, z, r, theta;
9b48326f 465 geom->Detector2XYZ(digit->Detector(), digit->Ring(), digit->Sector(),
bf000c32 466 digit->Strip(), x, y, z);
467 // Correct for vertex offset.
468 z += fCurrentVertex;
469 phi = TMath::ATan2(y, x);
470 r = TMath::Sqrt(y * y + x * x);
471 theta = TMath::ATan2(r, z);
472 eta = -TMath::Log(TMath::Tan(theta / 2));
8f6ee336 473}
474
475
476
477//____________________________________________________________________
4347b38f 478void
1a1fdef7 479AliFMDReconstructor::FillESD(TTree* /* digitsTree */,
480 TTree* /* clusterTree */,
af885e0f 481 AliESDEvent* esd) const
121a60bd 482{
42403906 483 // nothing to be done
69b696b9 484 // FIXME: The vertex may not be known when Reconstruct is executed,
485 // so we may have to move some of that member function here.
f95a63c4 486 AliFMDDebug(2, ("Calling FillESD with two trees and one ESD"));
8f6ee336 487 // fESDObj->Print();
488
489 if (esd) {
f95a63c4 490 AliFMDDebug(2, ("Writing FMD data to ESD tree"));
8f6ee336 491 esd->SetFMDData(fESDObj);
a3537838 492 }
9684be2f 493
494 if (!fDiagnostics || !esd) return;
495 static bool first = true;
31fd97b2 496 Int_t evno = esd->GetEventNumberInFile(); // This is most likely NOT the event number you'd like to use. It has nothing to do with the 'real' event number.
f95a63c4 497 AliFMDDebug(1, ("Writing diagnostics histograms to FMD.Diag.root/%03d",
9684be2f 498 evno));
499 TFile f("FMD.Diag.root", (first ? "RECREATE" : "UPDATE"));
500 first = false;
501 f.cd();
502 TDirectory* d = f.mkdir(Form("%03d", evno),
503 Form("Diagnostics histograms for event # %d", evno));
504 d->cd();
505 if (fDiagStep1) fDiagStep1->Write();
506 if (fDiagStep2) fDiagStep2->Write();
507 if (fDiagStep3) fDiagStep3->Write();
508 if (fDiagStep4) fDiagStep4->Write();
509 if (fDiagAll) fDiagAll->Write();
510 d->Write();
511 f.Write();
512 f.Close();
513
514 if (fDiagStep1) fDiagStep1->Reset();
515 if (fDiagStep2) fDiagStep2->Reset();
516 if (fDiagStep3) fDiagStep3->Reset();
517 if (fDiagStep4) fDiagStep4->Reset();
518 if (fDiagAll) fDiagAll->Reset();
121a60bd 519}
520
54e415a8 521
522//____________________________________________________________________
523void
524AliFMDReconstructor::Reconstruct(AliRawReader*,TTree*) const
525{
526 // Cannot be used. See member function with same name but with 2
527 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 528 AliFMDDebug(2, ("Calling FillESD with loader and tree"));
54e415a8 529 AliError("MayNotUse");
530}
531//____________________________________________________________________
532void
533AliFMDReconstructor::Reconstruct(AliRunLoader*) const
534{
535 // Cannot be used. See member function with same name but with 2
536 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 537 AliFMDDebug(2, ("Calling FillESD with loader"));
54e415a8 538 AliError("MayNotUse");
539}
540//____________________________________________________________________
541void
542AliFMDReconstructor::Reconstruct(AliRunLoader*, AliRawReader*) const
543{
544 // Cannot be used. See member function with same name but with 2
545 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 546 AliFMDDebug(2, ("Calling FillESD with loader and raw reader"));
54e415a8 547 AliError("MayNotUse");
548}
549//____________________________________________________________________
550void
af885e0f 551AliFMDReconstructor::FillESD(AliRawReader*,TTree*,AliESDEvent*) const
54e415a8 552{
553 // Cannot be used. See member function with same name but with 2
554 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 555 AliFMDDebug(2, ("Calling FillESD with raw reader, tree, and ESD"));
54e415a8 556 AliError("MayNotUse");
557}
558//____________________________________________________________________
559void
af885e0f 560AliFMDReconstructor::FillESD(AliRunLoader*,AliESDEvent*) const
54e415a8 561{
562 // Cannot be used. See member function with same name but with 2
563 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 564 AliFMDDebug(2, ("Calling FillESD with loader and ESD"));
54e415a8 565 AliError("MayNotUse");
566}
567//____________________________________________________________________
568void
af885e0f 569AliFMDReconstructor::FillESD(AliRunLoader*,AliRawReader*,AliESDEvent*) const
54e415a8 570{
571 // Cannot be used. See member function with same name but with 2
572 // TTree arguments. Make sure you do local reconstrucion
f95a63c4 573 AliFMDDebug(2, ("Calling FillESD with loader, raw reader, and ESD"));
54e415a8 574 AliError("MayNotUse");
575}
576
42403906 577//____________________________________________________________________
578//
579// EOF
580//