PHOS Reconstruction

Definitions and objectives

What is PHOS ?

PHOS is the association of two sub-detectors: Both sub-detectors are physically organized in modules. A separate document is dedicated to the description of the PHOS geometry.

Digitization

After the simulation, the file galice.root contains the hits stored in the tree TreeH. Then, the hits should be digitized, i.e. for each detector element (EMCA crystal, PPSD gas cell or CPV pad) the digitized amplitude should be obtained. The digitization is performed by the function AliPHOSv1::Hit2Digit(). The obtained digits are of the class AliPHOSDigit which is derives from AliDigit. It stores the digits as the pair (absolute Id, energy) in the array AliPHOSv1::fDigits which is written to the root file in the branch PHOS of the digit tree TreeD.

Reconstruction design

The reconstruction software consists in making particles from the digits generated after the simulation.

Use Case of reconstruction

Figure 1.: Use Case for reconstruction. Only PHOS is implemented. Other detectors are ignored. Click on image for full scale

The PHOS reconstruction is delegated by the detector to a reconstructioner algorithmic object. It proceeds along three steps:
Clusterization
The reconstruction of the showers (total energy loss and incident direction) induced in the EM calorimeter (AliPHOSEmcRecPoint) and of the pads hit in the pre-shower detector (AliPHOSPpsdRecPoint). The clustering is delegated to the algorithmic object AliPHOSClusterizerv1 which derives from the interface AliPHOSClusterizer. The clusterization is controlled by several parameters.
Sub-Track construction
Sub tracks are obtained by combining reconstructed points in EMCA and PPSD to form a track segment (AliPHOSTrackSegment). This construction is delegated to the algorithmic object AliPHOSTrackSegmentMakerv1 which derives from the interface AliPHOSTrackSegmentMaker.
Particle identification
From each track segment a reconstructed particle (AliPHOSRecParticle) is made the identification being based on the information delivered by the PPSD (which PPSD stage is hit) and by EMCA (shower profile and topology).The identification is delegated to the algorithmic object AliPHOSPIDv1 which derives from the interface AliPHOSPID. Several parameters control the identification procedure.
The first step is detector specific, i.e. does not require any information about other detectors. The two last steps may (and probably will) depend on other detector reconstructions (connecting sub tracks or combining particle identification delivered by various sub-detectors).

How to use it ? (CVS server /PHOS)


0. Starting point. I assume that the reader is able to performe an aliroot simulation of the PHOS detector. This implies that a root file has been created containing the information of digits (in the digits tree, fTreeD of AliRun) and, depending on the user choice, the information about hits (in the hit tree, fTreeH of AliRun). At this level of the simulation (version header of Aliroot in Mai 00), the reconstructed tree, fTreeR is empty.

1. Description of the reconstruction As it is explained if Fig.1, one will need to do the reconstruction proccess. The reconstruction proccess is done calling to the AliPHOSv1::Reconstruction() function. However, before call it, the user must define the reconstruction procedures, which is managed by the class AliPHOSReconstructioner. The latest contains the three reconstruction levels :
  1. Clusterizer procedure (AliPHOSClusterizerv1)
  2. TrackSegment procedure (AliPHOSTrackSegmentMakerv1)
  3. Particle identification (AliPHOSPIDv1)

Durign each stage of the reconstruction, reconstruction objects are stored in the reconstruction tree, fTreeR:
  1. List of emc and ppsd reconstructed points (AliPHOSEmcRecPoint and AliPHOSPpsdRecPoint) in the branches PHOSEmcRP and PHOSPpsdRP
  2. List of TrackSegments (AliPHOSTrackSegment) in the branch PHOSTS
  3. List of Particles (AliPHOSParticle) in the branch PHOSRP

2. Doing Reconstruction The reconstruction procedure has been compacted by the AliPHOSAnalyze class (see in the (CVS server) ). The following macro is enough to performed the reconstruction and fill the reconstruction tree in the root file. 
 Int_t ievent;
 Int_t nevents = 10; 
 AliPHOSAnalyze reconstruction("galice.root");
  for(ievent=0; i nevents; ievent++) {    
    reconstruction.AnalyzeOneEvent(ievent);    
  }

