[u/mrichter/AliRoot.git] / PWGCF / Correlations / JCORRAN / AliJPhoton.cxx

// $Id: AliJPhoton.cxx,v 1.4 2008/05/08 13:44:45 djkim Exp $

////////////////////////////////////////////////////
/*!
  \file AliJPhoton.cxx
  \brief
  \author J. Rak, D.J.Kim, R.Diaz (University of Jyvaskyla)
  \email: djkim@jyu.fi
  \version $Revision: 1.4 $
  \date $Date: 2008/05/08 13:44:45 $
  */
////////////////////////////////////////////////////

#include "AliJBaseTrack.h"
#include "AliJPhoton.h"
#include "TF1.h"

ClassImp(AliJPhoton);
//______________________________________________________________________________
AliJPhoton::AliJPhoton() : 
    AliJBaseTrack(),
    fChi2(-999),              
    fTof(-999),                   
    fX(-999),               
    fY(-999),                  
    fZ(-999),
    fProbPhot(-999), 
    fTrackDx(-1), 
    fTrackDz(-1), 
  fCaloType(kUndef),
  fDistToBadChannel(-999),
  fDispersion(-999),
  fM20(-999),
  fM02(-999),
  fEmcCpvDist(-999),
  fNCells(0),
  fSuperModuleId(-999),
  fCellsAbsId(0x0),
  fCellsAmpFraction(0x0)

{
  // default constructor
  for(Int_t i=0;i<kUnknownAli+1;i++) fCaloPID[i]=-1;

  SetPID((Double_t*)NULL);
}

//_____________________________________________________________________________
AliJPhoton::AliJPhoton(const AliJPhoton& a) : 
    AliJBaseTrack(a),
    fChi2(a.fChi2),
    fTof(a.fTof),
    fX(a.fX),
    fY(a.fY),
    fZ(a.fZ),
    fProbPhot(a.fProbPhot), 
    fTrackDx(a.fTrackDx), 
    fTrackDz(a.fTrackDz), 
  fCaloType(a.fCaloType),
  fDistToBadChannel(a.fDistToBadChannel),
  fDispersion(a.fDispersion),
  fM20(a.fM20),
  fM02(a.fM02),
  fEmcCpvDist(a.fEmcCpvDist),
  fNCells(a.fNCells),
  fSuperModuleId(a.fSuperModuleId),
  fCellsAbsId(NULL),
  fCellsAmpFraction(NULL)

{
  //copy constructor
  for(Int_t i=0;i<kUnknownAli+1;i++) fCaloPID[i] = a.fCaloPID[i];
  SetCellsAbsId( a.fCellsAbsId );
  SetCellsAmplitudeFraction( a.fCellsAmpFraction );
}

//_____________________________________________________________________________
AliJPhoton& AliJPhoton::operator=(const AliJPhoton& photon){
  //operator=    
  if(this != &photon){
    AliJBaseTrack::operator=(photon);
    fChi2 = photon.fChi2;
    fTof  = photon.fTof;
    fX    = photon.fX;
    fY    = photon.fY;
    fZ    = photon.fZ;
    fProbPhot = photon.fProbPhot;
    fTrackDx  = photon.fTrackDx;
    fTrackDz  = photon.fTrackDz;
    fCaloType = photon.fCaloType;
    for(Int_t i=0; i<kUnknownAli+1; i++){
      fCaloPID[i] = photon.fCaloPID[i];
    }
    fDistToBadChannel = photon.fDistToBadChannel;
    fDispersion    = photon.fDispersion;
    fM20           = photon.fM20;
    fM02           = photon.fM02;
    fEmcCpvDist    = photon.fEmcCpvDist;
    fNCells        = photon.fNCells;
    fSuperModuleId =  photon.fSuperModuleId;
    SetCellsAbsId( photon.fCellsAbsId );
    SetCellsAmplitudeFraction( photon.fCellsAmpFraction );

  }

  return *this;
}
//______________________________________________________________________________
void  AliJPhoton::SetCellsAbsId(const UShort_t *array)
{
    //  Set the array of cell absId numbers 
    if (fNCells) {
        if(fCellsAbsId){ delete [] fCellsAbsId; fCellsAbsId = NULL; }
	fCellsAbsId = new  UShort_t[fNCells];
	for (Int_t i = 0; i < fNCells; i++) fCellsAbsId[i] = array[i];
    }
}

//_______________________________________________________________________
void  AliJPhoton::SetCellsAmplitudeFraction(const Double32_t *array)
{
    //  Set the array of cell amplitude fraction
    if (fNCells) {
	if(fCellsAmpFraction){ delete [] fCellsAmpFraction; fCellsAmpFraction = NULL;}
	fCellsAmpFraction = new  Double32_t[fNCells];
	for (Int_t i = 0; i < fNCells; i++) fCellsAmpFraction[i] = array[i];
    }
}

