#include "TTree.h"
#include "TBRIK.h"
#include "TNode.h"
+#include "TMath.h"
// --- Standard library ---
#include <stdio.h>
AliPHOS::~AliPHOS(void)
{
+ delete fHits; // 28.12.1998
+ delete fTreePHOS; // 28.12.1998
fCradles->Delete();
delete fCradles;
}
//______________________________________________________________________________
-AliPHOSCradle *AliPHOS::GetCradleOfTheParticle(const Hep3Vector &p,const Hep3Vector &v) const
+AliPHOSCradle *AliPHOS::GetCradleOfTheParticle(const TVector3 &p,const TVector3 &v) const
{
// For a given direction 'p' and source point 'v' returns pointer to AliPHOSCradle
// in that direction or NULL if AliPHOSCradle was not found.
const float d = cradle->GetRadius()-cradle->GetCPV_PHOS_Distance()-cradle->GetCPV_Thikness();
cradle->GetXY(p,v,d,x,y,l);
- if( l>0 && fabs(x)<cradle->GetNz ()*cradle->GetCellSideSize()/2
- && fabs(y)<cradle->GetNphi()*cradle->GetCellSideSize()/2 )
+ if( l>0 && TMath::Abs(x)<cradle->GetNz ()*cradle->GetCellSideSize()/2
+ && TMath::Abs(y)<cradle->GetNphi()*cradle->GetCellSideSize()/2 )
return cradle;
}
//______________________________________________________________________________
+AliPHOSCradle::~AliPHOSCradle(void) // 28.12.1998
+{
+ fGammasReconstructed.Delete();
+ fParticles .Delete();
+}
+
+//______________________________________________________________________________
+
void AliPHOSCradle::Clear(Option_t *)
{
// Clear digit. information.
//______________________________________________________________________________
-void AliPHOSCradle::GetXY(const Hep3Vector &p,const Hep3Vector &v,float R,float &x,float &y,float &l) const
+void AliPHOSCradle::GetXY(const TVector3 &p,const TVector3 &v,float R,float &x,float &y,float &l) const
{
// This function calculates hit position (x,y) in the CRADLE cells plain from particle in
// the direction given by 'p' (not required to be normalized) and start point
//
// Example:
// AliPHOSCradle cradle(......);
-// Hep3Vector p(....), v(....);
+// TVector3 p(....), v(....);
// Float_t x,y,l;
// cradle.GetXY(p,v,x,y,l);
-// if( l<0 || fabs(x)>cradle.GetNz() *cradle.GetCellSideSize()/2
-// || fabs(y)>cradle.GetNphi()*cradle.GetCellSideSize()/2 )
+// if( l<0 || TMath::Abs(x)>cradle.GetNz() *cradle.GetCellSideSize()/2
+// || TMath::Abs(y)>cradle.GetNphi()*cradle.GetCellSideSize()/2 )
// cout << "Outside the CRADLE.\n";
// We have to create three vectors:
// n2 - second vector in CRADLE plain
// This three vectors are orthonormalized.
- double phi = fPhi/180*M_PI;
- Hep3Vector n1( 0 , 0 , 1 ), // Z direction (X)
+ double phi = fPhi/180*TMath::Pi();
+ TVector3 n1( 0.0 , 0.0 , 1.0 ), // Z direction (X)
n2( -sin(phi) , cos(phi) , 0 ), // around beam (Y)
s ( R*cos(phi) , R*sin(phi) , 0 ); // central point
const double l1_min = 1e-2;
double l1,
- p_n1 = p.dot(n1), // dot() - scalar product.
- p_n2 = p.dot(n2),
- v_n1 = v.dot(n1),
- v_n2 = v.dot(n2),
- s_n1 = s.dot(n1), // 0
- s_n2 = s.dot(n2); // 0
+ p_n1 = p*n1, // * - scalar product.
+ p_n2 = p*n2,
+ v_n1 = v*n1,
+ v_n2 = v*n2,
+ s_n1 = s*n1, // 0
+ s_n2 = s*n2; // 0
- if ( fabs(l1=p.x()-n1.x()*p_n1-n2.x()*p_n2)>l1_min )
- { l = (-v.x()+s.x()+n1.x()*(v_n1-s_n1)+n2.x()*(v_n2-s_n2))/l1; }
- else if ( fabs(l1=p.y()-n1.y()*p_n1-n2.y()*p_n2)>l1_min )
- { l = (-v.y()+s.y()+n1.y()*(v_n1-s_n1)+n2.y()*(v_n2-s_n2))/l1; }
- else if ( fabs(l1=p.z()-n1.z()*p_n1-n2.z()*p_n2)>l1_min )
- { l = (-v.z()+s.z()+n1.z()*(v_n1-s_n1)+n2.z()*(v_n2-s_n2))/l1; }
-
-// double lx = (-v.x()+s.x()+n1.x()*(v.dot(n1)-s.dot(n1))+n2.x()*(v.dot(n2)-s.dot(n2)))/
-// (p.x()-n1.x()*p.dot(n1)-n2.x()*p.dot(n2)),
-// ly = (-v.y()+s.y()+n1.y()*(v.dot(n1)-s.dot(n1))+n2.y()*(v.dot(n2)-s.dot(n2)))/
-// (p.y()-n1.y()*p.dot(n1)-n2.y()*p.dot(n2)),
-// lz = (-v.z()+s.z()+n1.z()*(v.dot(n1)-s.dot(n1))+n2.z()*(v.dot(n2)-s.dot(n2)))/
-// (p.z()-n1.z()*p.dot(n1)-n2.z()*p.dot(n2));
-// cout.form("x: %g %g %g %g\n",lx,-v.x()+s.x()+n1.x()*(v.dot(n1)-s.dot(n1))+n2.x()*(v.dot(n2)-s.dot(n2)),p.x()-n1.x()*p.dot(n1)-n2.x()*p.dot(n2));
-// cout.form("y: %g %g %g %g\n",lx,-v.y()+s.y()+n1.y()*(v.dot(n1)-s.dot(n1))+n2.y()*(v.dot(n2)-s.dot(n2)),p.y()-n1.y()*p.dot(n1)-n2.y()*p.dot(n2));
-// cout.form("z: %g %g %g %g\n",lx,-v.z()+s.z()+n1.z()*(v.dot(n1)-s.dot(n1))+n2.z()*(v.dot(n2)-s.dot(n2)),p.z()-n1.z()*p.dot(n1)-n2.z()*p.dot(n2));
+ if ( TMath::Abs(l1=p.X()-n1.X()*p_n1-n2.X()*p_n2)>l1_min )
+ { l = (-v.X()+s.X()+n1.X()*(v_n1-s_n1)+n2.X()*(v_n2-s_n2))/l1; }
+ else if ( TMath::Abs(l1=p.Y()-n1.Y()*p_n1-n2.Y()*p_n2)>l1_min )
+ { l = (-v.Y()+s.Y()+n1.Y()*(v_n1-s_n1)+n2.Y()*(v_n2-s_n2))/l1; }
+ else if ( TMath::Abs(l1=p.Z()-n1.Z()*p_n1-n2.Z()*p_n2)>l1_min )
+ { l = (-v.Z()+s.Z()+n1.Z()*(v_n1-s_n1)+n2.Z()*(v_n2-s_n2))/l1; }
+
+// double lx = (-v.X()+s.X()+n1.X()*(v.dot(n1)-s.dot(n1))+n2.X()*(v.dot(n2)-s.dot(n2)))/
+// (p.X()-n1.X()*p.dot(n1)-n2.X()*p.dot(n2)),
+// ly = (-v.Y()+s.Y()+n1.Y()*(v.dot(n1)-s.dot(n1))+n2.Y()*(v.dot(n2)-s.dot(n2)))/
+// (p.Y()-n1.Y()*p.dot(n1)-n2.Y()*p.dot(n2)),
+// lz = (-v.Z()+s.Z()+n1.Z()*(v.dot(n1)-s.dot(n1))+n2.Z()*(v.dot(n2)-s.dot(n2)))/
+// (p.Z()-n1.Z()*p.dot(n1)-n2.Z()*p.dot(n2));
+// cout.form("x: %g %g %g %g\n",lx,-v.X()+s.X()+n1.X()*(v.dot(n1)-s.dot(n1))+n2.X()*(v.dot(n2)-s.dot(n2)),p.X()-n1.X()*p.dot(n1)-n2.X()*p.dot(n2));
+// cout.form("y: %g %g %g %g\n",lx,-v.Y()+s.Y()+n1.Y()*(v.dot(n1)-s.dot(n1))+n2.Y()*(v.dot(n2)-s.dot(n2)),p.Y()-n1.Y()*p.dot(n1)-n2.Y()*p.dot(n2));
+// cout.form("z: %g %g %g %g\n",lx,-v.Z()+s.Z()+n1.Z()*(v.dot(n1)-s.dot(n1))+n2.Z()*(v.dot(n2)-s.dot(n2)),p.Z()-n1.Z()*p.dot(n1)-n2.Z()*p.dot(n2));
// cout.form("lx,ly,lz = %g,%g,%g\n",lx,ly,lz);
x = p_n1*l + v_n1 - s_n1;
// signal_step=0.001 GeV (1MeV)
// min_signal_reject = 15 or 30 MeV
+
common_for_event_storing.event_number = 0; // We do not know event number?
common_for_event_storing.crystals_matrix_amount_PHOS = 1;
common_for_event_storing.crystal_matrix_type = 1; // 1 - rectangular
const float stochastic_term = 0.03, // per cents over sqrt(E); E in GeV
electronic_noise = 0.01; // GeV
reconsfirst(stochastic_term,electronic_noise); // Call of reconstruction program.
-
for( int i=0; i<rcgamma.recons_gammas_amount; i++ )
{
Float_t thetta, alpha, betta, R=fRadius+rcgamma.Zgdev[i]/10;
g.fX = rcgamma.YW[i]/10;
- g.fXsigma = 0.2;
g.fY = rcgamma.XW[i]/10;
- g.fYsigma = 0.2;
g.fE = rcgamma.E [i];
- g.fEsigma = 0.01*sqrt(rcgamma.E[i])+0.05;
thetta = atan(g.fX/R);
alpha = atan(g.fY/R);
- betta = fPhi/180*M_PI + alpha;
+ betta = fPhi/180*TMath::Pi() + alpha;
g.fPx = g.fE * cos(thetta) * cos(betta);
g.fPy = g.fE * cos(thetta) * sin(betta);
else
mass = -sqrt(-mass);
- printf("XY=(%+7.2f,%+7.2f) (%+7.2f,%+7.2f,%+7.2f;%7.2f) mass=%8.4f\n",
- fX,fY,fPx,fPy,fPz,fE,mass);
+ printf("XY=(%+7.2f,%+7.2f) (%+7.2f,%+7.2f,%+7.2f;%7.2f) mass=%8.4f Ipart=%2d\n",
+ fX,fY,fPx,fPy,fPz,fE,mass,fIpart);
}
//______________________________________________________________________________
AliPHOSgamma &AliPHOSgamma::operator=(const AliPHOSgamma &g)
{
fX = g.fX;
- fXsigma = g.fXsigma;
fY = g.fY;
- fYsigma = g.fYsigma;
fE = g.fE;
- fEsigma = g.fEsigma;
fPx = g.fPx;
fPy = g.fPy;
fPz = g.fPz;
+ fIpart = g.fIpart;
return *this;
}
// Manager and hits classes for set:PHOS //
////////////////////////////////////////////////
-// --- CLHEP ---
-#include <CLHEP/Vector/ThreeVector.h>
-#include <CLHEP/Vector/LorentzVector.h>
-
// --- ROOT system ---
#include <TArray.h>
#include <TRandom.h>
#include <TH2.h>
+#include <TVector3.h>
// --- galice header files ---
#include "AliDetector.h"
virtual ~AliPHOSgamma(void) {}
AliPHOSgamma(void) {}
AliPHOSgamma(const AliPHOSgamma &g) { *this=g; }
- AliPHOSgamma(Float_t X, Float_t Xsigma,
- Float_t Y, Float_t Ysigma,
- Float_t E, Float_t Esigma,
- Float_t Px, Float_t Py, Float_t Pz) :
- fX(X), fXsigma(Xsigma),
- fY(Y), fYsigma(Ysigma),
- fE(E), fEsigma(Esigma),
- fPx(Px), fPy(Py), fPz(Pz)
+ AliPHOSgamma(Float_t X, Float_t Y, Float_t E,
+ Float_t Px, Float_t Py, Float_t Pz,
+ Int_t Ipart) :
+ fX(X), fY(Y), fE(E),
+ fPx(Px), fPy(Py), fPz(Pz),
+ fIpart(Ipart)
{}
Float_t fX; // cm. x-coordinate (in beam direction)
- Float_t fXsigma; // cm. x-coordinate error
-
Float_t fY; // cm. y-coordinate (around beam)
- Float_t fYsigma; // cm. y-coordinate error
Float_t fE; // GeV. energy
- Float_t fEsigma; // GeV. energy error
Float_t fPx; // GeV. Gamma momentum Px
Float_t fPy; // GeV. Gamma momentum Py
Float_t fPz; // GeV. Gamma momentum Pz
+ Int_t fIpart; // Current particle number (GEANT particle code)
+
void Print(Option_t *options=NULL);
AliPHOSgamma &operator=(const AliPHOSgamma &g);
public:
- virtual ~AliPHOSCradle(void) {}
+ virtual ~AliPHOSCradle(void);
AliPHOSCradle(void);
AliPHOSCradle(int Geometry ,
float CrystalSideSize ,
void Reconstruction(Float_t signal_step, UInt_t min_signal_reject);
- void GetXY(const Hep3Vector &p,const Hep3Vector &v,float R,float &x,float &y,float &l) const;
+ void GetXY(const TVector3 &p,const TVector3 &v,float R,float &x,float &y,float &l) const;
TObjArray &GetGammasReconstructed (void) {return fGammasReconstructed;}
TObjArray &GetParticles (void) {return fParticles;}
void FinishRun(void);
void ResetDigits(void);
void Print(Option_t *opt="");
- AliPHOSCradle *GetCradleOfTheParticle(const Hep3Vector &p,const Hep3Vector &v) const;
+ AliPHOSCradle *GetCradleOfTheParticle(const TVector3 &p,const TVector3 &v) const;
AliPHOSCradle &GetCradle(int n) {return *(AliPHOSCradle*)fCradles->operator[](n);}
// AliPHOSCradle &GetCradle(int n) {return *((AliPHOSCradle*) (*fCradles)[n]) ;}
void Reconstruction(Float_t signal_step, UInt_t min_signal_reject);
// --- galice header files ---
#include "AliPHOSv0.h"
#include "AliRun.h"
-#include "AliConst.h"
#include "AliMC.h"
ClassImp(AliPHOSv0)
yp1 = -r * TMath::Cos(pphi * 3.);
xp2 = -r * TMath::Sin(pphi);
yp2 = -r * TMath::Cos(pphi);
- pphi *= kRaddeg;
+ pphi *= 180/kPI;
AliMatrix(idrotm[0], 90.,-3*pphi, 90., 90-3*pphi, 0., 0.);
AliMatrix(idrotm[1], 90., -pphi, 90., 90-pphi, 0., 0.);
AliMatrix(idrotm[2], 90., pphi, 90., 90+pphi, 0., 0.);
#include "TTree.h"
#include "TBRIK.h"
#include "TNode.h"
+#include "TMath.h"
// --- galice header files ---
#include "AliPHOSv1.h"
const float cell_length = GetCrystalLength()+GetAirThickness()+GetWrapThickness()+GetPIN_Length(),
cell_side_size = GetCrystalSideSize()+2*GetAirThickness()+2*GetWrapThickness(),
-// cell_angle = 180/M_PI * 2 * atan(cell_side_size/2 / GetRadius()), // radians
+// cell_angle = 180/kPI * 2 * atan(cell_side_size/2 / GetRadius()), // radians
cradle_thikness = cell_length + GetCPV_Thickness() + GetCPV_PHOS_Distance(),
distance_to_CPV = GetRadius() - GetCPV_Thickness() - GetCPV_PHOS_Distance();
{
float c = distance_to_CPV, // Distance to CPV
l = cell_side_size*GetNphi()/2, // Cradle half size around beam (for rect. geom.)
- cradle_angle = 360/M_PI*atan(l/c),
+ cradle_angle = 360/kPI*atan(l/c),
cradle_angle_pos = -90+(i-(GetCradlesAmount()-1)/2.) * (cradle_angle+GetAngleBetweenCradles());
// Cradles are numerated in clock reversed order. (general way of angle increment)
float r = GetRadius() + cradle_thikness/2;
- x = r*cos(cradle_angle_pos*M_PI/180);
- y = r*sin(cradle_angle_pos*M_PI/180);
+ x = r*cos(cradle_angle_pos*kPI/180);
+ y = r*sin(cradle_angle_pos*kPI/180);
z = 0;
AliMatrix(rotation_matrix_number, 0,0 , 90,90+cradle_angle_pos , 90,180+cradle_angle_pos);
pMC->Gspos("PHOS",i+1,"ALIC",x,y,z,rotation_matrix_number,"ONLY");
AliMC *MC = AliMC::GetMC();
Int_t copy;
-// cout << form("PHOS STEP: %d %d %d %d\n",
-// MC->GetMedium(),cvolu->nlevel,ctrak->inwvol,inwold);
-
int cradle_number, cell_Z, cell_Phi; // Variables that describe cell position.
if( MC->GetMedium() == GetPHOS_IDTMED_PIN() && (MC->TrackInside() || MC->TrackExiting()==2) && inwold && MC->TrackCharge()!=0 )
//////////////////////////////////////////////////////////////////////////////
- if( MC->GetMedium() == GetPHOS_IDTMED_PbWO4() ) //&& gckine_.CHARGE )
+ if( MC->GetMedium() == GetPHOS_IDTMED_PbWO4() )
{
// GEANT particle into crystal.
Float_t xyz[3];
MC->TrackPosition(xyz);
- Hep3Vector p(xyz[0],xyz[1],xyz[2]),v;
+ TVector3 p(xyz[0],xyz[1],xyz[2]),v;
float x,y,l;
float R = cradle.GetRadius() - cradle.GetCPV_PHOS_Distance() - cradle.GetCPV_Thikness();
float Px = pmom[0] * pmom[3],
Py = pmom[1] * pmom[3],
Pz = pmom[2] * pmom[3];
+ Int_t Ipart = MC->TrackPid();
// TClonesArray &P=cradle.GetParticles();
// new( P[P.GetEntries()] ) AliPHOSgamma(x,0,y,0,ctrak->getot,0,Px,Py,Pz);
- cradle.GetParticles().Add(new AliPHOSgamma(x,0,y,0,MC->Etot(),0,Px,Py,Pz));
+ cradle.GetParticles().Add(new AliPHOSgamma(x,y,MC->Etot(),Px,Py,Pz,Ipart));
if( MC->TrackCharge()!=0 )
cradle.AddCPVHit(x,y);
--- /dev/null
+/////////////////////////////////////////////////////////
+// Manager and hits classes for set:PHOS version 3 //
+/////////////////////////////////////////////////////////
+
+// --- ROOT system ---
+#include "TH1.h"
+#include "TRandom.h"
+#include "TFile.h"
+#include "TTree.h"
+#include "TBRIK.h"
+#include "TNode.h"
+
+// --- galice header files ---
+#include "AliPHOSv3.h"
+#include "AliRun.h"
+#include "TGeant3.h"
+
+ClassImp(AliPHOSv3)
+
+//______________________________________________________________________________
+
+
+AliPHOSv3::AliPHOSv3() : AliPHOS()
+{
+}
+
+//______________________________________________________________________________
+
+AliPHOSv3::AliPHOSv3(const char *name, const char *title)
+ : AliPHOS(name, title)
+{
+}
+
+//___________________________________________
+void AliPHOSv3::CreateGeometry()
+{
+
+ cout << "AliPHOSv3::CreateGeometry() PHOS creation\n";
+
+ AliMC* pMC = AliMC::GetMC();
+
+ AliPHOS *PHOS_tmp = (AliPHOS*)gAlice->GetModule("PHOS");
+ if( NULL==PHOS_tmp )
+ {
+ printf("There isn't PHOS detector!\n");
+ return;
+ }
+// AliPHOS &PHOS = *PHOS_tmp;
+
+ //////////////////////////////////////////////////////////////////////////////
+
+ Int_t rotation_matrix_number=0;
+ Float_t par[11],
+ x,y,z;
+
+ const float cell_length = GetCrystalLength()+GetAirThickness()+GetWrapThickness()+GetPIN_Length(),
+ cell_side_size = GetCrystalSideSize()+2*GetAirThickness()+2*GetWrapThickness(),
+// cell_angle = 180/kPI * 2 * atan(cell_side_size/2 / GetRadius()), // radians
+ cradle_thikness = cell_length + GetCPV_Thickness() + GetCPV_PHOS_Distance(),
+ distance_to_CPV = GetRadius() - GetCPV_Thickness() - GetCPV_PHOS_Distance();
+
+ //////////////////////////////////////////////////////////////////////////////
+ // CELL volume and subvolumes creation
+ //////////////////////////////////////////////////////////////////////////////
+
+ par[0] = GetCrystalSideSize()/2 + GetWrapThickness();
+ par[1] = GetCrystalSideSize()/2 + GetWrapThickness();
+ par[2] = GetCrystalLength() /2 + GetWrapThickness()/2;
+ pMC->Gsvolu("WRAP","BOX ",GetPHOS_IDTMED_Tyvek(),par,3);
+
+ par[0] = GetCrystalSideSize()/2;
+ par[1] = GetCrystalSideSize()/2;
+ par[2] = GetCrystalLength()/2;
+ pMC->Gsvolu("CRST","BOX ",GetPHOS_IDTMED_PbWO4(),par,3);
+
+ // PIN
+ par[0] = GetPIN_SideSize()/2;
+ par[1] = GetPIN_SideSize()/2;
+ par[2] = GetPIN_Length()/2;
+ pMC->Gsvolu("PIN ","BOX ",GetPHOS_IDTMED_PIN(),par,3);
+
+ //////////////////////////////////////////////////////////////////////////////
+ // CRADLE,CPV creation.
+ //////////////////////////////////////////////////////////////////////////////
+
+ par[0] = cell_side_size/2 * GetNz();
+ par[1] = cell_side_size/2 * GetNphi();
+ par[2] = cradle_thikness/2;
+ pMC->Gsvolu("PHOS","BOX ",GetPHOS_IDTMED_AIR(),par,3);
+
+//par[0] : the same as above
+//par[1] : the same as above
+ par[2] = GetCPV_Thickness()/2;
+ pMC->Gsvolu("CPV ","BOX ",GetPHOS_IDTMED_CPV(),par,3);
+
+ x = 0;
+ y = 0;
+ z = (cell_length+GetCPV_PHOS_Distance())/2;
+ pMC->Gspos("CPV ",1,"PHOS",x,y,z,0,"ONLY");
+
+ par[0] = cell_side_size/2 * GetNz();
+ par[1] = cell_side_size/2 * GetNphi();
+ par[2] = cell_length/2;
+ pMC->Gsvolu("CRS0","BOX ",GetPHOS_IDTMED_AIR(),par,3);
+
+ x = 0;
+ y = 0;
+ z = -(cradle_thikness-cell_length)/2;
+ pMC->Gspos("CRS0",1,"PHOS",x,y,z,0,"ONLY");
+
+ pMC->Gsdvn("CRS1","CRS0",GetNphi(),2);
+ pMC->Gsdvn("CELL","CRS1",GetNz() ,1);
+
+ //////////////////////////////////////////////////////////////////////////////
+ // CELL creation
+ //////////////////////////////////////////////////////////////////////////////
+
+ x = 0;
+ y = 0;
+ z = -GetWrapThickness()/2;
+ pMC->Gspos("CRST",1,"WRAP",x,y,z,0,"ONLY");
+
+ x = 0;
+ y = 0;
+ z = GetPIN_Length()/2;
+ pMC->Gspos("WRAP",1,"CELL",x,y,z,0,"ONLY");
+
+ x = 0;
+ y = 0;
+ z = -GetCrystalLength()/2-GetWrapThickness()/2;
+ pMC->Gspos("PIN ",1,"CELL",x,y,z,0,"ONLY");
+
+ //////////////////////////////////////////////////////////////////////////////
+ // CELL has been created.
+ //////////////////////////////////////////////////////////////////////////////
+
+// int n=0;
+// z = -(GetCPV_Thickness()+GetCPV_PHOS_Distance())/2;
+//
+// for( int iy=0; iy<GetNphi(); iy++ )
+// {
+// y = (iy-(GetNphi()-1)/2.)*cell_side_size;
+// for( int ix=0; ix<GetNz(); ix++ )
+// {
+// x = (ix-(GetNz()-1)/2.)*cell_side_size;
+// pMC->Gspos("CELL",++n,"PHOS",x,y,z,0,"ONLY");
+// }
+// }
+
+ //////////////////////////////////////////////////////////////////////////////
+ // End of CRADLE creation.
+ //////////////////////////////////////////////////////////////////////////////
+
+
+ //////////////////////////////////////////////////////////////////////////////
+ // PHOS creation
+ //////////////////////////////////////////////////////////////////////////////
+
+ for( int i=0; i<GetCradlesAmount(); i++ )
+ {
+ float c = distance_to_CPV, // Distance to CPV
+ l = cell_side_size*GetNphi()/2, // Cradle half size around beam (for rect. geom.)
+ cradle_angle = 360/kPI*atan(l/c),
+ cradle_angle_pos = -90+(i-(GetCradlesAmount()-1)/2.) * (cradle_angle+GetAngleBetweenCradles());
+ // Cradles are numerated in clock reversed order. (general way of angle increment)
+
+ float r = GetRadius() + cradle_thikness/2;
+ x = r*cos(cradle_angle_pos*kPI/180);
+ y = r*sin(cradle_angle_pos*kPI/180);
+ z = 0;
+ AliMatrix(rotation_matrix_number, 0,0 , 90,90+cradle_angle_pos , 90,180+cradle_angle_pos);
+ pMC->Gspos("PHOS",i+1,"ALIC",x,y,z,rotation_matrix_number,"ONLY");
+
+ GetCradleAngle(i) = cradle_angle_pos;
+//
+// int n = PHOS.fCradles->GetEntries();
+// PHOS.fCradles->Add(new AliPHOSCradle( 1, // geometry.
+// GetCrystalSideSize (),
+// GetCrystalLength (),
+// GetWrapThickness (),
+// GetAirThickness (),
+// GetPIN_SideSize (),
+// GetPIN_Length (),
+// GetRadius (),
+// GetCPV_Thickness (),
+// GetCPV_PHOS_Distance (),
+// GetNz (),
+// GetNphi (),
+// cradle_angle_pos ));
+//
+// if( n+1 != PHOS.fCradles->GetEntries() ||
+// NULL == PHOS.fCradles->At(n) )
+// {
+// cout << " Can not create or add AliPHOSCradle.\n";
+// exit(1);
+// }
+ }
+ AddPHOSCradles();
+
+ //////////////////////////////////////////////////////////////////////////////
+ // All is done.
+ // Print some information.
+ //////////////////////////////////////////////////////////////////////////////
+}
+
+void AliPHOSv3::StepManager()
+{
+ static Bool_t inwold=0; // Status of previous ctrak->inwvol
+ AliMC *MC = AliMC::GetMC();
+ Int_t copy;
+
+// if( MC->TrackEntering() ) {
+// Int_t Volume_ID = MC->CurrentVol(Volume_name, copy);
+// cout << "AliPHOSv3::StepManager() entered to PHOS to the volume " << Volume_name << "!\n";
+// }
+
+ int cradle_number, cell_Z, cell_Phi; // Variables that describe cell position.
+
+ if( MC->GetMedium()==GetPHOS_IDTMED_PIN() && MC->TrackEntering() && MC->TrackCharge()!=0 )
+ {
+ // GEANT particle just have entered into PIN diode.
+
+ AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
+
+ MC->CurrentVolOff(4,0,copy);
+ cradle_number = copy-1;
+ MC->CurrentVolOff(1,0,copy);
+ cell_Z = copy-1;
+ MC->CurrentVolOff(2,0,copy);
+ cell_Phi = copy-1;
+
+ TH2S &h = PHOS.GetCradle(cradle_number).fChargedTracksInPIN;
+ h.AddBinContent(h.GetBin(cell_Z,cell_Phi));
+
+// cout << "AliPHOSv3::StepManager() entered to PHOS pin diode\n";
+// cout << " cradle_nimber = " << cradle_number << endl;
+// cout << " cell_z = " << cell_Z << endl;
+// cout << " cell_Phi = " << cell_Phi << endl;
+ }
+
+ //////////////////////////////////////////////////////////////////////////////
+
+ if( MC->GetMedium() == GetPHOS_IDTMED_PbWO4() )
+ {
+ // GEANT particle into crystal.
+
+ AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
+
+ MC->CurrentVolOff(5,0,copy);
+ cradle_number = copy-1;
+ MC->CurrentVolOff(2,0,copy);
+ cell_Z = copy-1;
+ MC->CurrentVolOff(3,0,copy);
+ cell_Phi = copy-1;
+
+ TH2F &h = PHOS.GetCradle(cradle_number).fCellEnergy;
+ h.AddBinContent(h.GetBin(cell_Z,cell_Phi),MC->Edep());
+ }
+
+ //////////////////////////////////////////////////////////////////////////////
+
+ if( MC->GetMedium()==GetPHOS_IDTMED_CPV() && MC->TrackEntering() )
+ {
+ // GEANT particle just have entered into CPV detector.
+
+ AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
+
+ MC->CurrentVolOff(1,0,cradle_number);
+ cradle_number--;
+
+ // Save CPV x,y hits position of charged particles.
+
+ AliPHOSCradle &cradle = PHOS.GetCradle(cradle_number);
+
+ Float_t xyz[3];
+ MC->TrackPosition(xyz);
+ TVector3 p(xyz[0],xyz[1],xyz[2]),v;
+
+ float x,y,l;
+ float R = cradle.GetRadius() - cradle.GetCPV_PHOS_Distance() - cradle.GetCPV_Thikness();
+ cradle.GetXY(p,v,R,x,y,l);
+ if( PHOS.fDebugLevel>0 )
+ if( l<0 )
+ printf("PHOS_STEP: warning: negative distance to CPV!! %f\n", l);
+
+ // Store current particle in the list of Cradle particles.
+ Float_t pmom[4];
+ MC->TrackMomentum(pmom);
+ Float_t Px = pmom[0] * pmom[3],
+ Py = pmom[1] * pmom[3],
+ Pz = pmom[2] * pmom[3];
+ Float_t Getot = MC->Etot();
+ Int_t Ipart = MC->TrackPid();
+
+ cradle.GetParticles().Add(new AliPHOSgamma(x,y,Getot,Px,Py,Pz,Ipart));
+
+ // printf ("Cradle %i, x,y = %8.3f, %8.3f cm, E,Px,Py,Pz = %8.3f, %8.3f, %8.3f GeV, %8.3f, Ipart = %i\n",
+ // cradle_number,x,y,Getot,Px,Py,Pz,Ipart);
+
+ }
+
+ inwold=MC->TrackEntering(); // Save current status of GEANT variable.
+}
--- /dev/null
+#ifndef PHOSv3_H
+#define PHOSv3_H
+////////////////////////////////////////////////////////
+// Manager and hits classes for set:PHOS version 1 //
+////////////////////////////////////////////////////////
+
+// --- galice header files ---
+#include "AliPHOS.h"
+
+class AliPHOSv3 : public AliPHOS {
+
+ public:
+ AliPHOSv3();
+ AliPHOSv3(const char *name, const char *title);
+ virtual ~AliPHOSv3(){}
+ virtual void CreateGeometry();
+ virtual Int_t IsVersion() const {return 3;}
+ virtual void StepManager();
+
+ ClassDef(AliPHOSv3,1) //Hits manager for set:PHOS version 3
+};
+
+#endif
+
# C++ sources
-SRCS = AliPHOS.cxx AliPHOSv0.cxx AliPHOSv1.cxx AliPHOSv2.cxx
+SRCS = AliPHOS.cxx AliPHOSv0.cxx AliPHOSv1.cxx AliPHOSv2.cxx AliPHOSv3.cxx
# C++ Headers
#pragma link C++ class AliPHOSv0;
#pragma link C++ class AliPHOSv1;
#pragma link C++ class AliPHOSv2;
+#pragma link C++ class AliPHOSv3;
#pragma link C++ class AliPHOShitv2;
#pragma link C++ class AliPHOSCradle;
#pragma link C++ class AliPHOSgamma;