An the end of the stage, the galice root file will contain the information of the reconstruction in the treeR:
AliPHOSEmcRecPoint and AliPHOSPpsdRecPoint -- AliPHOSTrackSegment -- AliPHOSParticle

3. Analysing the Reconstructed tree The following code shows how one can analyse one event of the reconstructed tree. I this example, we get the primary particle in the Kine tree at the origin of the detected AliPHOSParticle. This example show the that full Fererico's parabola is connected and one can go back uption the initial point.
{
// Opening root file with reconstructed events
TFile * in = new TFile("galice.root");
// Getting gAlice 
gAlice = (AliRun*) in->Get("gAlice");
// Getting PHOS Detector, Geometry and Indexig    
AliPHOSv1            * fPHOS      = (AliPHOSv1 *)gAlice->GetDetector("PHOS")         ;
AliPHOSGeometry      * fGeom      = AliPHOSGeometry::GetInstance( fPHOS->GetGeometry()->GetName(), fPHOS->GetGeometry()->GetTitle() );
AliPHOSIndexToObject * fObjGetter = AliPHOSIndexToObject::GetInstance(fPHOS) ; 

// Clones Array of tracksegments and particles
TObjArray    * EmcRPList        = new  TObjArray(1000);
TObjArray    * PpsdRPList       = new  TObjArray(1000);
TClonesArray * ParticleList     = new  TClonesArray("AliPHOSRecParticle", 1000) ;
TClonesArray * TrackSegmentList = new  TClonesArray("AliPHOSTrackSegment", 1000) ;
TClonesArray * ArrayOfParticles = gAlice->Particles();    
     

// Getting Event
Int_t ievent = 8;
Int_t nprimaryparticles = gAlice->GetEvent(ievent); // Getting the physics event

cout << "FinalAnalysis >>> Event #" << ievent << " with " << nprimaryparticles << " primary particles." << endl;


// Getting Reconstrcuted Tree and Brach with TrackSegments and Reconstructed Particles
TTree     * Rec       = gAlice->TreeR();
TTree     * TK        = gAlice->TreeK();
TTree     * TD        = gAlice->TreeD(); 
TBranch   * RecBranch = Rec->GetBranch("PHOSRP");
TBranch   * TraBranch = Rec->GetBranch("PHOSTS");
TBranch   * EmcBranch = Rec->GetBranch("PHOSEmcRP");
TBranch   * PpsdBranch = Rec->GetBranch("PHOSPpsdRP");
EmcBranch-> SetAddress(&EmcRPList);
PpsdBranch->SetAddress(&PpsdRPList);
TraBranch-> SetAddress(&TrackSegmentList);
RecBranch-> SetAddress(&ParticleList);

TD->GetEvent(0);  //This tree only contains one tree-event with the List of Digits of the physicis event.
TK->GetEvent(0);  // This tree only contains one tree-event with the List of TParticles of the physicis event.
Int_t nbits = Rec->GetEvent(0); // This tree only contains one tree-event with the List of Reconstructed Objects (AliPHOSEmcRecPoint, AliPHOSPpsdRecPoint, AliPHOSTrackSegment and AliPHOSParticle)
cout << "FinalAnalysis >>> " << nbits << " bits in reconstructed tree envent 0." << endl;

Int_t nparticles =  ParticleList.GetEntries();
Int_t iparticle;
cout << "FinalAnalysis>>> Reconstructed particle list entries is "  << nparticles <At(iparticle);
  AliPHOSTrackSegment * ts = TrackSegmentList->At(rp->GetPHOSTrackSegmentIndex()); // TrackSegmentList->At(iparticle);
  AliPHOSEmcRecPoint  * emc = EmcRPList->At(ts->GetEmcRecPointIndex());
  TVector3 locpos; 
  emc->GetLocalPosition(locpos) ;
  Int_t * primaries; 
  Int_t nprimaries;
  Text_t particle[11];
  primaries = (emc->GetPrimaries(nprimaries));
  switch(rp->GetType())
    {
    case 0:
      strcpy( particle, "NEUTRAL_EM");
      break;
    case 1:
      strcpy(particle, "NEUTRAL_HA");
      break;
    case 2:
      strcpy(particle, "GAMMA");
      break ;
    case 3: 
      strcpy(particle, "GAMMA_H");
      break ;
    case 4:
      strcpy(particle, "ABSURD_EM") ;
      break ;
    case 5:
      strcpy(particle, "ABSURD_HA") ;
      break ;	
    case 6:
      strcpy(particle, "ELECTRON") ;
      break ;
    case 7:
      strcpy(particle, "CHARGED_HA") ;
      break ; 
     }
      
  cout << "FinalAnalysis>>> Partilce is " << particle <<  " and energy is " << rp->Energy() << " GeV" << endl ;
  cout << "FinalAnalysis>>> Index in list is " <<  rp->GetIndexInList() << " "  <<   endl ;
  cout << "FinalAnalysis>>> Position is: X="<>> Primaries are : " << primaries[0] << " " << primaries[1] << " " << primaries[2] << " " << primaries[3] << " " << endl;

  Int_t nprimaryparticle = ArrayOfParticles->GetEntriesFast();
  printf("FinalAnalysis >>> nprimaryparticle is %d \n",nprimaryparticle);
  
  Particle = (TParticle*)ArrayOfParticles->UncheckedAt(primaries[0]);
  Int_t   iPDGCode = Particle->GetPdgCode();
  cout << "FinalAnalysis >>> PDGCode of primary is " << iPDGCode  << endl;
}
}