//______________________________________________________________________________
void  AliJPhoton::SetPID(const Double_t *pid) {
   //set pid
    if(pid){
       for(Int_t i=0; i<kUnknownAli+1; ++i) fCaloPID[i]=pid[i];
       SetProbPhot(fCaloPID[kPhotonAli]);
    }else{
      for(Int_t i=0; i<kUnknownAli+1; fCaloPID[i++]=0.){} 
      fCaloPID[kUnknownAli]=1.;
    }
}


//______________________________________________________________________________
particleType AliJPhoton::GetParticleType() {

    // return the most problable particle type
    // Note: following the AliCaloPID implementation 
    
    //Init default weights 
    Float_t wPhoton = 0.75 ;
    Float_t wPi0 = 0.8 ;
    Float_t wElectron = 0.5 ;
    Float_t wCharged = 0.5 ;
    Float_t wNeutral = 0.5 ;

    Bool_t usePHOSweightFormula = kTRUE;
    //Formula to set the PID weight threshold for photon or pizero
    // ALICE-INT-2005-016 (2007)
    TFormula* wPhotonPHOSFormula = 
    new TFormula("photonWeight","0.75*(x<40)+ 0.68*(x>=100)+(x>=40 && x<100)*(0.98+x*(6e-3)-x*x*(2e-04)+x*x*x*(1.1e-06))");
    TFormula* wPi0PHOSFormula = 
    new TFormula("pi0Weight","0.80*(x<65)+ 0.915*(x>=100)+(x>=65 && x-x*(1.95e-3)-x*x*(4.31e-05)+x*x*x*(3.61e-07))");

    if(fCaloType == kPHOSCalo && usePHOSweightFormula){
       wPhoton  = wPhotonPHOSFormula->Eval(E()) ;
       wPi0 = wPi0PHOSFormula->Eval(E());
      // cout << "wPhotonF " << wPhoton << " wphoton " << fCaloPID[kPhotonAli] << " wPi0F " << wPi0 << " wpi0 " << fCaloPID[kPi0Ali] << endl;
    }
    
    if(fCaloType == kEMCALCalo){
       wPhoton   =  0.8 ;
       wPi0      =  0.5 ;
       wElectron =  0.8 ;
       wCharged  =  0.5 ;
       wNeutral  =  0.5 ;
    }
    
  
   particleType pdg = kHadron;

  //Select most probable ID
  if(fCaloType == kPHOSCalo){
    if(fCaloPID[kPhotonAli] > wPhoton) pdg = kPhoton ;
    else if(fCaloPID[kPi0Ali] > wPi0) pdg = kPizero ; 
    //else if(fCaloPID[kElectronAli] > wElectron)  pdg = electron ;
    //else if(fCaloPID[kEleConAli] >  wElectron) pdg = electronCon ;
    //else if(chargedHadronWeight > wCharged) pdg = chargedHadron ;  
    //else if(neutralHadronWeight > wNeutral) pdg = neutralHadron ; 
    //else if(allChargedWeight >  allNeutralWeight)
    // pdg = chargedUnknown ; 
    //else 
    //  pdg = neutralUnknown ;
  }
  else{//EMCAL
    //Temporal solution, electrons and photons not differenciated
    if(fCaloPID[kPhotonAli] + fCaloPID[kElectronAli]  > wPhoton) pdg = kPhoton ;
    else if(fCaloPID[kPi0Ali] > wPi0) pdg = kPizero ; 
    //else if(chargedHadronWeight + neutralHadronWeight > wCharged) pdg = chargedHadron ;  
    //else if(neutralHadronWeight + chargedHadronWeight > wNeutral) pdg = neutralHadron ; 
    //else pdg =  neutralUnknown ;

  }
  
  delete wPhotonPHOSFormula;
  delete wPi0PHOSFormula;
  
  //if(fCaloType == kPHOSCalo) cout << "pdg " << pdg << endl;
  
  return pdg ;

}
Commit	Line	Data
37dde34e	1	// $Id: AliJPhoton.cxx,v 1.4 2008/05/08 13:44:45 djkim Exp $
	2
	3	////////////////////////////////////////////////////
	4	/*!
	5	\file AliJPhoton.cxx
	6	\brief
	7	\author J. Rak, D.J.Kim, R.Diaz (University of Jyvaskyla)
	8	\email: djkim@jyu.fi
	9	\version $Revision: 1.4 $
	10	\date $Date: 2008/05/08 13:44:45 $
	11	*/
	12	////////////////////////////////////////////////////
	13
	14	#include "AliJBaseTrack.h"
	15	#include "AliJPhoton.h"
	16	#include "TF1.h"
	17
	18	ClassImp(AliJPhoton);
	19	//______________________________________________________________________________
	20	AliJPhoton::AliJPhoton() :
	21	AliJBaseTrack(),
	22	fChi2(-999),
	23	fTof(-999),
	24	fX(-999),
	25	fY(-999),
	26	fZ(-999),
	27	fProbPhot(-999),
	28	fTrackDx(-1),
	29	fTrackDz(-1),
	30	fCaloType(kUndef),
	31	fDistToBadChannel(-999),
	32	fDispersion(-999),
	33	fM20(-999),
	34	fM02(-999),
	35	fEmcCpvDist(-999),
	36	fNCells(0),
	37	fSuperModuleId(-999),
	38	fCellsAbsId(0x0),
	39	fCellsAmpFraction(0x0)
	40
	41	{
	42	// default constructor
	43	for(Int_t i=0;i<kUnknownAli+1;i++) fCaloPID[i]=-1;
	44
	45	SetPID((Double_t*)NULL);
	46	}
	47
	48	//_____________________________________________________________________________
	49	AliJPhoton::AliJPhoton(const AliJPhoton& a) :
	50	AliJBaseTrack(a),
	51	fChi2(a.fChi2),
	52	fTof(a.fTof),
	53	fX(a.fX),
	54	fY(a.fY),
	55	fZ(a.fZ),
	56	fProbPhot(a.fProbPhot),
	57	fTrackDx(a.fTrackDx),
	58	fTrackDz(a.fTrackDz),
	59	fCaloType(a.fCaloType),
	60	fDistToBadChannel(a.fDistToBadChannel),
	61	fDispersion(a.fDispersion),
	62	fM20(a.fM20),
	63	fM02(a.fM02),
	64	fEmcCpvDist(a.fEmcCpvDist),
