[u/mrichter/AliRoot.git] / FMD / AliFMDReconstructor.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

//____________________________________________________________________
//
// This is a class that constructs AliFMDMult (reconstructed
// multiplicity) from of Digits
//
// This class reads either digits from a TClonesArray or raw data from
// a DDL file (or similar), and stores the read ADC counts in an
// internal cache (fAdcs). 
//
// From the cached values it then calculates the number of particles
// that hit a region of the FMDs, as specified by the user. 
//
// The reconstruction can be done in two ways: Either via counting the
// number of empty strips (Poisson method), or by converting the ADC
// signal to an energy deposition, and then dividing by the typical
// energy loss of a particle.
// 
//      +---------------------+       +---------------------+
//      | AliFMDReconstructor |<>-----| AliFMDMultAlgorithm |
//      +---------------------+       +---------------------+
//                                               ^
//                                               |
//                                   +-----------+---------+
//                                   |                     |
//                         +-------------------+   +------------------+
//                         | AliFMDMultPoisson |   | AliFMDMultNaiive |
//                         +-------------------+   +------------------+
//
// AliFMDReconstructor acts as a manager class.  It contains a list of
// AliFMDMultAlgorithm objects.  The call graph looks something like 
//
//
//       +----------------------+            +----------------------+
//       | :AliFMDReconstructor |            | :AliFMDMultAlgorithm |
//       +----------------------+            +----------------------+
//                  |                                  |
//    Reconstruct  +-+                                 |
//    ------------>| |                         PreRun +-+
//                 | |------------------------------->| |   
//                 | |                                +-+
//                 | |-----+ (for each event)          |
//                 | |     | *ProcessEvent             |
//                 |+-+    |                           |
//                 || |<---+                 PreEvent +-+
//                 || |------------------------------>| |      
//                 || |                               +-+
//                 || |-----+                          |
//                 || |     | ProcessDigits            |
//                 ||+-+    |                          |
//                 ||| |<---+                          |
//                 ||| |         *ProcessDigit(digit) +-+
//                 ||| |----------------------------->| |
//                 ||| |                              +-+
//                 ||+-+                               |
//                 || |                     PostEvent +-+
//                 || |------------------------------>| |
//                 || |                               +-+
//                 |+-+                                |
//                 | |                        PostRun +-+
//                 | |------------------------------->| |
//                 | |                                +-+
//                 +-+                                 |
//                  |                                  |
//
//
// 
//-- Authors: Evgeny Karpechev(INR) and Alla Maevsksia
//  Latest changes by Christian Holm Christensen <cholm@nbi.dk>
//
//
//____________________________________________________________________

#include <AliLog.h>                        // ALILOG_H
#include <AliRun.h>                        // ALIRUN_H
#include <AliRunLoader.h>                  // ALIRUNLOADER_H
#include <AliLoader.h>                     // ALILOADER_H
#include <AliHeader.h>                     // ALIHEADER_H
#include <AliRawReader.h>                  // ALIRAWREADER_H
#include <AliGenEventHeader.h>             // ALIGENEVENTHEADER_H
#include "AliFMD.h"         		   // ALIFMD_H
#include "AliFMDGeometry.h"                // ALIFMDGEOMETRY_H
#include "AliFMDParameters.h"              // ALIFMDPARAMETERS_H
#include "AliFMDDetector.h"                // ALIFMDDETECTOR_H
#include "AliFMDRing.h"                    // ALIFMDRING_H
#include "AliFMDDigit.h"                   // ALIFMDDIGIT_H
#include "AliFMDReconstructor.h"           // ALIFMDRECONSTRUCTOR_H
#include "AliFMDRawStream.h"               // ALIFMDRAWSTREAM_H
#include "AliFMDRawReader.h"               // ALIFMDRAWREADER_H
#include "AliFMDMultAlgorithm.h" 	   // ALIFMDMULTALGORITHM_H
#include "AliFMDMultPoisson.h"		   // ALIFMDMULTPOISSON_H
#include "AliFMDMultNaiive.h"		   // ALIFMDMULTNAIIVE_H
#include "AliESD.h"			   // ALIESD_H