4. Examples of analysis
The class AliPHOSAnalyze contains several methods which performes physics analysis. The following code shows how one can analyse one event of the reconstructed tree.

Clusterization

A digit is a set of two numbers: the id of the elementary cell numbered from 1 to the maximum number of cells (EMCA+PPSD), and a digitized deposited energy. With the help of the geometry object the absolute id is converted into a relative numbering with a single EMCA or PPSD module, in terms of rows and columns ( ALIPHOSGeometry::AbsToRelNumbering() ). A clustering algorithm ( AliPHOSClusterizerv1::MakeClusters() ) groups neighbouring cells in EMCA and PPSD ( AliPHOSClusterizerv1::AreNeighbours() ) where two cells are defined as neighbours if they have a common edge or a common corner. The algorithm is controlled with alltogether 4 parameters: An ensemble of neighbouring digits is called a reconstructed point. The PHOS clusterizer produces EmcRecPoints and PpsdRecPoints deriving form the base class AliPHOSRecPoint. Each reconstructed point points, through a mechanism described later, to the digits at its origin. Since a digit remembers which primary particle(s) were at its origin, the list of primary particles at the origin of a reconstructed point is provided by the method AliPHOSRecPoint::GetPrimaries()).

The reconstructed points are contained in two TObjArray (mandatory because the size of a reconstructed point depends on the number of digits used to build it) which are stored in TreeR, branch PHOSEmcRP for EMCA and branch PHOSPpsdRP for PPSD.

Track Segments Maker

Each EMCA reconstructed point is either combined with one or two PPSD reconstructed points to form a track segment or constitutes a track segment by itself. As a first step the EMCA reconstructed points are unfolded in case there is more than one local maximum ( overlapping showers) to make two new reconstructed points which digit's energy is shared after a fitting procedure. The unfolding procedure can be switch on and off with the method: AliPHOSTrackSegmentMakerv1::(Un)SetUnfoldFlag(). EMCA reconstructed points and PPSD reconstructed points are then linked together to form pairs or triplets according to their relative distance. Each track segment points is linked, through a mechanism described later, to the reconstructed points at its origin.

The track segments are contained in a TClonesArray which are stored in TreeR, branch PHOSTS.

Particle Identification

Clusterization and track segments making are two activities that are detector-driven. This is not the case for the particle identification. Indeed a particle might be identified by combining the track segments calculated by several detectors. Up to know the implementation for such a mechanism is not yet done. Instead PHOS has its own stand alone particle identification algorithm which associates a reconstructed particle to each track segment. The particle identification is controlled by several parameters which can be set by the user. Each reconstructed particle is linked, through a mechanism described later, to the track segment at its origin. The various types of reconstructed particles are the following:
NEUTRAL_EM: a reconstructed point in EMCA, none in the two layers of the PPSD, and narrow lateral development;
NEUTRAL_HA : a reconstructed point in EMCA, none in the two layers of the PPSD, and wide lateral development;
GAMMA : a reconstructed point in EMCA, none in the first layer CPV, one in the second layer PC, narrow lateral development;
GAMMA_H: a reconstructed point in EMCA, none in the first layer CPV, one in the second layer PC, wide lateral development;
ABSURD_EM: a reconstructed point in EMCA, one in the first layer CPV, none in the second layer PC, narrow lateral development;
ABSURD_HA: a reconstructed point in EMCA, one in the first layer CPV, none in the second layer PC, wide lateral development;
ELECTRON: a reconstructed point in EMCA, one in both layers of the PPSD, and narrow lateral development;
CHARGED_HA : a reconstructed point in EMCA, one in both layers of the PPSD, and wide lateral development;