65	fNCells(a.fNCells),
66	fSuperModuleId(a.fSuperModuleId),
67	fCellsAbsId(NULL),
68	fCellsAmpFraction(NULL)
69
70	{
71	//copy constructor
72	for(Int_t i=0;i<kUnknownAli+1;i++) fCaloPID[i] = a.fCaloPID[i];
73	SetCellsAbsId( a.fCellsAbsId );
74	SetCellsAmplitudeFraction( a.fCellsAmpFraction );
75	}
76
77	//_____________________________________________________________________________
78	AliJPhoton& AliJPhoton::operator=(const AliJPhoton& photon){
79	//operator=
80	if(this != &photon){
81	AliJBaseTrack::operator=(photon);
82	fChi2 = photon.fChi2;
83	fTof = photon.fTof;
84	fX = photon.fX;
85	fY = photon.fY;
86	fZ = photon.fZ;
87	fProbPhot = photon.fProbPhot;
88	fTrackDx = photon.fTrackDx;
89	fTrackDz = photon.fTrackDz;
89e21aee	90	fCaloType = photon.fCaloType;
37dde34e	91	for(Int_t i=0; i<kUnknownAli+1; i++){
	92	fCaloPID[i] = photon.fCaloPID[i];
	93	}
	94	fDistToBadChannel = photon.fDistToBadChannel;
	95	fDispersion = photon.fDispersion;
	96	fM20 = photon.fM20;
	97	fM02 = photon.fM02;
	98	fEmcCpvDist = photon.fEmcCpvDist;
	99	fNCells = photon.fNCells;
	100	fSuperModuleId = photon.fSuperModuleId;
89e21aee	101	SetCellsAbsId( photon.fCellsAbsId );
89e21aee	102	SetCellsAmplitudeFraction( photon.fCellsAmpFraction );
37dde34e	103
	104	}
	105
	106	return *this;
	107	}
	108	//______________________________________________________________________________
	109	void AliJPhoton::SetCellsAbsId(const UShort_t *array)
	110	{
	111	// Set the array of cell absId numbers
	112	if (fNCells) {
	113	if(fCellsAbsId){ delete [] fCellsAbsId; fCellsAbsId = NULL; }
	114	fCellsAbsId = new UShort_t[fNCells];
	115	for (Int_t i = 0; i < fNCells; i++) fCellsAbsId[i] = array[i];
	116	}
	117	}
	118
	119	//_______________________________________________________________________
	120	void AliJPhoton::SetCellsAmplitudeFraction(const Double32_t *array)
	121	{
	122	// Set the array of cell amplitude fraction
	123	if (fNCells) {
	124	if(fCellsAmpFraction){ delete [] fCellsAmpFraction; fCellsAmpFraction = NULL;}
	125	fCellsAmpFraction = new Double32_t[fNCells];
	126	for (Int_t i = 0; i < fNCells; i++) fCellsAmpFraction[i] = array[i];
	127	}
	128	}
	129
	130	//______________________________________________________________________________