//____________________________________________________________________
ClassImp(AliFMDReconstructor)
#if 0
  ; // This is here to keep Emacs for indenting the next line
#endif

//____________________________________________________________________
AliFMDReconstructor::AliFMDReconstructor() 
  : AliReconstructor(),
    fPedestal(0), 
    fPedestalWidth(0),
    fPedestalFactor(0)
{
  // Make a new FMD reconstructor object - default CTOR.
  AliFMDParameters* pars = AliFMDParameters::Instance();
  fPedestal       = pars->GetPedestal();
  fPedestalWidth  = pars->GetPedestalWidth();
  fPedestalFactor = pars->GetPedestalFactor();
  
  fAlgorithms.Add(new AliFMDMultNaiive);
  fAlgorithms.Add(new AliFMDMultPoisson);

  TIter next(&fAlgorithms);
  AliFMDMultAlgorithm* algorithm = 0;
  while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
    algorithm->PreRun(0);
}
  

//____________________________________________________________________
AliFMDReconstructor::AliFMDReconstructor(const AliFMDReconstructor& other) 
  : AliReconstructor(),
    fPedestal(0), 
    fPedestalWidth(0),
    fPedestalFactor(0)
{
  // Copy constructor 
  fPedestal       = other.fPedestal;
  fPedestalWidth  = other.fPedestalWidth;
  fPedestalFactor = other.fPedestalFactor;

  
  fAlgorithms.Delete();
  TIter next(&(other.fAlgorithms));
  AliFMDMultAlgorithm* algorithm = 0;
  while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
    fAlgorithms.Add(algorithm);
  fAlgorithms.SetOwner(kFALSE);
}
  

//____________________________________________________________________
AliFMDReconstructor&
AliFMDReconstructor::operator=(const AliFMDReconstructor& other) 
{
  // Assignment operator
  fPedestal       = other.fPedestal;
  fPedestalWidth  = other.fPedestalWidth;
  fPedestalFactor = other.fPedestalFactor;

  fAlgorithms.Delete();
  TIter next(&(other.fAlgorithms));
  AliFMDMultAlgorithm* algorithm = 0;
  while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
    fAlgorithms.Add(algorithm);
  fAlgorithms.SetOwner(kFALSE);

  return *this;
}

//____________________________________________________________________
AliFMDReconstructor::~AliFMDReconstructor() 
{
  // Destructor 
  fAlgorithms.Delete();
}

//____________________________________________________________________
void 
AliFMDReconstructor::Init(AliRunLoader* runLoader) 
{
  // Initialize the reconstructor 
  AliDebug(1, Form("Init called with runloader 0x%x", runLoader));
  fCurrentVertex = 0;
  if (!runLoader) { 
    Warning("Init", "No run loader");
    return;
  }
  AliHeader* header = runLoader->GetHeader();
  if (!header) {
    Warning("Init", "No header");
    return;
  }
  AliGenEventHeader* eventHeader = header->GenEventHeader();
  if (eventHeader) {
    TArrayF vtx;
    eventHeader->PrimaryVertex(vtx);
    fCurrentVertex = vtx[2];
    AliDebug(1, Form("Primary vertex Z coordinate for event # %d/%d is %f", 
		     header->GetRun(), header->GetEvent(), fCurrentVertex));
    Warning("Init", "no generator event header");
  }
  else {
    Warning("Init", "No generator event header - "
	    "perhaps we get the vertex from ESD?");
  }
  // Get the ESD tree 
  SetESD(new AliESD);
}

//____________________________________________________________________
void 
AliFMDReconstructor::ConvertDigits(AliRawReader* reader, 
				   TTree* digitsTree) const
{
  // Convert Raw digits to AliFMDDigit's in a tree 
  AliFMDRawReader rawRead(reader, digitsTree);
  // rawRead.SetSampleRate(fFMD->GetSampleRate());
  rawRead.Exec();
}