Disk Storage Structure

All the three types of reconstructed objects( reconstructed point, track segment and reconstructed particle) are stored in the reconstruction Tree (TreeR). The container of the reconstructed points is a TObjArray and the container of the track segment and the reconstructed particles is a TClonesArray. The names of the branches are respectively: PHOSEmcRP, PHOSPpsdRP, PHOSTS and PHOSRP. The TreeR is best seen in the following example: 
      root [2] TFile myfile("junk.root")
      root [3] TTree * reconstructedtree = (TTree *)myfile.Get("TreeR0")
      root [5] reconstructedtree.Print()
      ******************************************************************************
      *Tree    :TreeR0    : Reconstruction                                         *
      *Entries :        1 : Total  Size =     10882 bytes  File  Size =       4024 *
      *        :          : Tree compression factor =   5.65                       *
      ******************************************************************************
      *Branch  :PHOSEmcRP : PHOSEmcRP                                              *
      *Entries :        1 : Total  Size =         0 bytes  File Size  =          0 *
      *Baskets :        0 : Basket Size =     16000 bytes  Compression=   1.00     *
      *............................................................................*
      *Branch  :PHOSPpsdRP : PHOSPpsdRP                                            *
      *Entries :        1 : Total  Size =      8334 bytes  File Size  =       1476 *
      *Baskets :        1 : Basket Size =     16000 bytes  Compression=   5.65     *
      *............................................................................*
      *Branch  :PHOSTS    : PHOSTS                                                 *
      *Entries :        1 : Total  Size =         0 bytes  File Size  =          0 *
      *Baskets :        0 : Basket Size =     16000 bytes  Compression=   1.00     *
      *............................................................................*
      *Branch  :PHOSRP    : PHOSRP                                                 *
      *Entries :        1 : Total  Size =         0 bytes  File Size  =          0 *
      *Entries :        1 : Total  Size =         0 bytes  File Size  =          0 *
      *Baskets :        0 : Basket Size =     16000 bytes  Compression=   1.00     *
      *............................................................................*
      root [6]
branches in reconstructed Tree

Figure 2: Branches in TreeR generated by the PHOS reconstruction. PHOSEmc(Ppsd)RP contains EMCA(PPSD) reconstructed points in TObjArray; PHOSTS contains track segments in a TClonesArray; PHOSRP contains reconstructed particles in a TClonesArray.

How are the links implemented ? ( The information bellow must be updated, so forget it for the moment...)

An algorithmic class, AliPHOSIndexToObject returns a pointer to the object of interest given the index of its storage in the array (TObjArray or TClonesArray). It is a singleton (only one instance of the object does exist at run time). It must be initialized once: 
      AliPHOS * phos = (AliPHOSv1 *)gAlice->GetDetector("PHOS") ;
      AliPHOSIndexToObject::GetInstance(phos) ;
It can then be used as follow, for example for retrieving a particle from the kine tree (TreeK): 
     AliPHOSIndexToObject * please = AliPHOSIndexToObject::GetInstance() ;
     Int_t index = 123 ; 
     TParticle * primaryparticle = please->GimePrimaryParticle(index) ;
We shall now detail how the various reconstruction objects are linked together.

Digits to primary particles