a10c1c2a	131	void AliJPhoton::SetPID(const Double_t *pid) {
37dde34e	132	//set pid
	133	if(pid){
	134	for(Int_t i=0; i<kUnknownAli+1; ++i) fCaloPID[i]=pid[i];
	135	SetProbPhot(fCaloPID[kPhotonAli]);
	136	}else{
	137	for(Int_t i=0; i<kUnknownAli+1; fCaloPID[i++]=0.){}
	138	fCaloPID[kUnknownAli]=1.;
	139	}
	140	}
	141
	142
	143	//______________________________________________________________________________
	144	particleType AliJPhoton::GetParticleType() {
	145
	146	// return the most problable particle type
	147	// Note: following the AliCaloPID implementation
	148
	149	//Init default weights
	150	Float_t wPhoton = 0.75 ;
	151	Float_t wPi0 = 0.8 ;
	152	Float_t wElectron = 0.5 ;
	153	Float_t wCharged = 0.5 ;
	154	Float_t wNeutral = 0.5 ;
	155
	156	Bool_t usePHOSweightFormula = kTRUE;
	157	//Formula to set the PID weight threshold for photon or pizero
	158	// ALICE-INT-2005-016 (2007)
	159	TFormula* wPhotonPHOSFormula =
	160	new TFormula("photonWeight","0.75(x<40)+ 0.68(x>=100)+(x>=40 && x<100)(0.98+x(6e-3)-xx(2e-04)+xxx*(1.1e-06))");
	161	TFormula* wPi0PHOSFormula =
	162	new TFormula("pi0Weight","0.80(x<65)+ 0.915(x>=100)+(x>=65 && x-x(1.95e-3)-xx(4.31e-05)+xxx(3.61e-07))");
	163
	164	if(fCaloType == kPHOSCalo && usePHOSweightFormula){
	165	wPhoton = wPhotonPHOSFormula->Eval(E()) ;
	166	wPi0 = wPi0PHOSFormula->Eval(E());
	167	// cout << "wPhotonF " << wPhoton << " wphoton " << fCaloPID[kPhotonAli] << " wPi0F " << wPi0 << " wpi0 " << fCaloPID[kPi0Ali] << endl;
	168	}
	169
	170	if(fCaloType == kEMCALCalo){
	171	wPhoton = 0.8 ;
	172	wPi0 = 0.5 ;
	173	wElectron = 0.8 ;
	174	wCharged = 0.5 ;
	175	wNeutral = 0.5 ;
	176	}
	177
	178
	179	particleType pdg = kHadron;
	180
	181	//Select most probable ID
	182	if(fCaloType == kPHOSCalo){
	183	if(fCaloPID[kPhotonAli] > wPhoton) pdg = kPhoton ;
	184	else if(fCaloPID[kPi0Ali] > wPi0) pdg = kPizero ;
	185	//else if(fCaloPID[kElectronAli] > wElectron) pdg = electron ;
	186	//else if(fCaloPID[kEleConAli] > wElectron) pdg = electronCon ;
	187	//else if(chargedHadronWeight > wCharged) pdg = chargedHadron ;
	188	//else if(neutralHadronWeight > wNeutral) pdg = neutralHadron ;
	189	//else if(allChargedWeight > allNeutralWeight)
	190	// pdg = chargedUnknown ;
	191	//else
	192	// pdg = neutralUnknown ;
	193	}
	194	else{//EMCAL
	195	//Temporal solution, electrons and photons not differenciated
196	if(fCaloPID[kPhotonAli] + fCaloPID[kElectronAli] > wPhoton) pdg = kPhoton ;
197	else if(fCaloPID[kPi0Ali] > wPi0) pdg = kPizero ;
198	//else if(chargedHadronWeight + neutralHadronWeight > wCharged) pdg = chargedHadron ;
199	//else if(neutralHadronWeight + chargedHadronWeight > wNeutral) pdg = neutralHadron ;
200	//else pdg = neutralUnknown ;
201
202	}
203
204	delete wPhotonPHOSFormula;
205	delete wPi0PHOSFormula;
206
207	//if(fCaloType == kPHOSCalo) cout << "pdg " << pdg << endl;
208
209	return pdg ;
210
211	}
212
213