//____________________________________________________________________
void 
AliFMDReconstructor::Reconstruct(TTree* digitsTree, 
				 TTree* clusterTree) const 
{
  // Reconstruct event from digits in tree 
  // Get the FMD branch holding the digits. 
  // FIXME: The vertex may not be known yet, so we may have to move
  // some of this to FillESD. 
  AliDebug(1, "Reconstructing from digits in a tree");
#if 0
  if (fESD) {
    // const AliESDVertex* vertex = fESD->GetVertex();
    // if (vertex) {
    //   AliDebug(1, Form("Got vertex from ESD: %f", vertex->GetZv()));
    //   fCurrentVertex = vertex->GetZv();
    // }
  }
#endif  
  TBranch *digitBranch = digitsTree->GetBranch("FMD");
  if (!digitBranch) {
    Error("Exec", "No digit branch for the FMD found");
    return;
  }
  TClonesArray* digits = new TClonesArray("AliFMDDigit");
  digitBranch->SetAddress(&digits);

  TIter next(&fAlgorithms);
  AliFMDMultAlgorithm* algorithm = 0;
  while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
    algorithm->PreEvent(clusterTree, fCurrentVertex);
  digitBranch->GetEntry(0);
  
  ProcessDigits(digits);

  next.Reset();
  algorithm = 0;
  while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
    algorithm->PostEvent();
  clusterTree->Fill();
  digits->Delete();
  delete digits;
}
 
//____________________________________________________________________
void
AliFMDReconstructor::ProcessDigits(TClonesArray* digits) const
{
  // For each digit, find the pseudo rapdity, azimuthal angle, and
  // number of corrected ADC counts, and pass it on to the algorithms
  // used. 
  Int_t nDigits = digits->GetEntries();
  AliDebug(1, Form("Got %d digits", nDigits));
  for (Int_t i = 0; i < nDigits; i++) {
    AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
    AliFMDGeometry* fmd = AliFMDGeometry::Instance();
    AliFMDDetector* subDetector = fmd->GetDetector(digit->Detector());
    if (!subDetector) { 
      Warning("ProcessDigits", "Unknown detector: FMD%d" , digit->Detector());
      continue;
    }
    
    AliFMDRing* ring  = subDetector->GetRing(digit->Ring());
    Float_t     ringZ = subDetector->GetRingZ(digit->Ring());
    if (!ring) {
      Warning("ProcessDigits", "Unknown ring: FMD%d%c", digit->Detector(), 
	      digit->Ring());
      break;
    }
    
    Float_t  realZ    = fCurrentVertex + ringZ;
    Float_t  stripR   = ((ring->GetHighR() - ring->GetLowR()) 
			 / ring->GetNStrips() * (digit->Strip() + .5) 
			 + ring->GetLowR());
    Float_t  theta    = TMath::ATan2(stripR, realZ);
    Float_t  phi      = (2 * TMath::Pi() / ring->GetNSectors() 
			 * (digit->Sector() + .5));
    Float_t  eta      = -TMath::Log(TMath::Tan(theta / 2));
    UShort_t counts   = SubtractPedestal(digit);
    
    TIter next(&fAlgorithms);
    AliFMDMultAlgorithm* algorithm = 0;
    while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next()))) 
      algorithm->ProcessDigit(digit, eta, phi, counts);
  }
}
      
//____________________________________________________________________
UShort_t
AliFMDReconstructor::SubtractPedestal(AliFMDDigit* digit) const
{
  // Member function to subtract the pedestal from a digit
  // This implementation does nothing, but a derived class could over
  // load this to subtract a pedestal that was given in a database or
  // something like that. 

  Int_t counts = 0;
  Float_t ped = fPedestal + fPedestalFactor * fPedestalWidth;
  if (digit->Count3() > 0)      counts = digit->Count3();
  else if (digit->Count2() > 0) counts = digit->Count2();
  else                          counts = digit->Count1();
  counts = TMath::Max(Int_t(counts - ped), 0);
  return  UShort_t(counts);
}