Digits have three additional data members (fPrimary1(2,3)) that gives the index of the primary particles (maximum 3) that have contributed to the formation of the digit. Remember that in PHOS a digit is obtained from all the hits accumulated in a single cell. This is not really elegant, but because the container of digits is a TClonesArray it is not possible to replace the three data members by an array of integers of variable length as all objects in a TClonesArray must be identical. To retrieve the particles that have contibuted to a digit, do: 
    // initialization
          // open root file
     TFile rootfile("junk.root") ;               
          // get AliRun object
     gAlice = (AliRun *)rootfile.Get("gAlice") ;
          // get detector object
     AliPHOSv1 * phos = (AliPHOSv1 *)gAlice->GetDetector("PHOS") ;
          // get the geometry associated with the detector
     AliPHOSGeometry::GetInstance( phos->GetGeometry()->GetName(), phos->GetGeometry()->GetTitle() ) ; 
          // initializes the index to object converter
     AliPHOSIndexToObject::GetInstance(phos) ;
          // get the list of digits
     TClonesArray * digitslist = phos->Digits() ; 
          // loop over the list of digits
     TIter nextdigit(digitslist) ; 
     AliPHOSDigit * digit ; 
     TParticle * primaryparticle[3] ; 
          // get the pointer of the index to object converter
     AliPHOSIndexToObject * please = AliPHOSIndexToObject::GetInstance() ;
     while ( digit = (AliPHOSDigit * )nextdigit() ) {
       for ( Int_t i = 0 ; i < digit->GetNprimary() ; i++ ) 
       primary[i] = please->GimePrimaryParticle(digit->GetPrimary(i) ) ;
     }

Reconstructed point to Digit

A digit has an additional data member which is the index of the object stored in the TObjArray. It is set by FinishEvent(). A reconstructed point containing a list of digits, it can retrieve the digit object with the help of AliPHOSIndexToObject: 
      // initialization
          // open root file
     TFile rootfile("junk.root") ;               
          // get AliRun object
     gAlice = (AliRun *)rootfile.Get("gAlice") ;
          // get detector object
     AliPHOSv1 * phos = (AliPHOSv1 *)gAlice->GetDetector("PHOS") ;
          // get the geometry associated with the detector
     AliPHOSGeometry::GetInstance( phos->GetGeometry()->GetName(), phos->GetGeometry()->GetTitle() ) ; 
          // initializes the index to object converter
     AliPHOSIndexToObject::GetInstance(phos) ;
          // get the reconstructed points list
     Int_t evt = 123 ; 
     TObjArray * recpointslist = phos->EmcRecPoints(evt) ; 
          // loop over reconstructed points 
     TIter nextrecpoint(recpointslist) ; 
          // get the pointer of the index to object converter
     AliPHOSIndexToObject * please = AliPHOSIndexToObject::GetInstance() ;
     AliPHOSEmcRecPoint * recpoint ;  
     while ( recpoint = (AliPHOSEmcRecPoint * )nextrecpoint() ) {
             // get the associated digits list
       Int_t * digitsindexeslist = recpoint->GetDigitsList() ;
             // loop over the list of digits
       for ( Int_t i = 0 ; i < recpoint->GetDigitsMultiplicity() ; i++ ) {
         AliPHOSDigit * digit = please->GimeDigit(digitsindexeslist[i] ) ;
               // get the primary particle associated with that digit
 	 Int_t numberofprimaries = 0 ; 
	 Int_t * prim = digit->GetPrimaries(numberofprimaries) ;
         for (Int_t i = 0 ; i < numberofprimaries ; i++ )
           // and print them 
	   please->GimePrimaryParticle( prim[i] )->Print()  ;
         } 
       } 
         // or get the primaries directly from the reconstructed point
       Int_t numberofprimaries = 0 ; 
       Int_t * prim = recpoint->GetPrimaries(numberofprimaries) ;  
     }

Track segment to reconstructed point