//____________________________________________________________________
void 
AliFMDReconstructor::FillESD(TTree*  /* digitsTree */, 
			     TTree*  /* clusterTree */,
			     AliESD* /* esd*/) const
{
  // nothing to be done
  // FIXME: The vertex may not be known when Reconstruct is executed,
  // so we may have to move some of that member function here. 
#if 0
  TClonesArray* multStrips  = 0;
  TClonesArray* multRegions = 0;
  TTree*        treeR  = fmdLoader->TreeR();
  TBranch*      branchRegions = treeR->GetBranch("FMDPoisson");
  TBranch*      branchStrips  = treeR->GetBranch("FMDNaiive");
  branchRegions->SetAddress(&multRegions);
  branchStrips->SetAddress(&multStrips);
  
  Int_t total = 0;
  Int_t nEntries  = clusterTree->GetEntries();
  for (Int_t entry = 0; entry < nEntries; entry++) {
    AliDebug(5, Form("Entry # %d in cluster tree", entry));
    treeR->GetEntry(entry);
    
    
    Int_t nMults = multRegions->GetLast();
    for (Int_t i = 0; i <= nMults; i++) {
      AliFMDMultRegion* multR =
	static_cast<AliFMDMultRegion*>(multRegions->UncheckedAt(i));
      Int_t nParticles=multR->Particles();
      if (i>=0 && i<=13)   hEtaPoissonI1->AddBinContent(i+1,nParticles);
      if (i>=14 && i<=27 ) hEtaPoissonI2->AddBinContent(i-13,nParticles);
      if (i>=28 && i<=33 );
      if (i>=34 && i<=47 ) hEtaPoissonI3->AddBinContent(48-i,nParticles);
      if (i>=48 && i<=53)  hEtaPoissonO3->AddBinContent(54-i,nParticles);
    }
  }
#endif   
}


//____________________________________________________________________
void 
AliFMDReconstructor::Reconstruct(AliRawReader*,TTree*) const 
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}
//____________________________________________________________________
void 
AliFMDReconstructor::Reconstruct(AliRunLoader*) const 
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}
//____________________________________________________________________
void 
AliFMDReconstructor::Reconstruct(AliRunLoader*, AliRawReader*) const 
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}
//____________________________________________________________________
void 
AliFMDReconstructor::FillESD(AliRawReader*,TTree*,AliESD*) const 
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}
//____________________________________________________________________
void 
AliFMDReconstructor::FillESD(AliRunLoader*,AliESD*) const
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}
//____________________________________________________________________
void 
AliFMDReconstructor::FillESD(AliRunLoader*,AliRawReader*,AliESD*) const 
{
  // Cannot be used.  See member function with same name but with 2
  // TTree arguments.   Make sure you do local reconstrucion 
  AliError("MayNotUse");
}