A reconstructed point has an additional data member which is the index of the object stored in the TClnesArray. It is set by the reconstructioner. A track segment contains the index of one EMCA reconstructed points and two PPSD reconstructed points, it can retrieve this reconstructed point objects with the help of AliPHOSIndexToObject: 
     // initialization
          // open root file
     TFile rootfile("junk.root") ;               
          // get AliRun object
     gAlice = (AliRun *)rootfile.Get("gAlice") ;
          // get detector object
     AliPHOSv1 * phos = (AliPHOSv1 *)gAlice->GetDetector("PHOS") ;
          // get the geometry associated with the detector
     AliPHOSGeometry::GetInstance( phos->GetGeometry()->GetName(), phos->GetGeometry()->GetTitle() ) ; 
          // initializes the index to object converter
     AliPHOSIndexToObject::GetInstance(phos) ;
          // get the track segments list
     Int_t evt = 123 ; 
     TClonesArray * tracksegmentslist = phos->TrackSegments(evt) ; 
          // loop over track segments
     TIter nexttracksegment(tracksegmentslist) ; 
     AliPHOSTrackSegment * tracksegment ;  
         // get the pointer of the index to object converter
     AliPHOSIndexToObject * please = AliPHOSIndexToObject::GetInstance() ;
     while ( tracksegment = (AliPHOSTrackSegment * )nexttracksegment() ) {
           // get the associated reconstructed points
         AliPHOSEmcRecPoint  * emcrecpoint     = please->GimeRecPoint(tracksegment->GetEmcRecPoint(),    "emc" ) ;
         AliPHOSPpsdRecPoint * ppsduprecpoint  = please->GimeRecPoint(tracksegment->GetPpsdUpRecPoint(), "ppsd" ) ;
         AliPHOSPpsdRecPoint * ppsdlowrecpoint = please->GimeRecPoint(tracksegment->GetPpsdLowRecPoint(),"ppsd" ) ;
          // get the primaries particles  
         Int_t numberofprimariestoemc = 0 ; 
         Int_t * primemc = tracksegment->GetPrimariesEmc(numberofprimariestoemc) ;
         Int_t numberofprimariestoppsdup = 0 ; 
         Int_t * primppsdup = tracksegment->GetPrimariesPpsdUp(numberofprimariestoppsdup) ;
         Int_t numberofprimariestoppsdlow = 0 ; 
         Int_t * primppsdlow = tracksegment->GetPrimariesPpsdLow(numberofprimariestoppsdlow) ;
          // print one as example  
         please->GimePrimaryParticle( primppsdlow[0] )->Print()  ;     
       }
     rootfile.Close() ;
     gAlice = 0 ; 
     phos = 0 ; 
     recparticleslist = 0 ;

Reconstructed particle to Track segment

A track segment has an additional data member which is the index of the object stored in the TClnesArray. It is set by the reconstructioner. A reconstructed particle contains the index of a track segment, it can retrieve this reconstructed point objects with the help of AliPHOSIndexToObject: 
     // initialization
          // open root file
     TFile rootfile("junk.root") ;               
          // get AliRun object
     gAlice = (AliRun *)rootfile.Get("gAlice") ;
          // get detector object
     AliPHOSv1 * phos = (AliPHOSv1 *)gAlice->GetDetector("PHOS") ;
          // get the geometry associated with the detector
     AliPHOSGeometry::GetInstance( phos->GetGeometry()->GetName(), phos->GetGeometry()->GetTitle() ) ; 
          // initializes the index to object converter
     AliPHOSIndexToObject::GetInstance(phos) ;
          // get the reconstructed particles list
     Int_t evt = 123 ; 
     TClonesArray * recparticleslist = phos->RecParticles(evt) ; 
     // loop over reconstructed particles
     TIter nextrecparticle(recparticleslist) ; 
     AliPHOSRecParticle * recparticle ;  
       // get the pointer of the index to object converter
     AliPHOSIndexToObject * please = AliPHOSIndexToObject::GetInstance() ;
     while ( recparticle = (AliPHOSRecParticle * )nextrecparticle() ) {
           // get the track segment ...
         AliPHOSTrackSegment  * tracksegment  = recparticle->GetPHOSTrackSegment()  ;
           // and print it
	 tracksegment->Print() ;
           // get the list of primaries ...
	 Int_t numberofprimaries = 0 ; 
	 Int_t * prim = recparticle->GetPrimaries(numberofprimaries) ;
         for (Int_t i = 0 ; i < numberofprimaries ; i++ )
           // and print them 
	   please->GimePrimaryParticle( prim[i] )->Print()  ; 
    }     
     rootfile.Close() ;
     gAlice = 0 ; 
     phos = 0 ; 
     recparticleslist = 0 ;
A short cut allows to access the primaries at the origin of the reconstructed particle: 
     Int_t numberofprimaries = 0 ; 
     Int_t * prim = recparticle->GetPrimariesPpsdLow(numberofprimaries) ;
© Groupe Photons Subatech [Go to the GPS Home Page]
Last modified: Mon Jul 31 11:36:07 CEST 2000