//____________________________________________________________________
//
// EOF
//
Commit	Line	Data
	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	**************************************************************************/
	15
	16	/* $Id$ */
	17
	18	//____________________________________________________________________
	19	//
	20	// This is a class that constructs AliFMDMult (reconstructed
	21	// multiplicity) from of Digits
	22	//
	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).
	26	//
	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.
	29	//
	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.
	34	//
	35	// +---------------------+ +---------------------+
	36	// \| AliFMDReconstructor \|<>-----\| AliFMDMultAlgorithm \|
	37	// +---------------------+ +---------------------+
	38	// ^
	39	// \|
	40	// +-----------+---------+
	41	// \| \|
	42	// +-------------------+ +------------------+
	43	// \| AliFMDMultPoisson \| \| AliFMDMultNaiive \|
	44	// +-------------------+ +------------------+
	45	//
	46	// AliFMDReconstructor acts as a manager class. It contains a list of
	47	// AliFMDMultAlgorithm objects. The call graph looks something like
	48	//
	49	//
	50	// +----------------------+ +----------------------+
	51	// \| :AliFMDReconstructor \| \| :AliFMDMultAlgorithm \|
	52	// +----------------------+ +----------------------+
	53	// \| \|
	54	// Reconstruct +-+ \|
	55	// ------------>\| \| PreRun +-+
	56	// \| \|------------------------------->\| \|
	57	// \| \| +-+
	58	// \| \|-----+ (for each event) \|
	59	// \| \| \| *ProcessEvent \|
	60	// \|+-+ \| \|
	61	// \|\| \|<---+ PreEvent +-+
	62	// \|\| \|------------------------------>\| \|
	63	// \|\| \| +-+
	64	// \|\| \|-----+ \|
	65	// \|\| \| \| ProcessDigits \|
	66	// \|\|+-+ \| \|
	67	// \|\|\| \|<---+ \|
	68	// \|\|\| \| *ProcessDigit(digit) +-+
	69	// \|\|\| \|----------------------------->\| \|
	70	// \|\|\| \| +-+
	71	// \|\|+-+ \|
	72	// \|\| \| PostEvent +-+
	73	// \|\| \|------------------------------>\| \|
	74	// \|\| \| +-+
	75	// \|+-+ \|
	76	// \| \| PostRun +-+
	77	// \| \|------------------------------->\| \|
	78	// \| \| +-+
	79	// +-+ \|
	80	// \| \|
	81	//
	82	//
	83	//
	84	//-- Authors: Evgeny Karpechev(INR) and Alla Maevsksia
	85	// Latest changes by Christian Holm Christensen <cholm@nbi.dk>
	86	//
	87	//
	88	//____________________________________________________________________
	89
	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 "AliFMDGeometry.h" // ALIFMDGEOMETRY_H
	99	#include "AliFMDParameters.h" // ALIFMDPARAMETERS_H
	100	#include "AliFMDDetector.h" // ALIFMDDETECTOR_H
	101	#include "AliFMDRing.h" // ALIFMDRING_H
	102	#include "AliFMDDigit.h" // ALIFMDDIGIT_H
	103	#include "AliFMDReconstructor.h" // ALIFMDRECONSTRUCTOR_H
	104	#include "AliFMDRawStream.h" // ALIFMDRAWSTREAM_H
	105	#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
	106	#include "AliFMDMultAlgorithm.h" // ALIFMDMULTALGORITHM_H
	107	#include "AliFMDMultPoisson.h" // ALIFMDMULTPOISSON_H
	108	#include "AliFMDMultNaiive.h" // ALIFMDMULTNAIIVE_H
	109	#include "AliESD.h" // ALIESD_H
	110
	111	//____________________________________________________________________
	112	ClassImp(AliFMDReconstructor)
	113	#if 0
	114	; // This is here to keep Emacs for indenting the next line
	115	#endif
	116
	117	//____________________________________________________________________
	118	AliFMDReconstructor::AliFMDReconstructor()
	119	: AliReconstructor(),
	120	fPedestal(0),
	121	fPedestalWidth(0),
	122	fPedestalFactor(0)
	123	{
	124	// Make a new FMD reconstructor object - default CTOR.
	125	AliFMDParameters* pars = AliFMDParameters::Instance();
	126	fPedestal = pars->GetPedestal();
	127	fPedestalWidth = pars->GetPedestalWidth();
	128	fPedestalFactor = pars->GetPedestalFactor();
	129
	130	fAlgorithms.Add(new AliFMDMultNaiive);
	131	fAlgorithms.Add(new AliFMDMultPoisson);
	132
	133	TIter next(&fAlgorithms);
	134	AliFMDMultAlgorithm* algorithm = 0;
	135	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	136	algorithm->PreRun(0);
	137	}
	138
	139
	140	//____________________________________________________________________
	141	AliFMDReconstructor::AliFMDReconstructor(const AliFMDReconstructor& other)
	142	: AliReconstructor(),
	143	fPedestal(0),
	144	fPedestalWidth(0),
	145	fPedestalFactor(0)
	146	{
	147	// Copy constructor
	148	fPedestal = other.fPedestal;
	149	fPedestalWidth = other.fPedestalWidth;
	150	fPedestalFactor = other.fPedestalFactor;
	151
	152
	153	fAlgorithms.Delete();
	154	TIter next(&(other.fAlgorithms));
	155	AliFMDMultAlgorithm* algorithm = 0;
	156	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	157	fAlgorithms.Add(algorithm);
	158	fAlgorithms.SetOwner(kFALSE);
	159	}
	160
	161
	162	//____________________________________________________________________
	163	AliFMDReconstructor&
	164	AliFMDReconstructor::operator=(const AliFMDReconstructor& other)
	165	{
	166	// Assignment operator
	167	fPedestal = other.fPedestal;
	168	fPedestalWidth = other.fPedestalWidth;
	169	fPedestalFactor = other.fPedestalFactor;
	170
	171	fAlgorithms.Delete();
	172	TIter next(&(other.fAlgorithms));
	173	AliFMDMultAlgorithm* algorithm = 0;
	174	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	175	fAlgorithms.Add(algorithm);
	176	fAlgorithms.SetOwner(kFALSE);
	177
	178	return *this;
	179	}
	180
	181	//____________________________________________________________________
	182	AliFMDReconstructor::~AliFMDReconstructor()
	183	{
	184	// Destructor
	185	fAlgorithms.Delete();
	186	}
	187
	188	//____________________________________________________________________
	189	void
	190	AliFMDReconstructor::Init(AliRunLoader* runLoader)
	191	{
	192	// Initialize the reconstructor
	193	AliDebug(1, Form("Init called with runloader 0x%x", runLoader));
	194	fCurrentVertex = 0;
	195	if (!runLoader) {
	196	Warning("Init", "No run loader");
	197	return;
	198	}
	199	AliHeader* header = runLoader->GetHeader();
	200	if (!header) {
	201	Warning("Init", "No header");
	202	return;
	203	}
	204	AliGenEventHeader* eventHeader = header->GenEventHeader();
	205	if (eventHeader) {
	206	TArrayF vtx;
	207	eventHeader->PrimaryVertex(vtx);
	208	fCurrentVertex = vtx[2];
	209	AliDebug(1, Form("Primary vertex Z coordinate for event # %d/%d is %f",
	210	header->GetRun(), header->GetEvent(), fCurrentVertex));
	211	Warning("Init", "no generator event header");
	212	}
	213	else {
	214	Warning("Init", "No generator event header - "
	215	"perhaps we get the vertex from ESD?");
	216	}
	217	// Get the ESD tree
	218	SetESD(new AliESD);
	219	}
	220
	221	//____________________________________________________________________
	222	void
	223	AliFMDReconstructor::ConvertDigits(AliRawReader* reader,
	224	TTree* digitsTree) const
	225	{
	226	// Convert Raw digits to AliFMDDigit's in a tree
	227	AliFMDRawReader rawRead(reader, digitsTree);
	228	// rawRead.SetSampleRate(fFMD->GetSampleRate());
	229	rawRead.Exec();
	230	}
	231
	232	//____________________________________________________________________
	233	void
	234	AliFMDReconstructor::Reconstruct(TTree* digitsTree,
	235	TTree* clusterTree) const
	236	{
	237	// Reconstruct event from digits in tree
	238	// Get the FMD branch holding the digits.
	239	// FIXME: The vertex may not be known yet, so we may have to move
	240	// some of this to FillESD.
	241	AliDebug(1, "Reconstructing from digits in a tree");
	242	#if 0
	243	if (fESD) {
	244	// const AliESDVertex* vertex = fESD->GetVertex();
	245	// if (vertex) {
	246	// AliDebug(1, Form("Got vertex from ESD: %f", vertex->GetZv()));
	247	// fCurrentVertex = vertex->GetZv();
	248	// }
	249	}
	250	#endif
	251	TBranch *digitBranch = digitsTree->GetBranch("FMD");
	252	if (!digitBranch) {
	253	Error("Exec", "No digit branch for the FMD found");
	254	return;
	255	}
	256	TClonesArray* digits = new TClonesArray("AliFMDDigit");
	257	digitBranch->SetAddress(&digits);
	258
	259	TIter next(&fAlgorithms);
	260	AliFMDMultAlgorithm* algorithm = 0;
	261	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	262	algorithm->PreEvent(clusterTree, fCurrentVertex);
	263	digitBranch->GetEntry(0);
	264
	265	ProcessDigits(digits);
	266
	267	next.Reset();
	268	algorithm = 0;
	269	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	270	algorithm->PostEvent();
	271	clusterTree->Fill();
	272	digits->Delete();
	273	delete digits;
	274	}
	275
	276	//____________________________________________________________________
	277	void
	278	AliFMDReconstructor::ProcessDigits(TClonesArray* digits) const
	279	{
	280	// For each digit, find the pseudo rapdity, azimuthal angle, and
	281	// number of corrected ADC counts, and pass it on to the algorithms
	282	// used.
	283	Int_t nDigits = digits->GetEntries();
	284	AliDebug(1, Form("Got %d digits", nDigits));
	285	for (Int_t i = 0; i < nDigits; i++) {
	286	AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->At(i));
	287	AliFMDGeometry* fmd = AliFMDGeometry::Instance();
	288	AliFMDDetector* subDetector = fmd->GetDetector(digit->Detector());
	289	if (!subDetector) {
	290	Warning("ProcessDigits", "Unknown detector: FMD%d" , digit->Detector());
	291	continue;
	292	}
	293
	294	AliFMDRing* ring = subDetector->GetRing(digit->Ring());
	295	Float_t ringZ = subDetector->GetRingZ(digit->Ring());
	296	if (!ring) {
	297	Warning("ProcessDigits", "Unknown ring: FMD%d%c", digit->Detector(),
	298	digit->Ring());
	299	break;
	300	}
	301
	302	Float_t realZ = fCurrentVertex + ringZ;
	303	Float_t stripR = ((ring->GetHighR() - ring->GetLowR())
	304	/ ring->GetNStrips() * (digit->Strip() + .5)
	305	+ ring->GetLowR());
	306	Float_t theta = TMath::ATan2(stripR, realZ);
	307	Float_t phi = (2 * TMath::Pi() / ring->GetNSectors()
	308	* (digit->Sector() + .5));
	309	Float_t eta = -TMath::Log(TMath::Tan(theta / 2));
	310	UShort_t counts = SubtractPedestal(digit);
	311
	312	TIter next(&fAlgorithms);
	313	AliFMDMultAlgorithm* algorithm = 0;
	314	while ((algorithm = static_cast<AliFMDMultAlgorithm*>(next())))
	315	algorithm->ProcessDigit(digit, eta, phi, counts);
	316	}
	317	}
	318
	319	//____________________________________________________________________
	320	UShort_t
	321	AliFMDReconstructor::SubtractPedestal(AliFMDDigit* digit) const
	322	{
	323	// Member function to subtract the pedestal from a digit
	324	// This implementation does nothing, but a derived class could over
	325	// load this to subtract a pedestal that was given in a database or
	326	// something like that.
	327
	328	Int_t counts = 0;
	329	Float_t ped = fPedestal + fPedestalFactor * fPedestalWidth;
	330	if (digit->Count3() > 0) counts = digit->Count3();
	331	else if (digit->Count2() > 0) counts = digit->Count2();
	332	else counts = digit->Count1();
	333	counts = TMath::Max(Int_t(counts - ped), 0);
	334	return UShort_t(counts);
	335	}
	336
	337	//____________________________________________________________________
	338	void
	339	AliFMDReconstructor::FillESD(TTree* /* digitsTree */,
	340	TTree* /* clusterTree */,
	341	AliESD* /* esd*/) const
	342	{
	343	// nothing to be done
	344	// FIXME: The vertex may not be known when Reconstruct is executed,
	345	// so we may have to move some of that member function here.
	346	#if 0
	347	TClonesArray* multStrips = 0;
	348	TClonesArray* multRegions = 0;
	349	TTree* treeR = fmdLoader->TreeR();
	350	TBranch* branchRegions = treeR->GetBranch("FMDPoisson");
	351	TBranch* branchStrips = treeR->GetBranch("FMDNaiive");
	352	branchRegions->SetAddress(&multRegions);
	353	branchStrips->SetAddress(&multStrips);
	354
	355	Int_t total = 0;
	356	Int_t nEntries = clusterTree->GetEntries();
	357	for (Int_t entry = 0; entry < nEntries; entry++) {
	358	AliDebug(5, Form("Entry # %d in cluster tree", entry));
	359	treeR->GetEntry(entry);
	360
	361
	362	Int_t nMults = multRegions->GetLast();
	363	for (Int_t i = 0; i <= nMults; i++) {
	364	AliFMDMultRegion* multR =
	365	static_cast<AliFMDMultRegion*>(multRegions->UncheckedAt(i));
	366	Int_t nParticles=multR->Particles();
	367	if (i>=0 && i<=13) hEtaPoissonI1->AddBinContent(i+1,nParticles);
	368	if (i>=14 && i<=27 ) hEtaPoissonI2->AddBinContent(i-13,nParticles);
	369	if (i>=28 && i<=33 );
	370	if (i>=34 && i<=47 ) hEtaPoissonI3->AddBinContent(48-i,nParticles);
	371	if (i>=48 && i<=53) hEtaPoissonO3->AddBinContent(54-i,nParticles);
	372	}
	373	}
	374	#endif
	375	}
	376
	377
	378	//____________________________________________________________________
	379	void
	380	AliFMDReconstructor::Reconstruct(AliRawReader,TTree) const
	381	{
	382	// Cannot be used. See member function with same name but with 2
	383	// TTree arguments. Make sure you do local reconstrucion
	384	AliError("MayNotUse");
	385	}
	386	//____________________________________________________________________
	387	void
	388	AliFMDReconstructor::Reconstruct(AliRunLoader*) const
	389	{
	390	// Cannot be used. See member function with same name but with 2
	391	// TTree arguments. Make sure you do local reconstrucion
	392	AliError("MayNotUse");
	393	}
	394	//____________________________________________________________________
	395	void
	396	AliFMDReconstructor::Reconstruct(AliRunLoader, AliRawReader) const
	397	{
	398	// Cannot be used. See member function with same name but with 2
	399	// TTree arguments. Make sure you do local reconstrucion
	400	AliError("MayNotUse");
	401	}
	402	//____________________________________________________________________
	403	void
	404	AliFMDReconstructor::FillESD(AliRawReader,TTree,AliESD*) const
	405	{
	406	// Cannot be used. See member function with same name but with 2
	407	// TTree arguments. Make sure you do local reconstrucion
	408	AliError("MayNotUse");
	409	}
	410	//____________________________________________________________________
	411	void
	412	AliFMDReconstructor::FillESD(AliRunLoader,AliESD) const
	413	{
	414	// Cannot be used. See member function with same name but with 2
	415	// TTree arguments. Make sure you do local reconstrucion
	416	AliError("MayNotUse");
	417	}
	418	//____________________________________________________________________
	419	void
	420	AliFMDReconstructor::FillESD(AliRunLoader,AliRawReader,AliESD*) const
	421	{
	422	// Cannot be used. See member function with same name but with 2
	423	// TTree arguments. Make sure you do local reconstrucion
	424	AliError("MayNotUse");
	425	}
	426
	427	//____________________________________________________________________
	428	//
	429	// EOF
	430	//