RICH lib splitted
authorkir <kir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 15 Dec 2004 08:51:41 +0000 (08:51 +0000)
committerkir <kir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 15 Dec 2004 08:51:41 +0000 (08:51 +0000)
20 files changed:
RICH/AliRICH.cxx
RICH/AliRICH.h
RICH/AliRICHChamber.cxx
RICH/AliRICHParam.cxx
RICH/AliRICHParam.h
RICH/AliRICHTracker.cxx
RICH/AliRICHv0.cxx
RICH/AliRICHv0.h
RICH/AliRICHv1.h
RICH/Opticals.h
RICH/RICHLinkDef.h [deleted file]
RICH/RICHbaseLinkDef.h [new file with mode: 0644]
RICH/RICHrecLinkDef.h [new file with mode: 0644]
RICH/RICHsimLinkDef.h [new file with mode: 0644]
RICH/RichMake
RICH/libRICH.pkg [deleted file]
RICH/libRICHbase.pkg [new file with mode: 0644]
RICH/libRICHrec.pkg [new file with mode: 0644]
RICH/libRICHsim.pkg [new file with mode: 0644]
RICH/menu.C

index bb0431b..7601834 100644 (file)
@@ -171,15 +171,15 @@ void AliRICH::BuildGeometry()
   TNode *node, *subnode, *top;
   top=gAlice->GetGeometry()->GetNode("alice");
 
-  Float_t widx=P()->SectorSizeX();
-  Float_t leny=P()->SectorSizeY();
-  Float_t dz = P()->Zfreon()+P()->Zwin()+P()->Pc2Win();
-  Float_t dead = P()->DeadZone();
+  Float_t widx =P()->SectorSizeX();
+  Float_t leny =P()->SectorSizeY();
+  Float_t dz   =P()->Zfreon()+P()->Zwin()+P()->Pc2Win();
+  Float_t dead =P()->DeadZone();
 
   new TBRIK("RICH","RICH","void",widx+dead/2,leny+leny/2+dead,dz+0.1); //RICH chamber
   new TBRIK("RPC" ,"RPC" ,"void",widx/2,leny/2,0.01);                  //RICH sector 
 
-  for(int i=1;i<=kNchambers;i++){
+  for(int i=1;i<=P()->Nchambers();i++){
     top->cd();
     node = new TNode(Form("RICH%i",i),Form("RICH%i",i),"RICH",C(i)->Center().X(),C(i)->Center().Y(),C(i)->Center().Z(),C(i)->RotMatrixName());
     node->SetLineColor(kRed);
@@ -207,22 +207,28 @@ void AliRICH::BuildGeometry()
 
   AliDebug(1,"Stop.");    
 }//void AliRICH::BuildGeometry()
-
-//______________________________________________________________________________
+//__________________________________________________________________________________________________
 void AliRICH::CreateMaterials()
 {
 // Definition of available RICH materials  
         
   Int_t   material=0; //tmp material id number
   Float_t a=0,z=0,den=0,radl=0,absl=0; //tmp material parameters
+  
   Float_t tmaxfd=-10.0, deemax=-0.2, stemax=-0.1,epsil=0.001, stmin=-0.001; 
   Int_t   isxfld = gAlice->Field()->Integ();
   Float_t sxmgmx = gAlice->Field()->Max();
     
-  Float_t aAir[4]={12.,14.,16.,36.};  Float_t zAir[4]={6.,7.,8.,18.}; Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};//total 0.9999999
+  Float_t aAir[4]={12,14,16,36};  Float_t zAir[4]={6,7,8,18}; Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};//total 0.9999999
   AliMixture(++material, "RichAir",aAir,zAir,den=0.00120479,4,wAir);                                          //1 (Air) 0.01% C 75% N  23% O 1% Ar
   AliMedium(kAir, "RichAir",material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
+  AliMixture(++material, "RichAerogel",aAir,zAir,den=P()->DenGel(),4,wAir);                     //Aerogel represented by Air
+  AliMedium(kGel, "RichAerogel",material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  
+  AliMixture(++material, "RichAerogelReflector",aAir,zAir,den=P()->DenGel(),4,wAir);           //Aerogel reflector represented by Air
+  AliMedium(kReflector, "RichAerogelReflector",material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  
   AliMaterial(++material, "RichRohacell", a=12.01,z=6.0, den=0.1,     radl=18.8,   absl=0);                   //2 Rohacell 51 C-equiv radl rad cover
   AliMedium(kRoha, "RichRohacell", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
@@ -230,68 +236,73 @@ void AliRICH::CreateMaterials()
   AliMixture(++material, "RichSiO2",aQuartz,zQuartz,den=2.64,-2, wQuartz);                                    //3 Quarz (SiO2) -trasparent rad window
   AliMedium(kSiO2, "RichSiO2",material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
-  Float_t  aFreon[2]={12,19};  Float_t  zFreon[2]={6,9};  Float_t wmatFreon[2]={6,14};
+  Float_t  aFreon[2]={12,19};  Float_t  zFreon[2]={6,9};  Float_t wmatFreon[2]={6,14};                        // C12-6 F19-9   
   AliMixture(++material, "RichC6F14",aFreon,zFreon,den=1.68,-2,wmatFreon);                                    //4 Freon (C6F14) 
   AliMedium(kC6F14, "RichC6F14",material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
   Float_t aMethane[2]={12.01,1}; Float_t zMethane[2]={6,1}; Float_t wMethane[2]={1,4};
   AliMixture (++material, "RichCH4", aMethane, zMethane, den=7.17e-4,-2, wMethane);                        //5,9 methane (CH4) normal and for Gap    
   AliMedium(kCH4, "RichCH4"   , material, 1, isxfld, sxmgmx, tmaxfd, stemax,  deemax, epsil,  stmin);  
-  AliMixture (++material, "RichCH4", aMethane, zMethane, den=7.17e-4,-2, wMethane);                        //5,9 methane (CH4) normal and for Gap    
-  AliMedium(kGap, "RichCH4gap$", material, 1, isxfld, sxmgmx, tmaxfd, 0.1   , -deemax, epsil, -stmin);
+  AliMixture (++material, "RichCH4gap", aMethane, zMethane, den=7.17e-4,-2, wMethane);                      //5,9 methane (CH4) normal and for Gap    
+  AliMedium(kGap, "RichCH4gap", material, 1, isxfld, sxmgmx, tmaxfd, 0.1   , -deemax, epsil, -stmin);
     
   AliMaterial(++material, "RichCsI",      a=12.01,z=6.0, den=0.1,     radl=18.8,   absl=0);                   //6 CsI-radl equivalent
-  AliMedium(kCsI, "RichCsI$", material, 1, isxfld, sxmgmx,tmaxfd, stemax, deemax, epsil, stmin);
+  AliMedium(kCsI, "RichCsI", material, 1, isxfld, sxmgmx,tmaxfd, stemax, deemax, epsil, stmin);
   
-  AliMaterial(++material, "GridCu",    a=63.54,z=29.0,den=8.96,    radl=1.43,   absl=0);                   //7 anode grid (Cu) 
-  AliMedium(kGridCu, "GRID$", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  AliMaterial(++material, "RichGridCu",    a=63.54,z=29.0,den=8.96,    radl=1.43,   absl=0);                   //7 anode grid (Cu) 
+  AliMedium(kGridCu, "RichGridCu", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+    
+  AliMaterial(++material, "RichPcbCu",     a=63.54,z=29.0,den=8.96,    radl=1.4,    absl=0);                   //12 Cu
+  AliMedium(kCu, "RichPcbCu", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
-  AliMixture (++material, "OpSiO2",aQuartz, zQuartz, den=2.64, -2, wQuartz);                             //8 Quarz (SiO2) - opaque
-  AliMedium(kOpSiO2, "QUARTZO$",material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  AliMixture (++material, "RichOpSiO2",aQuartz, zQuartz, den=2.64, -2, wQuartz);                             //8 Quarz (SiO2) - opaque
+  AliMedium(kOpSiO2, "RichOpSiO2",material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
-  AliMaterial(++material, "ALUM",     a=26.98,z=13.0,den=2.699,     radl=8.9,    absl=0);                 //10 aluminium sheet (Al)
-  AliMedium(kAl, "ALUMINUM$", material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  AliMaterial(++material, "RichAl",     a=26.98,z=13.0,den=2.699,     radl=8.9,    absl=0);                 //10 aluminium sheet (Al)
+  AliMedium(kAl, "RichAl", material, 1, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
   Float_t aGlass[5]={12.01,28.09,16,10.8,23}; Float_t zGlass[5]={6,14,8,5,11};  Float_t wGlass[5]={0.5,0.105,0.355,0.03,0.01};
-  AliMixture(++material,"GLASS",aGlass, zGlass, den=1.74, 5, wGlass);                                    //11 Glass 50%-C 10.5%-Si 35.5%-O 3%-B 1%-Na
-  AliMedium(kGlass, "GLASS", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
-  
-  AliMaterial(++material, "COPPER$",  a=63.54,z=29.0,den=8.96,    radl=1.4,    absl=0);                   //12 Cu
-  AliMedium(kCu, "PCB_COPPER", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
-  
+  AliMixture(++material,"RichGlass",aGlass, zGlass, den=1.74, 5, wGlass);                                    //11 Glass 50%-C 10.5%-Si 35.5%-O 3%-B 1%-Na
+  AliMedium(kGlass, "RichGlass", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
-  AliMaterial(++material, "W$",  a=183.84,z=74.0,den=19.3,    radl=0.35,    absl=185.0/den);              //13 W - anod wires
-  AliMedium(kW, "W", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+  AliMaterial(++material, "RichW",  a=183.84,z=74.0,den=19.3,    radl=0.35,    absl=185.0/den);              //13 W - anod wires
+  AliMedium(kW, "RichW", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
   
   if(P()->IsRadioSrc()){
     AliInfo("Special radioactive source materials");
-    AliMaterial(++material, "STEEL$",  a=63.54,z=29.0,den=8.96,    radl=1.4,    absl=0);                  //14 Steel
-    AliMedium(kSteel, "STEEL", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+    AliMaterial(++material, "RichSteel",  a=55.845,z=26.0,den=7.87,    radl=1.76,    absl=131.9/den);        //14 Steel (Fe)
+    AliMedium(kSteel, "RichSteel", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   
-    AliMaterial(++material, "PERPEX$",  a=63.54,z=29.0,den=8.96,    radl=1.4,    absl=0);                 //15 Perpex
-    AliMedium(kPerpex, "PERPEX", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+    AliMaterial(++material, "RichPerpex",  a=63.54,z=29.0,den=8.96,    radl=1.4,    absl=0);                 //15 Perpex
+    AliMedium(kPerpex, "RichPerpex", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
     
-    AliMaterial(++material, "STRONZIUM$",  a=87.62,z=38.0,den=2.54,    radl=4.24,    absl=0);             //16 Sr90
-    AliMedium(kSr90, "STRONZIUM", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+    AliMaterial(++material, "RichSr90",  a=87.62,z=38.0,den=2.54,    radl=4.24,    absl=0);                  //16 Sr90
+    AliMedium(kSr90, "RichSr90", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
+    
+    Float_t aMylar[5]={12.01,1,16}; Float_t zMylar[5]={6,1,8};  Float_t wMylar[5]={5,4,5};                  //17 Mylar C5H4O5
+    AliMixture(++material,"RichMylar",aMylar, zMylar, den=1.39, -3, wMylar); 
+    AliMedium(kMylar, "RichMylar", material, 0, isxfld, sxmgmx, tmaxfd, stemax, deemax, epsil, stmin);
   }
   
 //Optical properties:
 #include "Opticals.h"
-  gMC->SetCerenkov((*fIdtmed)[kAir]      , kNbins, aPckov, aAbsCH4,    aQeAll, aIdxCH4);       //1 Air
-  gMC->SetCerenkov((*fIdtmed)[kRoha]     , kNbins, aPckov, aAbsCH4,    aQeAll, aIdxCH4);       //2 Honeycomb  
-  gMC->SetCerenkov((*fIdtmed)[kSiO2]     , kNbins, aPckov, aAbsSiO2,   aQeAll, aIdxSiO2);      //3 Quartz SiO2 
-  gMC->SetCerenkov((*fIdtmed)[kC6F14]    , kNbins, aPckov, aAbsC6F14,  aQeAll, aIdxC6F14);     //4 Freon C6F14
-  gMC->SetCerenkov((*fIdtmed)[kCH4]      , kNbins, aPckov, aAbsCH4,    aQeAll, aIdxCH4);       //5 Methane CH4 
-  gMC->SetCerenkov((*fIdtmed)[kCsI]      , kNbins, aPckov, aAbsCsI,    aQeCsI, aIdxCH4);       //6 CsI
-  gMC->SetCerenkov((*fIdtmed)[kGridCu]   , kNbins, aPckov, aAbsGrid,   aQeAll, aIdxGrid);      //7 grid Cu
-  gMC->SetCerenkov((*fIdtmed)[kOpSiO2]   , kNbins, aPckov, aAbsOpSiO2, aQeAll, aIdxOpSiO2);    //8 Opaque quartz SiO2
-  gMC->SetCerenkov((*fIdtmed)[kGap]      , kNbins, aPckov, aAbsCH4,    aQeAll, aIdxCH4);       //9 Special methane gap
-  gMC->SetCerenkov((*fIdtmed)[kAl]       , kNbins, aPckov, aAbsGrid,   aQeAll, aIdxGrid);      //10 Aluminium
-  gMC->SetCerenkov((*fIdtmed)[kGlass]    , kNbins, aPckov, aAbsOpSiO2, aQeAll, aIdxOpSiO2);    //11 Glass    
+  gMC->SetCerenkov((*fIdtmed)[kAir]      , kNbins, aPckov, aAbsCH4    , aQeAll, aIdxCH4);       //1 Air
+  gMC->SetCerenkov((*fIdtmed)[kRoha]     , kNbins, aPckov, aAbsCH4    , aQeAll, aIdxCH4);       //2 Honeycomb  
+  gMC->SetCerenkov((*fIdtmed)[kSiO2]     , kNbins, aPckov, aAbsSiO2   , aQeAll, aIdxSiO2);      //3 Quartz SiO2 
+  gMC->SetCerenkov((*fIdtmed)[kC6F14]    , kNbins, aPckov, aAbsC6F14  , aQeAll, aIdxC6F14);     //4 Freon C6F14
+  gMC->SetCerenkov((*fIdtmed)[kCH4]      , kNbins, aPckov, aAbsCH4    , aQeAll, aIdxCH4);       //5 Methane CH4 
+  gMC->SetCerenkov((*fIdtmed)[kCsI]      , kNbins, aPckov, aAbsCsI    , aQeCsI, aIdxCH4);       //6 CsI
+  gMC->SetCerenkov((*fIdtmed)[kGridCu]   , kNbins, aPckov, aAbsGrid   , aQeAll, aIdxMetal);     //7 grid Cu
+  gMC->SetCerenkov((*fIdtmed)[kOpSiO2]   , kNbins, aPckov, aAbsOpSiO2 , aQeAll, aIdxMetal);     //8 Opaque quartz SiO2
+  gMC->SetCerenkov((*fIdtmed)[kGap]      , kNbins, aPckov, aAbsCH4    , aQeAll, aIdxCH4);       //9 Special methane gap
+  gMC->SetCerenkov((*fIdtmed)[kAl]       , kNbins, aPckov, aAbsGrid   , aQeAll, aIdxMetal);     //10 Aluminium
+  gMC->SetCerenkov((*fIdtmed)[kGlass]    , kNbins, aPckov, aAbsOpSiO2 , aQeAll, aIdxMetal);     //11 Glass    
+  gMC->SetCerenkov((*fIdtmed)[kGel]      , kNbins, aPckov, aAbsGel    , aQeAll, aIdxGel);       //12 Aerogel
+  gMC->SetCerenkov((*fIdtmed)[kReflector], kNbins, aPckov, aAbsRef    , aQeAll, aIdxMetal);     //13 Aerogel reflector
 }//void AliRICH::CreateMaterials()
 //__________________________________________________________________________________________________
-Float_t AliRICH::Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)const
+Float_t AliRICH::Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)
 {
 
     //ENE(EV), PDOTI=COS(INC.ANG.), PDOTR=COS(POL.PLANE ROT.ANG.)
@@ -348,38 +359,6 @@ Float_t AliRICH::Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)const
     return(fresn);
 }//Fresnel()
 //__________________________________________________________________________________________________
-Float_t AliRICH::AbsoCH4(Float_t x)const
-{
-//Evaluate the absorbtion lenght of CH4
-  Float_t sch4[9] = {.12,.16,.23,.38,.86,2.8,7.9,28.,80.};              //MB X 10^22
-  Float_t em[9] = {8.1,8.158,8.212,8.267,8.322,8.378,8.435,8.493,8.55};
-  const Float_t kLoschmidt=2.686763e19;                                      // LOSCHMIDT NUMBER IN CM-3
-  const Float_t kPressure=750.,kTemperature=283.;                                      
-  const Float_t kPn=kPressure/760.;
-  const Float_t kTn=kTemperature/273.16;
-  const Float_t kC0=-1.655279e-1;
-  const Float_t kC1=6.307392e-2;
-  const Float_t kC2=-8.011441e-3;
-  const Float_t kC3=3.392126e-4;
-               
-  Float_t crossSection=0;                        
-  if (x<7.75) 
-    crossSection=.06e-22;
-  else if(x>=7.75 && x<=8.1){                 //------ METHANE CROSS SECTION cm-2 ASTROPH. J. 214, L47 (1978)                                               
-       crossSection=(kC0+kC1*x+kC2*x*x+kC3*x*x*x)*1.e-18;
-  }else if (x> 8.1){
-    Int_t j=0;
-    while (x<=em[j] || x>=em[j+1]){
-      j++;
-      Float_t a=(sch4[j+1]-sch4[j])/(em[j+1]-em[j]);
-      crossSection=(sch4[j]+a*(x-em[j]))*1e-22;
-    }
-  }//if
-    
-    Float_t density=kLoschmidt*kPn/kTn; //CH4 molecular density 1/cm-3
-    return 1./(density*crossSection);
-}//AbsoCH4()
-//__________________________________________________________________________________________________
 void AliRICH::MakeBranch(Option_t* option)
 {
 //Create Tree branches for the RICH.
@@ -460,105 +439,12 @@ void AliRICH::Print(Option_t *option)const
   fCounters.Print();
 }//void AliRICH::Print(Option_t *option)const
 //__________________________________________________________________________________________________
-void AliRICH::CreateGeometry()
-{
-//Creates detailed geometry simulation (currently GEANT volumes tree)         
-  AliDebug(1,"Start main.");
-  Double_t cm=1,mm=0.1*cm,mkm=0.001*cm;
-  Float_t par[3];
-  Int_t id=0;
-       
-//place chambers into mother volume ALIC
-  par[0]=(6*mm+1681*mm+6*mm)/2;par[1]=(6*mm+1466*mm+6*mm)/2;par[2]=(80*mm+40*mm)*2/2;gMC->Gsvolu("RICH","BOX ",(*fIdtmed)[kCH4],par,3);//2033P1  z p84 TDR
-  if(P()->IsRadioSrc()){
-    AliInfo("Special test beam geometry");
-    gMC->Gspos("RICH",1,"ALIC",0,0,0,0, "ONLY");  //single RICH chamber without rotation test beam geometry
-  }else
-    for(int i=1;i<=kNchambers;i++){ //full RICH with 7 chambers
-      AliMatrix(id,C(i)->ThetaXd(),C(i)->PhiXd(),  C(i)->ThetaYd(),C(i)->PhiYd(),  C(i)->ThetaZd(),C(i)->PhiZd());
-      gMC->Gspos("RICH",i,"ALIC",C(i)->Center().X(),C(i)->Center().Y(),C(i)->Center().Z(),id, "ONLY");
-    }
-//Pad Panel frame  6 sectors
-  par[0]=648*mm/2;par[1]=  411*mm/2;par[2]=40  *mm/2;gMC->Gsvolu("RPPF","BOX ",(*fIdtmed)[kAl]  ,par,3);//PPF 2001P2 inner size of the slab by 1mm more
-  par[0]=181*mm/2;par[1]=89.25*mm/2;par[2]=38.3*mm/2;gMC->Gsvolu("PPFL","BOX ",(*fIdtmed)[kAir] ,par,3);//large whole
-  par[0]=114*mm/2;par[1]=89.25*mm/2;par[2]=38.3*mm/2;gMC->Gsvolu("PPFS","BOX ",(*fIdtmed)[kAir] ,par,3);//small whole
-  par[0]=644*mm/2;par[1]=  407*mm/2;par[2]= 1.7*mm/2;gMC->Gsvolu("RPC ","BOX ",(*fIdtmed)[kCsI] ,par,3);//by 0.2 mm more then actual size (PCB 2006P1)
-  
-  gMC->Gspos("RPPF",1,"RICH",    -335*mm,      -433*mm,  8*cm+20*mm,  0,"ONLY");//F1 2040P1 z p.84 TDR
-  gMC->Gspos("RPPF",2,"RICH",    +335*mm,      -433*mm,  8*cm+20*mm,  0,"ONLY");
-  gMC->Gspos("RPPF",3,"RICH",    -335*mm,         0*mm,  8*cm+20*mm,  0,"ONLY");
-  gMC->Gspos("RPPF",4,"RICH",    +335*mm,         0*mm,  8*cm+20*mm,  0,"ONLY");
-  gMC->Gspos("RPPF",5,"RICH",    -335*mm,      +433*mm,  8*cm+20*mm,  0,"ONLY");
-  gMC->Gspos("RPPF",6,"RICH",    +335*mm,      +433*mm,  8*cm+20*mm,  0,"ONLY");  
-  gMC->Gspos("RPC ",1,"RPPF",       0*mm,         0*mm,   -19.15*mm,  0,"ONLY");//PPF 2001P2 
-  gMC->Gspos("PPFL",1,"RPPF",  -224.5*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",2,"RPPF",  -224.5*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",3,"RPPF",  -224.5*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",4,"RPPF",  -224.5*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",1,"RPPF",  - 65.0*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",2,"RPPF",  - 65.0*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",3,"RPPF",  - 65.0*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",4,"RPPF",  - 65.0*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",5,"RPPF",  + 65.0*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",6,"RPPF",  + 65.0*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",7,"RPPF",  + 65.0*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFS",8,"RPPF",  + 65.0*mm,  +151.875*mm,     0.85*mm,  0,"ONLY"); 
-  gMC->Gspos("PPFL",5,"RPPF",  +224.5*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",6,"RPPF",  +224.5*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",7,"RPPF",  +224.5*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
-  gMC->Gspos("PPFL",8,"RPPF",  +224.5*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
-//Gap - anod wires 6 copies to RICH
-  par[0]=648*mm/2;par[1]=  411*mm/2 ;par[2]=4.45*mm/2;gMC->Gsvolu("RGAP ","BOX ",(*fIdtmed)[kCH4] ,par,3);//xy as PPF 2001P2 z WP 2099P1
-  par[0]=  0*mm  ;par[1]=  20*mkm/2 ;par[2]= 648*mm/2;gMC->Gsvolu("RANO","TUBE",(*fIdtmed)[kW]   ,par,3);//WP 2099P1 z = gap x PPF 2001P2
-  AliMatrix(id,180,0, 90,90, 90,0); //wires along x
-  
-  gMC->Gspos("RGAP",1,"RICH",    -335*mm,      -433*mm,8*cm-2.225*mm, 0,"ONLY"); //F1 2040P1 z WP 2099P1
-  gMC->Gspos("RGAP",2,"RICH",    +335*mm,      -433*mm,8*cm-2.225*mm, 0,"ONLY"); 
-  gMC->Gspos("RGAP",3,"RICH",    -335*mm,         0*mm,8*cm-2.225*mm, 0,"ONLY"); 
-  gMC->Gspos("RGAP",4,"RICH",    +335*mm,         0*mm,8*cm-2.225*mm, 0,"ONLY"); 
-  gMC->Gspos("RGAP",5,"RICH",    -335*mm,      +433*mm,8*cm-2.225*mm, 0,"ONLY"); 
-  gMC->Gspos("RGAP",6,"RICH",    +335*mm,      +433*mm,8*cm-2.225*mm, 0,"ONLY"); 
-  for(int i=1;i<=96;i++)
-    gMC->Gspos("RANO",i,"RGAP",     0*mm, -411/2*mm+i*4*mm, 0.185*mm, id,"ONLY"); //WP 2099P1  
-//Radiator 3 modules
-  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=  24*mm/2;  gMC->Gsvolu("RRAD","BOX ",(*fIdtmed)[kC6F14]     ,par,3); // Rad 2011P1
-  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=   4*mm/2;  gMC->Gsvolu("RRFR","BOX ",(*fIdtmed)[kRoha]      ,par,3); //front 
-  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=   5*mm/2;  gMC->Gsvolu("RRWI","BOX ",(*fIdtmed)[kSiO2]      ,par,3); //window
-  par[0]=1330*mm/2 ;par[1]=   5*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRLO","BOX ",(*fIdtmed)[kRoha]      ,par,3); //long side  
-  par[0]=  10*mm/2 ;par[1]= 403*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRSH","BOX ",(*fIdtmed)[kRoha]      ,par,3); //short side 
-  par[0]=   0      ;par[1]=  10*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRSP","TUBE",(*fIdtmed)[kSiO2]      ,par,3); //spacer        
-    
-  gMC->Gspos("RRAD",1,"RICH",   0*mm,-434*mm,   -12*mm,  0,"ONLY"); //radiator to RICH
-  gMC->Gspos("RRAD",2,"RICH",   0*mm,   0*mm,   -12*mm,  0,"ONLY"); 
-  gMC->Gspos("RRAD",3,"RICH",   0*mm,+434*mm,   -12*mm,  0,"ONLY"); 
-    gMC->Gspos("RRFR",1,"RRAD",   0*mm,   0*mm, -10.0*mm,  0,"ONLY"); //front cover 
-    gMC->Gspos("RRWI",1,"RRAD",   0*mm,   0*mm,   9.5*mm,  0,"ONLY"); //quartz window (back cover)
-    gMC->Gspos("RRLO",1,"RRAD",   0*mm,-204*mm,  -0.5*mm,  0,"ONLY"); //long side
-    gMC->Gspos("RRLO",2,"RRAD",   0*mm,+204*mm,  -0.5*mm,  0,"ONLY"); //long side
-    gMC->Gspos("RRSH",1,"RRAD",-660*mm,   0*mm,  -0.5*mm,  0,"ONLY"); //short side
-    gMC->Gspos("RRSH",2,"RRAD",+660*mm,   0*mm,  -0.5*mm,  0,"ONLY"); //short side 
-    for(int i=0;i<3;i++)
-      for(int j=0;j<10;j++)
-        gMC->Gspos("RRSP",10*i+j,"RRAD",-1330*mm/2+116*mm+j*122*mm,(i-1)*105*mm,-0.5*mm,0,"ONLY");
-
-//Sandbox  
-  par[0]=1419*mm/2 ;par[1]=1378*mm/2;par[2]=50.5*mm/2; gMC->Gsvolu("RSNB","BOX ",(*fIdtmed)[kAir]  ,par,3);  //2072P1   
-  par[0]=1419*mm/2 ;par[1]=1378*mm/2;par[2]= 0.5*mm/2; gMC->Gsvolu("RSCO","BOX ",(*fIdtmed)[kAl]   ,par,3);  //cover
-  par[0]=1359*mm/2 ;par[1]=1318*mm/2;par[2]=49.5*mm/2; gMC->Gsvolu("RSHO","BOX ",(*fIdtmed)[kRoha] ,par,3); //honeycomb structure 
-  
-  gMC->Gspos("RSNB",1,"RICH",   0*mm, 0*mm, -73.75*mm, 0,"ONLY"); //p.84 TDR sandbox to rich
-    gMC->Gspos("RSHO",1,"RSNB", 0*mm, 0*mm,      0*mm, 0,"ONLY"); //2072P1 honeycomv to sandbox
-    gMC->Gspos("RSCO",1,"RSNB", 0*mm, 0*mm,    +25*mm, 0,"ONLY"); //cover to sandbox
-    gMC->Gspos("RSCO",2,"RSNB", 0*mm, 0*mm,    -25*mm, 0,"ONLY"); //cover to sandbox
-             
-  AliDebug(1,"Stop main.");  
-}//void AliRICH::CreateGeometry()
-//__________________________________________________________________________________________________
 void AliRICH::ControlPlots()
 { 
 // Creates a set of hists to control the results of simulation. Hists are in file $HOME/RCP.root
      
-  TH1F *pHxD=0,*pHyD=0,*pNumClusH1=0,
+  TH1F             *pElecP=0 ,*pMuonP=0 ,*pPionP=0 ,*pKaonP=0 ,*pProtP=0,  //stack particles
+                   *pHxD=0,*pHyD=0,*pNumClusH1=0,
                    *pQdcH1=0,       *pSizeH1=0,
                    *pPureMipQdcH1=0,*pPureMipSizeH1=0,
                    *pPureCerQdcH1=0,*pPureCerSizeH1=0,
@@ -566,14 +452,22 @@ void AliRICH::ControlPlots()
                    *pMipQdcH1=0,    *pPhotQdcH1=0;  
   TH2F *pMapH2=0,*pPureMipMapH2=0,*pPureCerMapH2=0,*pPureFeeMapH2=0;
   
+  GetLoader()->GetRunLoader()->LoadHeader();  
+  GetLoader()->GetRunLoader()->LoadKinematics();  
+  
   Bool_t isDig =!GetLoader()->LoadDigits();
   Bool_t isClus=!GetLoader()->LoadRecPoints();
 
-  if(!isDig && !isClus){AliError("No digits and clusters! Nothing to do.");return;}
+//  if(!isDig && !isClus){AliError("No digits and clusters! Nothing to do.");return;}
   
   TStopwatch sw;TDatime time;
     
   TFile *pFile = new TFile("$(HOME)/RCP.root","RECREATE");   
+  pElecP=new TH1F("Pelec","e  versus momentum;P [GeV]",1000,-10,10); 
+  pMuonP=new TH1F("Pmuon","mu versus momentum;P [GeV]",1000,-10,10); 
+  pPionP=new TH1F("Ppion","pi versus momentum;P [GeV]",1000,-10,10); 
+  pKaonP=new TH1F("Pkaon","K  versus momentum;P [GeV]",1000,-10,10); 
+  pProtP=new TH1F("Pprot","p  versus momentum;P [GeV]",1000,-10,10); 
   
   if(isDig){
     AliInfo("Digits available");
@@ -582,36 +476,58 @@ void AliRICH::ControlPlots()
   }//isDig
   
   if(isClus){ 
-    cout<<"Clusters available\n";
+    AliInfo("Clusters available");
     pNumClusH1=new TH1F("NumClusPerEvent","Number of clusters per event;number",50,0,49);
     
-    pQdcH1        =new TH1F("ClusQdc",   "Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
+    pQdcH1        =new TH1F("ClusQdc",   "Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
     pSizeH1       =new TH1F("ClusSize",  "Cluster size all chambers;size [number of pads in cluster]",100,0,100);
-    pMapH2        =new TH2F("ClusMap",   "Cluster map;x [cm];y [cm]",1000,0,R()->P()->PcSizeX(),1000,0,R()->P()->PcSizeY());
+    pMapH2        =new TH2F("ClusMap",   "Cluster map;x [cm];y [cm]",1000,0,P()->PcSizeX(),1000,0,P()->PcSizeY());
   
-    pMipQdcH1     =new TH1F("QdcMip"      ,"MIP Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
-    pPhotQdcH1    =new TH1F("QdcPhot"     ,"Cer+Fee Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
+    pMipQdcH1     =new TH1F("QdcMip"      ,"MIP Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
+    pPhotQdcH1    =new TH1F("QdcPhot"     ,"Cer+Fee Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
         
-    pPureMipQdcH1 =new TH1F("QdcPureMip"  ,"MIP only Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
+    pPureMipQdcH1 =new TH1F("QdcPureMip"  ,"MIP only Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
     pPureMipSizeH1=new TH1F("SizePureMip" ,"MIP only Cluster size all chambers;size [number of pads in cluster]",100,0,100);
-    pPureMipMapH2 =new TH2F("MapPureMip"  ,"MIP only Cluster map;x [cm];y [cm]",1000,0,R()->P()->PcSizeX(),1000,0,R()->P()->PcSizeY());
+    pPureMipMapH2 =new TH2F("MapPureMip"  ,"MIP only Cluster map;x [cm];y [cm]",1000,0,P()->PcSizeX(),1000,0,P()->PcSizeY());
   
-    pPureCerQdcH1 =new TH1F("QdcPureCer"  ,"Cerenkov only Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
+    pPureCerQdcH1 =new TH1F("QdcPureCer"  ,"Cerenkov only Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
     pPureCerSizeH1=new TH1F("SizePureCer" ,"Cernekov only Cluster size all chambers;size [number of pads in cluster]",100,0,100);
-    pPureCerMapH2 =new TH2F("MapPureCer"  ,"Cerenkov only Cluster map;x [cm];y [cm]",1000,0,R()->P()->PcSizeX(),1000,0,R()->P()->PcSizeY());
+    pPureCerMapH2 =new TH2F("MapPureCer"  ,"Cerenkov only Cluster map;x [cm];y [cm]",1000,0,P()->PcSizeX(),1000,0,P()->PcSizeY());
     
-    pPureFeeQdcH1 =new TH1F("QdcPureFee"  ,"Feedback only Cluster Charge all chambers;q [QDC]",R()->P()->MaxQdc(),0,R()->P()->MaxQdc());
+    pPureFeeQdcH1 =new TH1F("QdcPureFee"  ,"Feedback only Cluster Charge all chambers;q [QDC]",P()->MaxQdc(),0,P()->MaxQdc());
     pPureFeeSizeH1=new TH1F("SizePureFee" ,"Feedback only Cluster size all chambers;size [number of pads in cluster]",100,0,100);
-    pPureFeeMapH2 =new TH2F("MapPureFee"  ,"Feedback only Cluster map;x [cm];y [cm]",1000,0,R()->P()->PcSizeX(),1000,0,R()->P()->PcSizeY());
+    pPureFeeMapH2 =new TH2F("MapPureFee"  ,"Feedback only Cluster map;x [cm];y [cm]",1000,0,P()->PcSizeX(),1000,0,P()->PcSizeY());
   }//isClus
   
   for(Int_t iEvtN=0;iEvtN < GetLoader()->GetRunLoader()->GetAliRun()->GetEventsPerRun();iEvtN++){//events loop
     GetLoader()->GetRunLoader()->GetEvent(iEvtN);    //gets current event
     
     if(!GetLoader()->TreeH()) GetLoader()->LoadHits();
-    for(Int_t iPrimN=0;iPrimN < GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop
-       GetLoader()->TreeH()->GetEntry(iPrimN);      
-    }
+    for(Int_t iPrimN=0;iPrimN < GetLoader()->TreeH()->GetEntries();iPrimN++){//hit tree loop
+      GetLoader()->TreeH()->GetEntry(iPrimN);      
+      for(Int_t j=0;j<Hits()->GetEntries();j++){//hits loop for a given primary
+      AliRICHhit *pHit = (AliRICHhit*)Hits()->At(j);
+        TParticle *pParticle = GetLoader()->GetRunLoader()->Stack()->Particle(pHit->GetTrack());
+        switch(pParticle->GetPdgCode()){
+          case kPositron : pElecP->Fill( pParticle->P()); break;
+          case kElectron : pElecP->Fill(-pParticle->P()); break;
+          
+          case kMuonPlus : pMuonP->Fill( pParticle->P()); break;
+          case kMuonMinus: pMuonP->Fill(-pParticle->P()); break;
+                    
+          case kPiPlus   : pPionP->Fill( pParticle->P()); break;
+          case kPiMinus  : pPionP->Fill(-pParticle->P()); break;
+          
+          case kKPlus    : pKaonP->Fill( pParticle->P()); break;
+          case kKMinus   : pKaonP->Fill(-pParticle->P()); break;
+          
+          case kProton   : pProtP->Fill( pParticle->P()); break;
+          case kProtonBar: pProtP->Fill(-pParticle->P()); break;
+              
+        }//switch PdgCode
+            
+      }//hits loop
+    }//hit tree loop
     
     if(isClus) GetLoader()->TreeR()->GetEntry(0);
     if(isDig)  GetLoader()->TreeD()->GetEntry(0);  
@@ -642,10 +558,10 @@ void AliRICH::ControlPlots()
         }//clusters loop
       }//isClus
       if(isDig){
-        for(Int_t iDigN=0;iDigN<R()->Digits(iChamN)->GetEntries();iDigN++){//digits loop
-          AliRICHdigit *pDig=(AliRICHdigit*)R()->Digits(iChamN)->At(iDigN);
+        for(Int_t iDigN=0;iDigN<Digits(iChamN)->GetEntries();iDigN++){//digits loop
+          AliRICHdigit *pDig=(AliRICHdigit*)Digits(iChamN)->At(iDigN);
           AliRICHhit   *pHit=Hit(pDig->GetTrack(0));//get first hit of this digit
-          TVector2 hitV2=R()->C(iChamN)->Mrs2Pc(pHit->OutX3()); TVector2 digV2=R()->P()->Pad2Loc(pDig->Pad());//center of pad for digit
+          TVector2 hitV2=C(iChamN)->Mrs2Pc(pHit->OutX3()); TVector2 digV2=P()->Pad2Loc(pDig->Pad());//center of pad for digit
           pHxD->Fill(hitV2.X()-digV2.X()); pHyD->Fill(hitV2.Y()-digV2.Y());
         }//digits loop
       }//isDig
@@ -656,6 +572,9 @@ void AliRICH::ControlPlots()
   if(isDig)  GetLoader()->UnloadDigits();
   if(isClus) GetLoader()->UnloadRecPoints();
   
+  GetLoader()->GetRunLoader()->UnloadHeader();  
+  GetLoader()->GetRunLoader()->UnloadKinematics();  
+  
   pFile->Write(); delete pFile;
   sw.Print();time.Print();
 }//ControlPlots()
@@ -663,15 +582,15 @@ void AliRICH::ControlPlots()
 AliRICHhit* AliRICH::Hit(Int_t tid)
 {
 //defines which hit provided by given tid for the currently loaded event
-  R()->GetLoader()->LoadHits();
-  for(Int_t iPrimN=0;iPrimN<R()->GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop      
-    R()->GetLoader()->TreeH()->GetEntry(iPrimN);
-    for(Int_t iHitN=0;iHitN<R()->Hits()->GetEntries();iHitN++){
-      AliRICHhit *pHit=(AliRICHhit*)R()->Hits()->At(iHitN);
-      if(tid==pHit->Track()) {R()->GetLoader()->UnloadHits();return pHit;}
+  GetLoader()->LoadHits();
+  for(Int_t iPrimN=0;iPrimN<GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop      
+    GetLoader()->TreeH()->GetEntry(iPrimN);
+    for(Int_t iHitN=0;iHitN<Hits()->GetEntries();iHitN++){
+      AliRICHhit *pHit=(AliRICHhit*)Hits()->At(iHitN);
+      if(tid==pHit->Track()) {GetLoader()->UnloadHits();return pHit;}
     }//hits
   }//prims loop
-  R()->GetLoader()->UnloadHits();
+  GetLoader()->UnloadHits();
   return 0;
 }
 //__________________________________________________________________________________________________
@@ -679,16 +598,16 @@ void AliRICH::PrintHits(Int_t iEvtN)
 {
 //Prints a list of RICH hits for a given event. Default is event number 0.
   AliInfo(Form("List of RICH hits for event %i",iEvtN));
-  R()->GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(R()->GetLoader()->LoadHits()) return;
+  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
+  if(GetLoader()->LoadHits()) return;
   
   Int_t iTotalHits=0;
-  for(Int_t iPrimN=0;iPrimN<R()->GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop
-    R()->GetLoader()->TreeH()->GetEntry(iPrimN);      
-    R()->Hits()->Print();
-    iTotalHits+=R()->Hits()->GetEntries();
+  for(Int_t iPrimN=0;iPrimN<GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop
+    GetLoader()->TreeH()->GetEntry(iPrimN);      
+    Hits()->Print();
+    iTotalHits+=Hits()->GetEntries();
   }
-  R()->GetLoader()->UnloadHits();
+  GetLoader()->UnloadHits();
   AliInfo(Form("totally %i hits",iTotalHits));
 }
 //__________________________________________________________________________________________________
@@ -696,29 +615,29 @@ void AliRICH::PrintSDigits(Int_t iEvtN)
 {
 //prints a list of RICH sdigits  for a given event
   Info("PrintSDigits","List of RICH sdigits for event %i",iEvtN);
-  R()->GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(R()->GetLoader()->LoadSDigits()) return;
+  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
+  if(GetLoader()->LoadSDigits()) return;
   
-  R()->GetLoader()->TreeS()->GetEntry(0);
-  R()->SDigits()->Print();
-  R()->GetLoader()->UnloadSDigits();
-  Info("PrintSDigits","totally %i sdigits",R()->SDigits()->GetEntries());
+  GetLoader()->TreeS()->GetEntry(0);
+  SDigits()->Print();
+  GetLoader()->UnloadSDigits();
+  Info("PrintSDigits","totally %i sdigits",SDigits()->GetEntries());
 }
 //__________________________________________________________________________________________________
 void AliRICH::PrintDigits(Int_t iEvtN)
 {
 //prints a list of RICH digits  for a given event
   Info("PrintDigits","List of RICH digits for event %i",iEvtN);
-  R()->GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(R()->GetLoader()->LoadDigits()) return;
+  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
+  if(GetLoader()->LoadDigits()) return;
   
   Int_t iTotalDigits=0;
-  R()->GetLoader()->TreeD()->GetEntry(0);
+  GetLoader()->TreeD()->GetEntry(0);
   for(Int_t iChamber=1;iChamber<=kNchambers;iChamber++){
-    R()->Digits(iChamber)->Print();
-    iTotalDigits+=R()->Digits(iChamber)->GetEntries();
+    Digits(iChamber)->Print();
+    iTotalDigits+=Digits(iChamber)->GetEntries();
   }
-  R()->GetLoader()->UnloadDigits();
+  GetLoader()->UnloadDigits();
   Info("PrintDigits","totally %i Digits",iTotalDigits);
 }
 //__________________________________________________________________________________________________
@@ -726,26 +645,227 @@ void AliRICH::PrintClusters(Int_t iEvtN)
 {
 //prints a list of RICH clusters  for a given event
   Info("PrintClusters","List of RICH clusters for event %i",iEvtN);
-  R()->GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(R()->GetLoader()->LoadRecPoints()) return;
+  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
+  if(GetLoader()->LoadRecPoints()) return;
   
   Int_t iTotalClusters=0;
-  R()->GetLoader()->TreeR()->GetEntry(0);
+  GetLoader()->TreeR()->GetEntry(0);
   for(Int_t iChamber=1;iChamber<=kNchambers;iChamber++){
-    R()->Clusters(iChamber)->Print();
-    iTotalClusters+=R()->Clusters(iChamber)->GetEntries();
+    Clusters(iChamber)->Print();
+    iTotalClusters+=Clusters(iChamber)->GetEntries();
   }
-  R()->GetLoader()->UnloadRecPoints();
+  GetLoader()->UnloadRecPoints();
   Info("PrintClusters","totally %i clusters",iTotalClusters);
 }
 //__________________________________________________________________________________________________
 void AliRICH::PrintTracks(Int_t iEvtN)
 {
-//prints a list of tracks (including secondaru) for a given event
-  Info("PrintTracks","List of all tracks for event %i",iEvtN);
-  R()->GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(R()->GetLoader()->GetRunLoader()->LoadHeader()) return;
-  if(R()->GetLoader()->GetRunLoader()->LoadKinematics()) return;
-  Int_t iTracksCounter=0;
-  Info("PrintTracks","totally %i tracks",iTracksCounter);
+//prints a list of tracks (including secondary) for a given event
+  AliInfo(Form("List of all tracks for event %i",iEvtN));
+  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
+  if(GetLoader()->GetRunLoader()->LoadHeader()) return;
+  if(GetLoader()->GetRunLoader()->LoadKinematics()) return;
+  AliStack *pStack=GetLoader()->GetRunLoader()->Stack();
+  
+  for(Int_t i=0;i<pStack->GetNtrack();i++){
+    pStack->Particle(i)->Print();
+  }
+  
+  AliInfo(Form("totally %i tracks including %i primaries",pStack->GetNtrack(),pStack->GetNprimary()));
+  GetLoader()->GetRunLoader()->UnloadHeader();
+  GetLoader()->GetRunLoader()->UnloadKinematics();
 }
+//__________________________________________________________________________________________________
+Int_t AliRICH::Nparticles(Int_t iPartID,Int_t iEvtN,AliRunLoader *pRL)
+{
+//counts total number of particles of given type (including secondary) for a given event
+  pRL->GetEvent(iEvtN);    
+  if(pRL->LoadHeader()) return 0;
+  if(pRL->LoadKinematics()) return 0;
+  AliStack *pStack=pRL->Stack();
+  
+  Int_t iCounter=0;
+  for(Int_t i=0;i<pStack->GetNtrack();i++){
+    if(pStack->Particle(i)->GetPdgCode()==iPartID) iCounter++;
+  }
+  
+  pRL->UnloadHeader();
+  pRL->UnloadKinematics();
+  return iCounter;
+}
+
+//__________________________________________________________________________________________________
+void AliRICH::GeomPadPanelFrame()
+{
+//Pad Panel frame  6 sectors
+  Double_t cm=1,mm=0.1*cm;//default is cm
+  Float_t par[3];
+  
+  par[0]=648*mm/2;par[1]=  411*mm/2;par[2]=40  *mm/2;gMC->Gsvolu("RPPF","BOX ",(*fIdtmed)[kAl]  ,par,3);//PPF 2001P2 inner size of the slab by 1mm more
+  par[0]=181*mm/2;par[1]=89.25*mm/2;par[2]=38.3*mm/2;gMC->Gsvolu("PPFL","BOX ",(*fIdtmed)[kAir] ,par,3);//large whole
+  par[0]=114*mm/2;par[1]=89.25*mm/2;par[2]=38.3*mm/2;gMC->Gsvolu("PPFS","BOX ",(*fIdtmed)[kAir] ,par,3);//small whole
+  par[0]=644*mm/2;par[1]=  407*mm/2;par[2]= 1.7*mm/2;gMC->Gsvolu("RPC ","BOX ",(*fIdtmed)[kCsI] ,par,3);//by 0.2 mm more then actual size (PCB 2006P1)
+  
+  gMC->Gspos("RPPF",1,"RICH",    -335*mm,      -433*mm,  8*cm+20*mm,  0,"ONLY");//F1 2040P1 z p.84 TDR
+  gMC->Gspos("RPPF",2,"RICH",    +335*mm,      -433*mm,  8*cm+20*mm,  0,"ONLY");
+  gMC->Gspos("RPPF",3,"RICH",    -335*mm,         0*mm,  8*cm+20*mm,  0,"ONLY");
+  gMC->Gspos("RPPF",4,"RICH",    +335*mm,         0*mm,  8*cm+20*mm,  0,"ONLY");
+  gMC->Gspos("RPPF",5,"RICH",    -335*mm,      +433*mm,  8*cm+20*mm,  0,"ONLY");
+  gMC->Gspos("RPPF",6,"RICH",    +335*mm,      +433*mm,  8*cm+20*mm,  0,"ONLY");  
+    gMC->Gspos("RPC ",1,"RPPF",       0*mm,         0*mm,   -19.15*mm,  0,"ONLY");//PPF 2001P2 
+    gMC->Gspos("PPFL",1,"RPPF",  -224.5*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",2,"RPPF",  -224.5*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",3,"RPPF",  -224.5*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",4,"RPPF",  -224.5*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",1,"RPPF",  - 65.0*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",2,"RPPF",  - 65.0*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",3,"RPPF",  - 65.0*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",4,"RPPF",  - 65.0*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",5,"RPPF",  + 65.0*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",6,"RPPF",  + 65.0*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",7,"RPPF",  + 65.0*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFS",8,"RPPF",  + 65.0*mm,  +151.875*mm,     0.85*mm,  0,"ONLY"); 
+    gMC->Gspos("PPFL",5,"RPPF",  +224.5*mm,  -151.875*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",6,"RPPF",  +224.5*mm,  - 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",7,"RPPF",  +224.5*mm,  + 50.625*mm,     0.85*mm,  0,"ONLY");
+    gMC->Gspos("PPFL",8,"RPPF",  +224.5*mm,  +151.875*mm,     0.85*mm,  0,"ONLY");
+}//GeomPadPanelFrame()
+//__________________________________________________________________________________________________
+void AliRICH::GeomAmpGap()
+{
+//Gap - anod wires 6 copies to RICH
+  Double_t cm=1,mm=0.1*cm,mkm=0.001*mm;//default is cm
+  Int_t matrixIdReturn=0; //matrix id returned by AliMatrix
+  Float_t par[3];
+  
+  par[0]=648*mm/2;par[1]=  411*mm/2 ;par[2]=4.45*mm/2;gMC->Gsvolu("RGAP ","BOX ",(*fIdtmed)[kCH4] ,par,3);//xy as PPF 2001P2 z WP 2099P1
+  par[0]=  0*mm  ;par[1]=  20*mkm/2 ;par[2]= 648*mm/2;gMC->Gsvolu("RANO","TUBE",(*fIdtmed)[kW]   ,par,3);//WP 2099P1 z = gap x PPF 2001P2
+  AliMatrix(matrixIdReturn,180,0, 90,90, 90,0); //wires along x
+  
+  gMC->Gspos("RGAP",1,"RICH",    -335*mm,      -433*mm,8*cm-2.225*mm, 0,"ONLY"); //F1 2040P1 z WP 2099P1
+  gMC->Gspos("RGAP",2,"RICH",    +335*mm,      -433*mm,8*cm-2.225*mm, 0,"ONLY"); 
+  gMC->Gspos("RGAP",3,"RICH",    -335*mm,         0*mm,8*cm-2.225*mm, 0,"ONLY"); 
+  gMC->Gspos("RGAP",4,"RICH",    +335*mm,         0*mm,8*cm-2.225*mm, 0,"ONLY"); 
+  gMC->Gspos("RGAP",5,"RICH",    -335*mm,      +433*mm,8*cm-2.225*mm, 0,"ONLY"); 
+  gMC->Gspos("RGAP",6,"RICH",    +335*mm,      +433*mm,8*cm-2.225*mm, 0,"ONLY"); 
+  for(int i=1;i<=96;i++)
+    gMC->Gspos("RANO",i,"RGAP",     0*mm, -411/2*mm+i*4*mm, 0.185*mm, matrixIdReturn,"ONLY"); //WP 2099P1  
+}//GeomAmpGap()
+//__________________________________________________________________________________________________
+void AliRICH::GeomRadiators()
+{
+//Defines radiators geometry  
+  Double_t mm=0.1;//default is cm
+  Float_t par[3];
+  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=  24*mm/2;  gMC->Gsvolu("RRAD","BOX ",(*fIdtmed)[kC6F14]     ,par,3); // Rad 2011P1
+  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=   4*mm/2;  gMC->Gsvolu("RRFR","BOX ",(*fIdtmed)[kRoha]      ,par,3); //front 
+  par[0]=1330*mm/2 ;par[1]= 413*mm/2  ;par[2]=   5*mm/2;  gMC->Gsvolu("RRWI","BOX ",(*fIdtmed)[kSiO2]      ,par,3); //window
+  par[0]=1330*mm/2 ;par[1]=   5*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRLO","BOX ",(*fIdtmed)[kRoha]      ,par,3); //long side  
+  par[0]=  10*mm/2 ;par[1]= 403*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRSH","BOX ",(*fIdtmed)[kRoha]      ,par,3); //short side 
+  par[0]=   0      ;par[1]=  10*mm/2  ;par[2]=  15*mm/2;  gMC->Gsvolu("RRSP","TUBE",(*fIdtmed)[kSiO2]      ,par,3); //spacer        
+    
+  gMC->Gspos("RRAD",1,"RICH",   0*mm,-434*mm,   -12*mm,  0,"ONLY"); //3 radiators to RICH
+  gMC->Gspos("RRAD",2,"RICH",   0*mm,   0*mm,   -12*mm,  0,"ONLY"); 
+  gMC->Gspos("RRAD",3,"RICH",   0*mm,+434*mm,   -12*mm,  0,"ONLY"); 
+    gMC->Gspos("RRFR",1,"RRAD",   0*mm,   0*mm, -10.0*mm,  0,"ONLY"); //front cover 
+    gMC->Gspos("RRWI",1,"RRAD",   0*mm,   0*mm,   9.5*mm,  0,"ONLY"); //quartz window (back cover)
+    gMC->Gspos("RRLO",1,"RRAD",   0*mm,-204*mm,  -0.5*mm,  0,"ONLY"); //long side
+    gMC->Gspos("RRLO",2,"RRAD",   0*mm,+204*mm,  -0.5*mm,  0,"ONLY"); //long side
+    gMC->Gspos("RRSH",1,"RRAD",-660*mm,   0*mm,  -0.5*mm,  0,"ONLY"); //short side
+    gMC->Gspos("RRSH",2,"RRAD",+660*mm,   0*mm,  -0.5*mm,  0,"ONLY"); //short side 
+    for(int i=0;i<3;i++)
+      for(int j=0;j<10;j++)
+        gMC->Gspos("RRSP",10*i+j,"RRAD",-1330*mm/2+116*mm+j*122*mm,(i-1)*105*mm,-0.5*mm,0,"ONLY");//spacers
+}//GeomRadiators()
+//__________________________________________________________________________________________________
+void AliRICH::GeomSandBox()
+{
+//Defines SandBox geometry
+  Double_t mm=0.1;//default is cm
+  Float_t par[3];
+  par[0]=1419*mm/2 ;par[1]=1378*mm/2;par[2]=50.5*mm/2; gMC->Gsvolu("RSNB","BOX ",(*fIdtmed)[kAir]  ,par,3);  //2072P1   
+  par[0]=1419*mm/2 ;par[1]=1378*mm/2;par[2]= 0.5*mm/2; gMC->Gsvolu("RSCO","BOX ",(*fIdtmed)[kAl]   ,par,3);  //cover
+  par[0]=1359*mm/2 ;par[1]=1318*mm/2;par[2]=49.5*mm/2; gMC->Gsvolu("RSHO","BOX ",(*fIdtmed)[kRoha] ,par,3); //honeycomb structure 
+  
+  gMC->Gspos("RSNB",1,"RICH",   0*mm, 0*mm, -73.75*mm, 0,"ONLY"); //p.84 TDR sandbox to rich
+    gMC->Gspos("RSHO",1,"RSNB", 0*mm, 0*mm,      0*mm, 0,"ONLY"); //2072P1 honeycomv to sandbox
+    gMC->Gspos("RSCO",1,"RSNB", 0*mm, 0*mm,    +25*mm, 0,"ONLY"); //cover to sandbox
+    gMC->Gspos("RSCO",2,"RSNB", 0*mm, 0*mm,    -25*mm, 0,"ONLY"); //cover to sandbox
+}//GeomSandBox()
+//__________________________________________________________________________________________________
+void AliRICH::GeomRadioSrc()
+{
+// Defines geometry for radioactive source  
+  Double_t cm=1,mm=0.1*cm,mkm=0.001*cm;
+  Float_t par[3];
+  
+  par[0]=0 ;par[1]= 70*mm/2  ;par[2]=  30*mm/2;      gMC->Gsvolu("RSRC","TUBE",(*fIdtmed)[kCH4]    ,par,3); //top src container
+    par[0]=0 ;par[1]= 38*mm/2  ;par[2]=  21.8*mm/2;  gMC->Gsvolu("RSAG","TUBE",(*fIdtmed)[kAl]     ,par,3); //Al glass
+      par[0]=0 ;par[1]= 34*mm/2  ;par[2]=  20*mm/2;  gMC->Gsvolu("RSPP","TUBE",(*fIdtmed)[kPerpex] ,par,3); //perpex plug
+        par[0]=0 ;par[1]= 5*mm/2  ;par[2]=  15*mm/2; gMC->Gsvolu("RSSC","TUBE",(*fIdtmed)[kSteel]  ,par,3); //steel screw in center of perpex
+        par[0]=0 ;par[1]= 2*mm/2  ;par[2]=  10*mm/2; gMC->Gsvolu("RSSS","TUBE",(*fIdtmed)[kSteel]  ,par,3); //Steel screw to support Sr90 
+          par[0]=0 ;par[1]= 1*mm/2  ;par[2]= 1*mm/2; gMC->Gsvolu("RSSR","TUBE",(*fIdtmed)[kSr90]   ,par,3); //Sr90 source
+        par[0]=0 ;par[1]= 4*mm/2  ;par[2]= 10*mm/2;  gMC->Gsvolu("RSWP","TUBE",(*fIdtmed)[kAir]    ,par,3); //Air hole in perpex plug      
+      par[0]=0 ;par[1]= 5*mm/2  ;par[2]= 1.8*mm/2;   gMC->Gsvolu("RSWA","TUBE",(*fIdtmed)[kAir]    ,par,3); //Air hole in Al glass bottom
+    par[0]=0 ;par[1]= 30*mm/2  ;par[2]= 50*mkm/2;    gMC->Gsvolu("RSMF","TUBE",(*fIdtmed)[kMylar]  ,par,3); //Mylar foil                
+    
+  gMC->Gspos("RSRC",1,"RICH",       30*cm,        0,     1*cm, 0,"ONLY"); //source to RICH
+    gMC->Gspos("RSMF",1,"RSRC",         0,        0,21.8*mm/2+50*mkm/2, 0,"ONLY");//mylar foil to top src volume
+    gMC->Gspos("RSAG",1,"RSRC",         0,        0,        0, 0,"ONLY");//Al glass to fake Src volume 
+      gMC->Gspos("RSWA",1,"RSAG",    6*mm,        0,   -10*mm, 0,"ONLY");//air whole in al glass bottom
+      gMC->Gspos("RSPP",1,"RSAG",       0,        0,   0.9*mm, 0,"ONLY");//perpex plug to Al glass
+        gMC->Gspos("RSWP",1,"RSPP",  6*mm,        0,    -5*mm, 0,"ONLY");//air whole in perpex plug
+        gMC->Gspos("RSSC",1,"RSPP",     0,        0,   2.5*mm, 0,"ONLY");//steel screw in center of perpex plug
+        gMC->Gspos("RSSS",1,"RSPP",  6*mm,        0,     5*mm, 0,"ONLY");//steel screw to support Sr90  in perpex plug
+          gMC->Gspos("RSSR",1,"RSSS",   0,        0,  -4.5*mm, 0,"ONLY");//Sr90  in support steel screw
+}//GeomSr90()
+//__________________________________________________________________________________________________
+void AliRICH::GeomAerogel()
+{
+//Creates detailed geometry for aerogel study.
+  AliDebug(1,"Start.");
+  Double_t cm=1;
+  Float_t par[3]; //tmp array for volume dimentions
+       
+  par[0]=10.1*cm/2;par[1]=10.1*cm/2;par[2]=10.1*cm/2;
+  gMC->Gsvolu("RREF","BOX ",(*fIdtmed)[kReflector],par,3);//reflector box
+  gMC->Gspos("RREF",1,"RICH",0,0,0,0, "ONLY");            //put it to RICH volume
+  
+  par[0]=10*cm/2;par[1]=10*cm/2;par[2]=10*cm/2;
+  gMC->Gsvolu("RGEL","BOX ",(*fIdtmed)[kGel],par,3);//10x10x10 cm^3 cubic of aerogel
+  gMC->Gspos("RGEL",1,"RREF",0,0,0,0,"ONLY");//put gel cell to reflector
+  AliDebug(1,"Stop.");  
+}//GeomAerogel()
+//__________________________________________________________________________________________________
+void AliRICH::CreateGeometry()
+{
+//Creates detailed geometry simulation (currently GEANT volumes tree)         
+  AliDebug(1,"Start main.");
+  Double_t mm=0.1;//default is cm
+  Float_t par[3];
+  Int_t matrixIdReturn=0; //matrix id returned by AliMatrix
+       
+//place chambers into mother volume ALIC
+  par[0]=(6*mm+1681*mm+6*mm)/2;par[1]=(6*mm+1466*mm+6*mm)/2;par[2]=(80*mm+40*mm)*2/2;
+  gMC->Gsvolu("RICH","BOX ",(*fIdtmed)[kCH4],par,3);//2033P1  z p84 TDR
+  for(int i=1;i<=P()->Nchambers();i++){ //test configuration with single chamber is taken into account automaticaly in AliRICHParam
+    AliMatrix(matrixIdReturn,
+                   C(i)->ThetaXd(),C(i)->PhiXd(),  
+                   C(i)->ThetaYd(),C(i)->PhiYd(),  
+                   C(i)->ThetaZd(),C(i)->PhiZd());
+    gMC->Gspos("RICH",i,"ALIC",C(i)->Center().X(),
+                               C(i)->Center().Y(),
+                               C(i)->Center().Z(),matrixIdReturn, "ONLY");
+  }
+  
+  if(P()->IsAerogel()) 
+    GeomAerogel();
+  else{
+    GeomPadPanelFrame();
+    GeomAmpGap();
+    if(P()->IsRadioSrc())    GeomRadioSrc(); else GeomRadiators(); 
+    GeomSandBox();           
+  }
+  AliDebug(1,"Stop main.");  
+}//CreateGeometry()
+//__________________________________________________________________________________________________
index 8ef6568..bd42271 100644 (file)
@@ -75,7 +75,7 @@ protected:
 class AliRICHcluster :public TObject
 {
 public:
-  enum ClusterStatus {kEdge,kShape,kSize,kRaw,kResolved,kEmpty=kBad};
+  enum ClusterStatus {kEdge,kShape,kSize,kRaw,kResolved,kEmpty}; 
                     AliRICHcluster():TObject(),fCFM(0),fSize(0),fShape(0),fQdc(0),fChamber(0),fX(0),fY(0),fStatus(kEmpty),fDigits(0) {}    
   virtual          ~AliRICHcluster()                 {AliDebug(1,"Start");/*Reset();*/}  
          void       Reset()                          {DeleteDigits();fCFM=fSize=fShape=fQdc=fChamber=0;fX=fY=0;fStatus=kEmpty;} //cleans the cluster
@@ -167,6 +167,7 @@ class AliESD;
 class AliRICH : public AliDetector 
 {
 public:
+//ctor & dtor    
             AliRICH();                                            
             AliRICH(const char *name, const char *title);
             AliRICH(const AliRICH& RICH):AliDetector(RICH) {;}  //copy ctor 
@@ -178,16 +179,19 @@ public:
   virtual void          StepManager()                               =0;                                  //interface from AliMC
   virtual void          Hits2SDigits();                                                                  //interface from AliSimulation
   virtual AliDigitizer* CreateDigitizer(AliRunDigitizer* man) const {return new AliRICHDigitizer(man);}  //interface from AliSimulation
-//  virtual void          Reconstruct()                         const;                                     //interface from AliReconstruction
-//  virtual void          FillESD(AliESD *pESD)                 const;                                     //interface from AliReconstruction          
   virtual void          SetTreeAddress();                                                                //interface from AliLoader
   virtual void          MakeBranch(Option_t *opt=" ");                                                   //interface from AliLoader
   virtual void          CreateMaterials();                                                               //interface from AliMC
   virtual void          CreateGeometry();                                                                //interface from AliMC
+          void          GeomPadPanelFrame();                                                             //defines PPF geometry
+          void          GeomAmpGap();                                                                    //defines gap geometry + anod wires
+          void          GeomRadiators();                                                                 //defines radiators geometry
+          void          GeomSandBox();                                                                   //defines sandbox geometry
+          void          GeomRadioSrc();                                                                  //defines radio source geometry
+          void          GeomAerogel();                                                                   //defines aerogel geometry
   virtual void          BuildGeometry();                                                                 //interface 
 //private part  
-          Float_t AbsoCH4(Float_t x)const;                               //calculates absorption length for methane
-          Float_t Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)const;  //deals with Fresnel absorption
+  static  Float_t Fresnel(Float_t ene,Float_t pdoti, Bool_t pola);       //deals with Fresnel absorption
   inline  void    CreateHits();                                          //create hits container as a simple list
   inline  void    CreateSDigits();                                       //create sdigits container as a simple list
   inline  void    CreateDigits();                                        //create digits container as 7 lists, one per chamber
@@ -199,42 +203,44 @@ public:
   TClonesArray*   SDigits()             const{return fSdigits;}
   TClonesArray*   Digits(Int_t iC)      const{if(fDigitsNew) return (TClonesArray *)fDigitsNew->At(iC-1);else return 0;}
   TClonesArray*   Clusters(Int_t iC)    const{if(fClusters)  return (TClonesArray *)fClusters->At(iC-1);else return 0;}
-  AliRICHChamber* C(Int_t iC)           const{return fpParam->C(iC);}                       //provides pointer to a given chamber
-  AliRICHParam*   P()                   const{return fpParam;}                              //provides pointer to a RICH params
-  AliRICH*        R()                        {return this;}                                 //provides pointer to RICH main object
-  TVector         Counters()            const{return fCounters;}                            //provides a set of counters
-  void            ControlPlots();                                                           //utility
-  virtual void    Print(Option_t *option="")               const;                           //prints current RICH status
-  void            PrintHits    (Int_t iEvent=0);                                            //utility
-  void            PrintSDigits (Int_t iEvent=0);                                            //utility
-  void            PrintDigits  (Int_t iEvent=0);                                            //utility
-  void            PrintClusters(Int_t iEvent=0);                                            //utility
-  void            PrintTracks  (Int_t iEvent=0);                                            //utility
+  AliRICHChamber* C(Int_t iC)           const{return fpParam->C(iC);}   //provides pointer to a given chamber
+  AliRICHParam*   P()                   const{return fpParam;}          //provides pointer to a RICH params
+  AliRICH*        R()                        {return this;}             //provides pointer to RICH main object
+  TVector         Counters()            const{return fCounters;}        //provides a set of counters
+  void            ControlPlots();                                       //creates ~/RCP.root with a set of QA plots
+  virtual void    Print(Option_t *option="")               const;       //prints current RICH status
+  void            PrintHits    (Int_t iEvent=0);                        //prints a list of hits for a given event
+  void            PrintSDigits (Int_t iEvent=0);                        //prints a list of sdigits for a given event
+  void            PrintDigits  (Int_t iEvent=0);                        //prints a list of digits for a given event
+  void            PrintClusters(Int_t iEvent=0);                        //prints a list of clusters for a given event
+  void            PrintTracks  (Int_t iEvent=0);                        //prints a list of tracks for a given event
+  static Int_t    Nparticles   (Int_t iPID,Int_t iEvent=0,AliRunLoader *p=0); //returns a number of particles 
+  AliRICHhit*     Hit(Int_t tid);                                       //returns pointer of the first RICH hit created by a given particle 
             
-  void AddHit(Int_t c,Int_t tid,TVector3 i3,TVector3 o3,Double_t eloss=0){TClonesArray &tmp=*fHits;new(tmp[fNhits++])AliRICHhit(c,tid,i3,o3,eloss);}
-  inline void AddSDigit(Int_t c,TVector pad,Double_t q,Int_t pid,Int_t tid); 
-  void AddDigit(int c,TVector pad,int q,int cfm,int *tid)//Add simulated digit
-       {TClonesArray &tmp=*((TClonesArray*)fDigitsNew->At(c-1));new(tmp[fNdigitsNew[c-1]++])AliRICHdigit(c,pad,q,cfm,tid[0],tid[1],tid[2]);}  
-  void AddDigit(Int_t c,TVector pad,Int_t q)//for real data digits
-       {TClonesArray &tmp=*((TClonesArray*)fDigitsNew->At(0));new(tmp[fNdigitsNew[0]++])AliRICHdigit(c,pad,q,0,-1,-1,-1);}  
-  void AddCluster(AliRICHcluster &cl)                     
-       {Int_t c=cl.C()-1;TClonesArray &tmp=*((TClonesArray*)fClusters->At(c));new(tmp[fNclusters[c]++])AliRICHcluster(cl);}
-  AliRICHhit* Hit(Int_t tid);           //returns pointer ot RICH hit for a given tid
+  inline void AddHit(Int_t chamber,Int_t tid,TVector3 in3,TVector3 out3,Double_t eloss=0);           //add new hit
+  inline void AddSDigit(Int_t c,TVector pad,Double_t q,Int_t pid,Int_t tid);                         //add new sdigit
+  inline void AddDigit(int c,TVector pad,int q,int cfm,int *tid);                                    //add new digit
+//  void AddDigit(Int_t c,TVector pad,Int_t q)//for real data digits
+//       {TClonesArray &tmp=*((TClonesArray*)fDigitsNew->At(0));new(tmp[fNdigitsNew[0]++])AliRICHdigit(c,pad,q,0,-1,-1,-1);}  
+  inline void AddCluster(AliRICHcluster &cl);                                                        //add new cluster                        
 protected:  
-  enum                  {kAir=1,kRoha,kSiO2,kC6F14,kCH4,kCsI,kGridCu,kOpSiO2,kGap,kAl,kGlass,kCu,kW,kSteel,kPerpex,kSr90};
-  AliRICHParam         *fpParam;             //main RICH parametrization     
-                                             //fHits and fDigits belong to AliDetector
-  TClonesArray         *fSdigits;            //! list of sdigits  
-  Int_t                 fNsdigits;           //! current number of sdigits
+  enum                  EMedia {kAir=1,kRoha,kSiO2,kC6F14,kCH4,kCsI,kGridCu,kOpSiO2,kGap,kAl,kGlass,kCu,kW
+                                      ,kSteel,kPerpex,kSr90,kMylar,kGel,kReflector};
+  enum                  ECounters {kStepManager=0,kCerProdTot,kCerProdRad,kCerKillTot,kCerKillRad,kCerKillRef,kEleProdTot};
+  AliRICHParam         *fpParam;                   //main RICH parametrization     
+                                                   //fHits and fDigits belong to AliDetector
+  TClonesArray         *fSdigits;                  //! list of sdigits  
+  Int_t                 fNsdigits;                 //! current number of sdigits
   
-  TObjArray            *fDigitsNew;          //! each chamber holds it's one lists of digits
-  Int_t                 fNdigitsNew[kNchambers];   //! array of current numbers of digits
+  TObjArray            *fDigitsNew;                //! each chamber holds it's one lists of digits
+  Int_t                 fNdigitsNew[7];            //! array of current numbers of digits
   
-  TObjArray            *fClusters;           //! each chamber holds it's one lists of clusters 
-  Int_t                 fNclusters[kNchambers];    //! array of current numbers of raw clusters
+  TObjArray            *fClusters;                 //! each chamber holds it's one lists of clusters 
+  Int_t                 fNclusters[7];             //! array of current numbers of raw clusters
   
-  TVector               fCounters;           //Photon history conters, explanation in StepManager() 
-  ClassDef(AliRICH,7)                        //Main RICH class 
+  TVector               fCounters;                 //Particle history counters, explanation in StepManager() 
+  
+  ClassDef(AliRICH,7)                              //Main RICH class 
 };//class AliRICH  
 
 //__________________________________________________________________________________________________
@@ -268,6 +274,13 @@ void AliRICH::CreateClusters()
   for(Int_t i=0;i<kNchambers;i++) {fClusters->AddAt(new TClonesArray("AliRICHcluster",10000), i); fNclusters[i]=0;}
 }
 //__________________________________________________________________________________________________
+void AliRICH::AddHit(Int_t c,Int_t tid,TVector3 i3,TVector3 o3,Double_t eloss)
+{
+//add new RICH hit to the list of hits  
+  TClonesArray &tmp=*fHits;
+  new(tmp[fNhits++])AliRICHhit(c,tid,i3,o3,eloss);
+}//AddHit()
+//__________________________________________________________________________________________________
 void AliRICH::AddSDigit(Int_t c,TVector pad,Double_t q,Int_t pid,Int_t tid) 
 { 
   Int_t cfm;  
@@ -279,4 +292,15 @@ void AliRICH::AddSDigit(Int_t c,TVector pad,Double_t q,Int_t pid,Int_t tid)
   TClonesArray &tmp=*fSdigits; new(tmp[fNsdigits++])AliRICHdigit(c,pad,q,cfm,tid,kBad,kBad);
 }
 //__________________________________________________________________________________________________
+void AliRICH::AddDigit(int c,TVector pad,int q,int cfm,int *tid)
+{
+  TClonesArray &tmp=*((TClonesArray*)fDigitsNew->At(c-1));
+  new(tmp[fNdigitsNew[c-1]++])AliRICHdigit(c,pad,q,cfm,tid[0],tid[1],tid[2]);
+}    
+//__________________________________________________________________________________________________
+void AliRICH::AddCluster(AliRICHcluster &cl)                     
+{
+  Int_t c=cl.C()-1;TClonesArray &tmp=*((TClonesArray*)fClusters->At(c));
+  new(tmp[fNclusters[c]++])AliRICHcluster(cl);
+}
 #endif//#ifndef AliRICH_h
index ee0ccc5..4221e04 100644 (file)
@@ -16,6 +16,7 @@
 #include "AliRICHChamber.h"
 #include <TRotMatrix.h>
 #include <AliLog.h>
+#include "AliRICHParam.h"
 
 ClassImp(AliRICHChamber)       
 //______________________________________________________________________________
@@ -26,39 +27,42 @@ AliRICHChamber::AliRICHChamber(Int_t iChamber):TNamed()
 // 5 4 3   |      
 //   2 1   |
 // <-----z y   . x     
-//  horizontal angle between chambers  19.5 grad
-//  vertical angle between chambers    20   grad     
-  RotY(90);//rotate around y
-  fCenterX3.SetXYZ(490,0,0);fPcX3.SetXYZ(490+8-0.2,0,0);fRadX3.SetXYZ(490-2,0,0); //shift center along x by 490 cm
+  if(iChamber!=0){ //if iChamber=0 then special test config: chamber without any shift and rotation
+    const Double_t kAngHor=19.5; //  horizontal angle between chambers  19.5 grad
+    const Double_t kAngVer=20;   //  vertical angle between chambers    20   grad     
+    const Double_t kAngCom=30;   //  common RICH rotation with respect to x axis  20   grad     
   
-  switch(iChamber){
-    case 0:                    //special test beam configuration without rotation.
-      break;
-    case 1:        
-      RotY(19.5); RotZ(-20);   //right and down 
-      break;      
-    case 2:
-      RotZ(-20);              //down
-      break;      
-    case 3:
-      RotY(19.5);             //right 
-      break;      
-    case 4:          
-      break;                  //no rotation
-    case 5:
-      RotY(-19.5);            //left   
-      break;      
-    case 6:
-      RotZ(20);               //up
-      break;      
-    case 7:
-      RotY(-19.5); RotZ(20);  //left and up 
-      break;      
-    default:
-      Fatal("named ctor","Wrong chamber number %i, check CreateChamber ctor",iChamber);
-  }//switch(iModuleN)
+    RotY(90);//rotate around y since initila position is in XY plane
+    fCenterX3.SetXYZ(490,0,0);  //shift center along x by 490 cm
+    fPcX3.SetXYZ(490+8-0.2,0,0);//position of the center of apmlification gap (anod wires plane)
+    fRadX3.SetXYZ(490-2,0,0);   //position of the entrance to freon
+    switch(iChamber){
+      case 1:        
+        RotY(kAngHor);  RotZ(-kAngVer);   //right and down 
+        break;      
+      case 2:
+                        RotZ(-kAngVer);   //down
+        break;      
+      case 3:
+        RotY(kAngHor);                   //right 
+        break;      
+      case 4:          
+        break;                           //no rotation
+      case 5:
+        RotY(-kAngHor);                  //left   
+        break;      
+      case 6:
+                        RotZ(kAngVer);     //up
+        break;      
+      case 7:
+        RotY(-kAngHor); RotZ(kAngVer);  //left and up 
+        break;      
+      default:
+        Fatal("named ctor","Wrong chamber number %i, check CreateChamber ctor",iChamber);
+    }//switch(iChamber)
+    RotZ(kAngCom);     //apply common rotation  
+  }//if(iChamber
   fName=Form("RICHc%i",iChamber);fTitle=Form("RICH chamber %i",iChamber);
-  RotZ(30);     //apply common rotation  
   fpRotMatrix=new TRotMatrix("rot"+fName,"rot"+fName, Rot().ThetaX()*TMath::RadToDeg(), Rot().PhiX()*TMath::RadToDeg(),
                                                       Rot().ThetaY()*TMath::RadToDeg(), Rot().PhiY()*TMath::RadToDeg(),
                                                       Rot().ThetaZ()*TMath::RadToDeg(), Rot().PhiZ()*TMath::RadToDeg());
@@ -67,6 +71,7 @@ AliRICHChamber::AliRICHChamber(Int_t iChamber):TNamed()
 void AliRICHChamber::Print(Option_t *opt) const
 {
 // Debug printout
-  ToAliInfo(fCenterX3.Print(opt));
+  TNamed::Print(opt);
+  fCenterX3.Print(opt);
 }//Print()
 //__________________________________________________________________________________________________
index 0271244..4145510 100644 (file)
 //  **************************************************************************
 #include "AliRICHParam.h"
 #include "AliRICHChamber.h"
+#include "AliRICHDisplFast.h"
+#include <TCanvas.h>
+#include <TLatex.h>
+#include <THStack.h>
+#include <TLegend.h>
+#include <TView.h>
+#include <TPolyMarker3D.h>
+#include <TPolyLine3D.h>
 
 ClassImp(AliRICHParam)
-Bool_t   AliRICHParam::fgIsWireSag            =kTRUE;
-Bool_t   AliRICHParam::fgIsResolveClusters    =kTRUE;
-Bool_t   AliRICHParam::fgIsRadioSrc           =kFALSE;
-Double_t AliRICHParam::fgAngleRot             =-60;
+Bool_t   AliRICHParam::fgIsWireSag            =kTRUE;   //take ware sagita into account?
+Bool_t   AliRICHParam::fgIsResolveClusters    =kTRUE;   //do cluster resolving?
+Bool_t   AliRICHParam::fgIsRadioSrc           =kFALSE;  //put radioactive source instead of radiators?
+Bool_t   AliRICHParam::fgIsAerogel            =kFALSE;  //special aerogel configuration
+Bool_t   AliRICHParam::fgIsTestBeam           =kFALSE;  //special test beam configuration
+
 Int_t    AliRICHParam::fgHV[kNsectors]        ={2050,2050,2050,2050,2050,2050};
 Int_t    AliRICHParam::fgNsigmaTh             =4;
 Float_t  AliRICHParam::fgSigmaThMean          =1.132; //QDC 
@@ -29,14 +39,16 @@ Float_t  AliRICHParam::fgSigmaThSpread        =0.035; //
 void AliRICHParam::Print(Option_t*)
 {
   AliInfo(Form("Pads in chamber (%3i,%3i) in sector (%2i,%2i)",NpadsX(),NpadsY(),NpadsXsec(),NpadsYsec()));
-  ToAliInfo(fpChambers->Print());
+  AliInfo(Form("Resolve clusters %i sagita %i Radio source %i Aerogel %i TestBeam %i",
+                IsResolveClusters(),IsWireSag(),IsRadioSrc(),IsAerogel(),IsTestBeam())); 
+  fpChambers->Print();
 }//Print()
 //__________________________________________________________________________________________________
 void AliRICHParam::CreateChambers()
 {
 //Create all RICH Chambers on each call. Previous chambers deleted.
   if(fpChambers) delete fpChambers;
-  if(IsRadioSrc()){ 
+  if(fgIsTestBeam){ 
     fpChambers=new TObjArray(1);//test beam configuration 1 chamber
     fpChambers->AddAt(new AliRICHChamber(0),0);  
   }else{ 
@@ -45,3 +57,196 @@ void AliRICHParam::CreateChambers()
   }
   fpChambers->SetOwner();
 }//CreateChambers()
+//__________________________________________________________________________________________________
+Float_t AliRICHParam::AbsCH4(Float_t eV)
+{
+//Evaluate the absorbtion lenght of CH4 for a photon of energy eV in electron-volts
+  const Float_t kLoschmidt=2.686763e19;                                      // LOSCHMIDT NUMBER IN CM-3
+  const Float_t kPressure=750.0;          //mm of Hg
+  const Float_t kTemperature=283.0;       //K (10 grad C)                               
+  const Float_t kPn=kPressure/760.;
+  const Float_t kTn=kTemperature/273.16;
+  const Float_t kC0=-1.655279e-1;
+  const Float_t kC1= 6.307392e-2;
+  const Float_t kC2=-8.011441e-3;
+  const Float_t kC3= 3.392126e-4;
+               
+  Float_t crossSection=0;                        
+  if (eV<7.75) 
+    crossSection=0.06e-22;
+  else                 //------ METHANE CROSS SECTION cm-2 ASTROPH. J. 214, L47 (1978)                                               
+    crossSection=(kC0+kC1*eV+kC2*eV*eV+kC3*eV*eV*eV)*1.e-18;
+    
+    Float_t density=kLoschmidt*kPn/kTn; //CH4 molecular concentration (cm^-3)
+    return 1.0/(density*crossSection);
+}//AbsoCH4()
+//__________________________________________________________________________________________________
+void AliRICHParam::TestSeg()
+{
+  
+  new TCanvas("name","PC segmentation");
+  gPad->Range(-20,-20,200,150);
+  AliRICHDisplFast::DrawSectors();
+  
+  TLatex t; t.SetTextSize(0.02);
+  t.DrawText(0,140,"View from interaction point");
+  t.DrawLatex(PcSizeX()+10,120,Form("Pc  %6.2fx%6.2fcm %3ix%3ipads",PcSizeX()    ,PcSizeY(),    NpadsX()   ,NpadsY()));
+  t.DrawLatex(PcSizeX()+10,115,Form("Sec %6.2fx%5.2fcm %3ix%2ipads",SectorSizeX(),SectorSizeY(),NpadsXsec(),NpadsYsec()));
+  t.DrawLatex(PcSizeX()+10,110,Form("Pad %6.2fx%4.2fcm DeadZone %6.2fcm",PadSizeX()   ,PadSizeY(),DeadZone()));
+  
+  TVector2 v2;
+  t.SetTextAlign(12);
+  v2=AliRICHParam::Pad2Loc( 1,24);  t.DrawText(v2.X(),v2.Y(),"sec 1");  
+  v2=AliRICHParam::Pad2Loc(81,24);  t.DrawText(v2.X(),v2.Y(),"sec 2");  
+  v2=AliRICHParam::Pad2Loc( 1,70);  t.DrawText(v2.X(),v2.Y(),"sec 3");  
+  v2=AliRICHParam::Pad2Loc(81,70);  t.DrawText(v2.X(),v2.Y(),"sec 4");  
+  v2=AliRICHParam::Pad2Loc( 1,120); t.DrawText(v2.X(),v2.Y(),"sec 5");  
+  v2=AliRICHParam::Pad2Loc(81,120); t.DrawText(v2.X(),v2.Y(),"sec 6");  
+  
+//  TGaxis *pAx=new TGaxis(0,0,140,  0,0,140,510,"-="); pAx->SetTitle("x, cm"); pAx->SetTextSize(0.05); pAx->Draw();
+//  TGaxis *pAy=new TGaxis(0,0,  0,140,0,140,510,"-="); pAy->SetTitle("y, cm"); pAy->SetTextSize(0.05); pAy->Draw();
+   
+  
+  t.SetTextColor(kBlue);  
+
+  Double_t margin=5;
+  Double_t x0=0; Double_t x1=SectorSizeX(); Double_t x2=SectorSizeX()+DeadZone(); Double_t x3=PcSizeX();
+  Double_t y0=0; Double_t y1=SectorSizeY(); Double_t y2=SectorSizeY()+DeadZone(); 
+  Double_t y3=2*SectorSizeY()+DeadZone(); Double_t y4=PcSizeY()-SectorSizeY();
+  Double_t y5=PcSizeY();
+   
+//write pad numbers along x  
+  t.SetTextAlign(13); t.DrawText(x0,y0,"1");  
+  t.SetTextAlign(33); t.DrawText(x1,y0,"80");
+  t.SetTextAlign(13); t.DrawText(x2,y0,"81");   
+  t.SetTextAlign(33); t.DrawText(x3,y0,"160");
+//write pad numbers along y    
+  t.SetTextAlign(31); t.DrawText(x0,y0,"1");  
+  t.SetTextAlign(33); t.DrawText(x0,y1,"48"); 
+  t.SetTextAlign(31); t.DrawText(x0,y2,"49");
+  t.SetTextAlign(33); t.DrawText(x0,y3,"96"); 
+  t.SetTextAlign(31); t.DrawText(x0,y4,"97");
+  t.SetTextAlign(33); t.DrawText(x0,y5,"144");        
+  
+   
+//positions along x   
+  t.SetTextColor(kGreen);
+  t.SetTextAlign(11);t.DrawText(x0,y0-margin,Form("%5.2f",x0));   
+  t.SetTextAlign(31);t.DrawText(x1,y0-margin,Form("%5.2f",x1));   
+  t.SetTextAlign(11);t.DrawText(x2,y0-margin,Form("%5.2f",x2));   
+  t.SetTextAlign(31);t.DrawText(x3,y0-margin,Form("%5.2f",x3));   
+//positions along y
+  t.SetTextAlign(31);t.DrawText(x0-margin,y0,Form("%5.2f",y0));   
+  t.SetTextAlign(33);t.DrawText(x0-margin,y1,Form("%5.2f",y1));   
+  t.SetTextAlign(31);t.DrawText(x0-margin,y2,Form("%5.2f",y2));   
+  t.SetTextAlign(33);t.DrawText(x0-margin,y3,Form("%5.2f",y3));   
+  t.SetTextAlign(31);t.DrawText(x0-margin,y4,Form("%5.2f",y4));   
+  t.SetTextAlign(33);t.DrawText(x0-margin,y5,Form("%5.2f",y5));   
+//coners      
+  t.SetTextColor(kRed);
+  t.SetTextSize(0.01);
+  TVector pad(2);
+//sector 1  
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y0)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y1)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y1)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y0)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+//secr 3
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y2)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y3)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y3)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y2)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+//secr 5   
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y4)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x0,y5)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y5)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x1,y4)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y4)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y5)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y5)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y4)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y2)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y3)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y3)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y2)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y0)));t.SetTextAlign(11);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x2,y1)));t.SetTextAlign(13);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y1)));t.SetTextAlign(33);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+  v2=Pad2Loc(Loc2Pad(TVector2(x3,y0)));t.SetTextAlign(31);t.DrawText(v2.X(),v2.Y(),Form("%5.2f,%5.2f",v2.X(),v2.Y()));
+}//TestSeg()
+//__________________________________________________________________________________________________
+void AliRICHParam::TestResp()
+{
+//test the response set of methodes  
+  TCanvas *pC=new TCanvas("c","Amplification test",900,800);
+  pC->Divide(1,2);
+  
+  
+  const Int_t nPoints=8;
+  THStack *pStackPhot=new THStack("StackPhot","photons");
+  THStack *pStackMip =new THStack("StackMip","mips");
+  TLegend *pLeg=new TLegend(0.6,0.2,0.9,0.5,"legend");    
+  TH1F *apHphot[nPoints];
+  TH1F *apHmip[nPoints];
+  
+  Double_t starty=0;
+  Double_t deltay=AliRICHParam::SectorSizeY()/nPoints;
+  
+  for(int i=0;i<nPoints;i++){
+    apHphot[i]=new TH1F(Form("hphot%i",i),"Qdc for Photon;QDC;Counts",500,0,500); apHphot[i]->SetLineColor(i);pStackPhot->Add(apHphot[i]);                 
+    apHmip[i] =new TH1F(Form("hmip%i",i),"Qdc for Mip;QDC;Counts",4000,0,4000);   apHmip[i]->SetLineColor(i);pStackMip->Add(apHmip[i]);                 
+    
+    pLeg->AddEntry(apHphot[i],Form("@(10,%5.2f->%5.2f)",starty+i*deltay,starty+i*deltay-SectorSizeY()/2));
+  }
+        
+  
+  TVector2 x2(0,0);  
+  for(Int_t i=0;i<10000;i++){//events loop
+    for(int j=0;j<nPoints;j++){
+      x2.Set(10,starty+j*deltay);
+      apHphot[j]->Fill(TotQdc(x2,0));
+      apHmip[j]->Fill(TotQdc(x2,gRandom->Landau(600,150)*1e-9));
+    }
+  }
+  
+  pC->cd(1);  pStackMip->Draw("nostack");
+  pC->cd(2);  pStackPhot->Draw("nostack"); pLeg->Draw();
+}//TestResp()
+//__________________________________________________________________________________________________
+void AliRICHParam::TestTrans()
+{
+//test the set of transformation methods
+  new TCanvas("trasform","Test LRS-MRS transform");
+  TLatex t; t.SetTextSize(0.02);
+  
+  TView *pView=new TView(1);
+  pView->SetRange(-600,-600,-600,600,600,600);
+  DrawAxis();  
+//Draw PC for all chambers by trasfering Pc plane using Pc2Mrs methode  
+  Int_t iNpointsX=50,iNpointsY=50;  
+  for(Int_t iChamberN=1;iChamberN<=7;iChamberN++){//chamber loop
+    TPolyMarker3D *pChamber=new TPolyMarker3D(iNpointsX*iNpointsY);
+    Int_t i=0;
+    for(Double_t x=0;x<PcSizeX();x+=PcSizeX()/iNpointsX)
+      for(Double_t y=0;y<PcSizeY();y+=PcSizeY()/iNpointsY){//step loop
+        TVector3 v3=C(iChamberN)->Pc2Mrs(TVector2(x,y));//from regular grid of local PC points to MRS presentation
+        pChamber->SetPoint(i++,v3.X(),v3.Y(),v3.Z());//Pc plane poing in MRS
+      }//step loop
+    pChamber->SetMarkerSize(1);
+    pChamber->SetMarkerColor(iChamberN);
+    pChamber->Draw();  
+    t.SetNDC();t.SetTextColor(iChamberN); t.DrawText(0.1,iChamberN*0.1,Form("Chamber %i",iChamberN));    
+  }//chamber loop   
+//  gPad->GetView()->RotateView(94,45);
+}//TestTrans()
+//__________________________________________________________________________________________________
+void AliRICHParam::DrawAxis()
+{
+  Double_t X[6]={0,0,0,300,0,0};  Double_t Y[6]={0,0,0,0,300,0};  Double_t Z[6]={0,0,0,0,0,300};  
+  TPolyLine3D *pXaxis=new TPolyLine3D(2,X);pXaxis->SetLineColor(kRed);   pXaxis->Draw();
+  TPolyLine3D *pYaxis=new TPolyLine3D(2,Y);pYaxis->SetLineColor(kGreen); pYaxis->Draw();
+  TPolyLine3D *pZaxis=new TPolyLine3D(2,Z);pZaxis->SetLineColor(kBlue);  pZaxis->Draw();  
+}
index b96b566..2fd04f1 100644 (file)
@@ -31,23 +31,45 @@ class AliRICHChamber;
 class AliRICHParam :public TObject  
 {
 public:
+//ctor&dtor    
                   AliRICHParam():TObject(),fpChambers(0)  {CreateChambers();}
   virtual        ~AliRICHParam()                          {delete fpChambers;}
-         void     CreateChambers();
+//test methodes  
+         void     Print(const Option_t *opt="");                                         //print current parametrization
+         void     Test()                            {TestSeg();TestTrans();TestResp();}  //test all groups of methodes
+         void     TestResp();                                                            //test the response group of methodes
+         void     TestSeg();                                                             //test the segmentation group of methodes
+         void     TestTrans();                                                           //test the transform group of methodes
+         void     DrawAxis();
+//flags staff         
+  static           void     SetAerogel(Bool_t a)                   {fgIsAerogel=a;}
+  static           Bool_t   IsAerogel()                            {return fgIsAerogel;}
+  static           void     SetRadioSrc(Bool_t a)                   {fgIsRadioSrc=a;}
+  static           Bool_t   IsRadioSrc()                            {return fgIsRadioSrc;}
+  static              void     SetTestBeam(Bool_t a)                {fgIsTestBeam=a;}
+  static              Bool_t   IsTestBeam()                         {return fgIsTestBeam;}
+  static                void     SetWireSag(Bool_t a)               {fgIsWireSag=a;}
+  static                Bool_t   IsWireSag()                        {return fgIsWireSag;}
+  static                   void     SetResolveClusters(Bool_t a)    {fgIsResolveClusters=a;}
+  static                   Bool_t   IsResolveClusters()             {return fgIsResolveClusters;}
+//Chambers manipulation  methodes 
+  void            CreateChambers();                                                      //form chamber structure  
   AliRICHChamber* C(Int_t i)                 {return (AliRICHChamber*)fpChambers->UncheckedAt(i-1);}      //returns pointer to chamber i
+  Int_t           Nchambers()                {return fpChambers->GetEntriesFast();}      //returns number of chambers 
+//Geometrical properties  
   static Int_t    NpadsX()                   {return kNpadsX;}                           //pads along X in chamber
   static Int_t    NpadsY()                   {return kNpadsY;}                           //pads along Y in chamber
   static Int_t    NpadsXsec()                {return NpadsX()/2;}                        //pads along X in sector
   static Int_t    NpadsYsec()                {return NpadsY()/3;}                        //pads along Y in sector
   static Double_t DeadZone()                 {return 2.6;}                               //dead zone size in cm  
-  static Double_t PadSizeX()                 {return 0.8;}                               //pad size x in cm 
-  static Double_t PadSizeY()                 {return 0.84;}                              //pad size y in cm   
+  static Double_t PadSizeX()                 {return 0.8;}                               //pad size x, cm 
+  static Double_t PadSizeY()                 {return 0.84;}                              //pad size y, cm   
   
-  static Double_t SectorSizeX()              {return NpadsX()*PadSizeX()/2;}             //sector size x in cm
-  static Double_t SectorSizeY()              {return NpadsY()*PadSizeY()/3;}             //sector size y in cm 
+  static Double_t SectorSizeX()              {return NpadsX()*PadSizeX()/2;}             //sector size x, cm
+  static Double_t SectorSizeY()              {return NpadsY()*PadSizeY()/3;}             //sector size y, cm 
   static Double_t PcSizeX()                  {return NpadsX()*PadSizeX()+DeadZone();}    //PC size x, cm
   static Double_t PcSizeY()                  {return NpadsY()*PadSizeY()+2*DeadZone();}  //PC size y, cm
-   
+//   
   static Double_t Zfreon()                   {return 1.5;}                               //freon thinkness, cm
   static Double_t Zwin()                     {return 0.5;}                               //radiator quartz window, cm   
   static Double_t Pc2Win()                   {return 8.0;}                               //cm between CsI PC and radiator quartz window
@@ -62,60 +84,61 @@ public:
   static Double_t IonisationPotential()      {return 26.0e-9;}                            //for CH4 in GeV taken from ????
   static TVector2 MathiesonDelta()           {return TVector2(5*0.18,5*0.18);}            //area of 5 sigmas of Mathieson distribution (cm)
   static Int_t    MaxQdc()                   {return 4095;}                               //QDC number of channels          
-    
-  static Bool_t   IsResolveClusters()         {return fgIsResolveClusters;}  //go after resolved clusters?
-  static Bool_t   IsWireSag()                 {return fgIsWireSag;}          //take wire sagita in account?
-  static Bool_t   IsRadioSrc()                {return fgIsRadioSrc;}         //add radioactive source inside CH4?
-  static Int_t    HV(Int_t sector)            {
+  
+  static Int_t    HV(Int_t sector)           {
     if (sector>=1 && sector <=6)
       return fgHV[sector-1];
     else {
-      ::Error("HV","Wrong sector %d",sector);
       return kBad;
     } 
   }       //high voltage for this sector
-  static void     SetDeclustering(Bool_t a)   {fgIsResolveClusters=a;}  
-  static void     SetRadioSrc(Bool_t a)       {fgIsRadioSrc=a;}  
-  static void     SetWireSag(Bool_t status)   {fgIsWireSag=status;}  
   static void     SetHV(Int_t sector,Int_t hv){fgHV[sector-1]=hv;}  
-  static void     SetAngleRot(Double_t rot)   {fgAngleRot =rot;}
-  static Double_t IndOfRefC6F14(Double_t eV)  {return eV*0.0172+1.177;}          // eV = photon energy in eV
-  static Double_t IndOfRefSiO2(Double_t eV)   {Double_t e1=10.666,e2=18.125,f1=46.411,f2= 228.71;
-                                     return TMath::Sqrt(1.+f1/(e1*e1-TMath::Power(eV,2))+f2/(e2*e2-TMath::Power(eV,2)));}//TDR p.35
-  static Double_t IndOfRefCH4()               {return 1.000444;}
-
-  inline static TVector  Loc2Area(TVector2 x2);                                                    //return area affected by hit x2
-  inline static TVector  Loc2Pad(TVector2 x2);                                                     //return pad containing given position
+//optical properties methodes  
+  static Float_t  PhotonEnergy(Int_t i)    {return 0.1*i+5.5;}             //photon energy (eV) for i-th point
+  static Float_t  AbsCH4(Float_t ev);                                      //CH4 absorption length (cm) for photon with given energy (eV)
+  static Float_t  AbsGel(Float_t)          {return 500;}                   //Aerogel absorption length (cm) for photon with given energy (eV)
+  static Float_t  RefIdxC6F14(Float_t eV)  {return eV*0.0172+1.177;}       //Freon ref index for photon with given energy (eV)
+  static Float_t  RefIdxCH4(Float_t)       {return 1.000444;}              //Methane ref index for photon with given energy (eV)
+  static Float_t  RefIdxSiO2(Float_t eV)   {Float_t e1=10.666,e2=18.125,f1=46.411,f2= 228.71;
+                                     return TMath::Sqrt(1.+f1/(e1*e1-eV*eV)+f2/(e2*e2-eV*eV));}//Quartz window ref index from TDR p.35
+  static Float_t  RefIdxGel(Float_t)       {return 1.05;}                  //aerogel ref index 
+  static Float_t  DenGel()                 {return (RefIdxGel(0)-1)/0.21;} //aerogel density gr/cm^3 parametrization by E.Nappi
+//trasformation methodes
+  inline static TVector  Loc2Area(const TVector2 &x2);                                             //return area affected by hit x2
+  inline static Int_t    Loc2Sec(const TVector2 &x2);                                              //return sector for given position
+  inline static TVector  Loc2Pad(const TVector2 &x2);                                              //return pad containing given position
   inline static TVector2 Pad2Loc(TVector pad);                                                     //return center of the pad
-         static TVector2 Pad2Loc(Int_t x,Int_t y) {TVector pad(2);pad[0]=x;pad[1]=y;return Pad2Loc(pad);}
+         static TVector2 Pad2Loc(Int_t x,Int_t y) {TVector pad(2);pad[0]=x;pad[1]=y;return Pad2Loc(pad);}//return center of the pad (x,y)
+  inline static Int_t    Pad2Sec(const TVector &pad);                                              //return sector of given pad
   inline static Int_t    PadNeighbours(Int_t iPadX,Int_t iPadY,Int_t aListX[4],Int_t aListY[4]);   //number of neighbours for this pad
-  
-  inline static Double_t Mathieson(Double_t x1,Double_t x2,Double_t y1,Double_t y2);               //Mathienson integral over these limits
+         static Bool_t   IsAccepted(const TVector2 &x2) {return ( x2.X()>=0 && x2.X()<=PcSizeX() && x2.Y()>=0 && x2.Y()<=PcSizeY() ) ? kTRUE:kFALSE;}
+//charge response methodes  
+  inline static Double_t Mathieson(Double_t x1,Double_t x2,Double_t y1,Double_t y2);               //Mathienson integral over given limits
   inline static Double_t GainSag(Double_t x,Int_t sector);                                         //gain variations in %
          static Double_t QdcSlope(Int_t sec){switch(sec){case kBad: return 0;  default:   return 33;}} //weight of electon in QDC channels
-         static Double_t Gain(TVector2 x2){if(IsWireSag()) return QdcSlope(Loc2Sec(x2))*(1+GainSag(x2.X(),Loc2Sec(x2))/100);else return QdcSlope(Loc2Sec(x2));}//gain for point in chamber RS 
-  inline static Double_t FracQdc(TVector2 x2,TVector pad);                                         //charge fraction to pad from hit
+         static Double_t Gain(const TVector2 &x2){//gives chamber gain in terms of QDC channels for given point in local ref system
+                          if(fgIsWireSag) return QdcSlope(Loc2Sec(x2))*(1+GainSag(x2.X(),Loc2Sec(x2))/100);
+                          else            return QdcSlope(Loc2Sec(x2));}
+  inline static Double_t FracQdc(const TVector2 &x2,const TVector &pad);                           //charge fraction to pad from hit
   inline static Int_t    TotQdc(TVector2 x2,Double_t eloss);                                       //total charge for hit eloss=0 for photons
   inline        Bool_t   IsOverTh(Int_t c,TVector pad,Double_t q);                                 //is QDC of the pad registered by FEE  
          static Int_t    NsigmaTh()                    {return fgNsigmaTh;}                        //
          static Float_t  SigmaThMean()                 {return fgSigmaThMean;}                     //QDC electronic noise mean
          static Float_t  SigmaThSpread()               {return fgSigmaThSpread;}                   //QDC electronic noise width
-                void     Print(const Option_t *opt="");                                            //virtual
                 
-  inline static Int_t    Loc2Sec(TVector2 &x2);             //return sector, x2->Sector RS
-  inline static Int_t    Pad2Sec(const TVector &pad);              //return sector
-         static Bool_t   IsAccepted(const TVector2 &x2) {return ( x2.X()>=0 && x2.X()<=PcSizeX() && x2.Y()>=0 && x2.Y()<=PcSizeY() ) 
-? kTRUE:kFALSE;}
-  inline static Double_t CogCorr(Double_t x) {return 3.31267e-2*TMath::Sin(2*TMath::Pi()/PadSizeX()*x)
+  inline static Double_t CogCorr(Double_t x) {return 3.31267e-2*TMath::Sin(2*TMath::Pi()/PadSizeX()*x) //correction of cluster CoG due to sinoidal
                                                     -2.66575e-3*TMath::Sin(4*TMath::Pi()/PadSizeX()*x)
                                                     +2.80553e-3*TMath::Sin(6*TMath::Pi()/PadSizeX()*x);}
+  
+  static Bool_t     fgIsAerogel;                            //aerogel geometry instead of normal RICH flag
 protected:
-         TObjArray *fpChambers;                             //list of chambers    
+  static Bool_t     fgIsRadioSrc;                           //radioactive source instead of radiators flag
+  static Bool_t     fgIsTestBeam;                           //test beam geometry instead of normal RICH flag
   static Bool_t     fgIsWireSag;                            //wire sagitta ON/OFF flag
   static Bool_t     fgIsResolveClusters;                    //declustering ON/OFF flag
-  static Bool_t     fgIsRadioSrc;                           //radioactive source ON/OFF flag
+
+         TObjArray *fpChambers;                             //list of chambers    
   static Int_t      fgHV[6];                                //HV applied to anod wires
-  static Double_t   fgAngleRot;                             //module rotation from up postion (0,0,490)cm
   static Int_t      fgNsigmaTh;                             //n. of sigmas to cut for zero suppression
   static Float_t    fgSigmaThMean;                          //sigma threshold value
   static Float_t    fgSigmaThSpread;                        //spread of sigma
@@ -124,8 +147,8 @@ protected:
 //__________________________________________________________________________________________________
 Int_t AliRICHParam::PadNeighbours(Int_t iPadX,Int_t iPadY,Int_t listX[4],Int_t listY[4])
 {
-// Determines all the neighbouring pads for the given one (iPadX,iPadY). Returns total number of these pads.
-// Dead zones are taken into account.    
+//Determines all the neighbouring pads for the given one (iPadX,iPadY). Returns total number of these pads.
+//Dead zones are taken into account, meaning pads from different sector are not taken. 
 //   1  
 // 2   3
 //   4     
@@ -138,10 +161,10 @@ Int_t AliRICHParam::PadNeighbours(Int_t iPadX,Int_t iPadY,Int_t listX[4],Int_t l
   return nPads;
 }//Pad2ClosePads()
 //__________________________________________________________________________________________________
-Int_t AliRICHParam::Loc2Sec(TVector2 &v2)
+Int_t AliRICHParam::Loc2Sec(const TVector2 &v2)
 {
-// Determines sector containing the given point and trasform this point to the local system of that sector.
-// Returns sector code:                       
+//Determines sector containing the given point.
+//Returns sector code:                       
 //y ^  5 6
 //  |  3 4
 //  |  1 2
@@ -151,44 +174,41 @@ Int_t AliRICHParam::Loc2Sec(TVector2 &v2)
   Double_t y4=PcSizeY()-SectorSizeY();      Double_t y5=PcSizeY();
   
   Int_t sector=kBad;  
-  Double_t x=v2.X(),y=v2.Y();  
-  if     (v2.X() >= x0 && v2.X() <= x1 )  {sector=1;}
-  else if(v2.X() >= x2 && v2.X() <= x3 )  {sector=2; x=v2.X()-x2;}
-  else                                    {return kBad;}
+  if     (v2.X() >= x0 && v2.X() <= x1 )  sector=1;
+  else if(v2.X() >= x2 && v2.X() <= x3 )  sector=2;
+  else                                    return kBad;
   
-  if     (v2.Y() >= y0 && v2.Y() <= y1 )  {}                                  //sectors 1 or 2 
-  else if(v2.Y() >= y2 && v2.Y() <= y3 )  {sector+=2; y=v2.Y()-y2;}           //sectors 3 or 4
-  else if(v2.Y() >= y4 && v2.Y() <= y5 )  {sector+=4; y=v2.Y()-y4;}           //sectors 5 or 6
-  else                                    {return kBad;}
-  v2.Set(x,y);
+  if     (v2.Y() >= y0 && v2.Y() <= y1 )  ;                    //sectors 1 or 2 
+  else if(v2.Y() >= y2 && v2.Y() <= y3 )  sector+=2;           //sectors 3 or 4
+  else if(v2.Y() >= y4 && v2.Y() <= y5 )  sector+=4;           //sectors 5 or 6
+  else                                    return kBad;
   return sector;
 }//Loc2Sec(Double_t x, Double_t y)
 //__________________________________________________________________________________________________
-TVector AliRICHParam::Loc2Pad(TVector2 x2)
+TVector AliRICHParam::Loc2Pad(const TVector2 &loc)
 {
-// Determines pad number TVector(padx,pady) containing the given point x2 defined the chamber RS.
-// Pad count starts in lower left corner from 1,1 to 144,160 in upper right corner of a chamber.
-// Returns sector number of the determined pad.      
+//Determines pad number TVector(padx,pady) containing the given point x2 defined in the chamber RS.
+//Pad count starts in lower left corner from 1,1 to 144,160 in upper right corner of a chamber.
 //y ^  5 6
 //  |  3 4
 //  |  1 2
 //   -------> x  
   TVector pad(2);
-  Int_t sector=Loc2Sec(x2);//trasforms x2 to sector reference system
-  if(sector==kBad) {pad[0]=pad[1]=kBad; return pad;}
-  
-  pad[0]=Int_t(x2.X()/PadSizeX())+1; if(pad[0]>NpadsXsec()) pad[0]= NpadsXsec();       
-  if(sector==2||sector==4||sector==6)   pad[0]+=  NpadsXsec();     
-
-  pad[1]=Int_t(x2.Y()/PadSizeY())+1; if(pad[1]>NpadsYsec()) pad[1]= NpadsYsec();
-  if(sector==3||sector==4)   pad[1]+=NpadsYsec();    
-  if(sector==5||sector==6)   pad[1]+=2*NpadsYsec();     
+  Int_t sec=Loc2Sec(loc);//trasforms x2 to sector reference system
+  if(sec==kBad) {pad[0]=pad[1]=kBad; return pad;}
+//first we deal with x  
+  if(sec==1||sec==3||sec==5)    pad[0]=           Int_t(            loc.X()   / PadSizeX() )+1; //sector 1 or 3 or 5
+  else                          pad[0]=NpadsX() - Int_t( (PcSizeX()-loc.X())  / PadSizeX() )  ; //sector 2 or 4 or 6
+//second deal with y
+       if(sec==1||sec==2)       pad[1]=Int_t(             loc.Y()                / PadSizeY() )+1;           //sector 1 or 2 
+  else if(sec==3||sec==4)       pad[1]=Int_t( (loc.Y()-SectorSizeY()-DeadZone()) / PadSizeY() )+NpadsYsec(); //sector 3 or 4    
+  else                          pad[1]=NpadsY() - Int_t( (PcSizeY()-loc.Y())  / PadSizeY() ); //sector 5 or 6        
   return pad;
 }
 //__________________________________________________________________________________________________
 Int_t AliRICHParam::Pad2Sec(const TVector &pad)
 {
-// Determines sector containing the given pad.
+//Determines sector containing the given pad.
   Int_t sector=kBad;      
   if     (pad[0] >= 1           && pad[0] <=   NpadsXsec() )    {sector=1;}
   else if(pad[0] >  NpadsXsec() && pad[0] <=   NpadsX()    )    {sector=2;} 
@@ -204,7 +224,7 @@ Int_t AliRICHParam::Pad2Sec(const TVector &pad)
 //__________________________________________________________________________________________________
 TVector2 AliRICHParam::Pad2Loc(TVector pad)
 {
-// Returns position of the center of the given pad in local system of the chamber (cm)    
+//Returns position of the center of the given pad in local system of the chamber (cm)    
 // y ^  5 6
 //   |  3 4        sector numbers
 //   |  1 2
@@ -231,9 +251,9 @@ TVector2 AliRICHParam::Pad2Loc(TVector pad)
 //__________________________________________________________________________________________________
 Double_t AliRICHParam::GainSag(Double_t x,Int_t sector)
 {
-// Returns % of gain variation due to wire sagita.
-// All curves are parametrized as per sector basis, so x must be apriory transformed to the Sector RS.    
-// Here x is a distance along wires.  
+//Returns % of gain variation due to wire sagita.
+//All curves are parametrized as per sector basis, so x must be apriory transformed to the Sector RS.    
+//Here x is a distance along wires.  
   x-=SectorSizeX()/2;
   if(x>SectorSizeX()) x-=SectorSizeX(); 
   switch(HV(sector)){
@@ -247,9 +267,9 @@ Double_t AliRICHParam::GainSag(Double_t x,Int_t sector)
 //__________________________________________________________________________________________________
 Int_t AliRICHParam::TotQdc(TVector2 x2,Double_t eloss)
 {
-// Calculates the total charge produced by the eloss in point x2 (Chamber RS).
-// Returns this change parametrised in QDC channels, or 0 if the hit in the dead zone.
-// eloss=0 means photon which produces 1 electron only eloss > 0 for Mip
+//Calculates the total charge produced by the eloss in point x2 (Chamber RS).
+//Returns this change parametrised in QDC channels, or 0 if the hit in the dead zone.
+//eloss=0 means photon which produces 1 electron only eloss > 0 for Mip
   if(Loc2Sec(x2)==kBad) return 0; //hit in the dead zone     
   Int_t iNelectrons=Int_t(eloss/IonisationPotential()); if(iNelectrons==0) iNelectrons=1;
   Double_t qdc=0;
@@ -257,10 +277,10 @@ Int_t AliRICHParam::TotQdc(TVector2 x2,Double_t eloss)
   return Int_t(qdc);
 }
 //__________________________________________________________________________________________________
-Double_t AliRICHParam::FracQdc(TVector2 x2,TVector pad)
+Double_t AliRICHParam::FracQdc(const TVector2 &x2,const TVector &pad)
 {
-// Calculates the charge fraction induced to given pad by the hit from the given point.
-// Integrated Mathieson distribution is used.  
+//Calculates the charge fraction induced to given pad by the hit from the given point.
+//Integrated Mathieson distribution is used.  
   TVector2 center2=Pad2Loc(pad);//gives center of requested pad
   Double_t normXmin=(x2.X()-center2.X()-PadSizeX()/2)  /Pc2Cath();//parametrise for Mathienson
   Double_t normXmax=(x2.X()-center2.X()+PadSizeX()/2)  /Pc2Cath();
@@ -273,8 +293,8 @@ Double_t AliRICHParam::FracQdc(TVector2 x2,TVector pad)
 //__________________________________________________________________________________________________
 Double_t AliRICHParam::Mathieson(Double_t xMin,Double_t yMin,Double_t xMax,Double_t yMax)
 {
-// All arguments are parametrised according to NIM A370(1988)602-603
-// Returns a charge fraction.   
+//All arguments are parametrised according to NIM A370(1988)602-603
+//Returns a charge fraction.   
   const Double_t kSqrtKx3=0.77459667;const Double_t kX2=0.962;const Double_t kX4=0.379;
   const Double_t kSqrtKy3=0.77459667;const Double_t kY2=0.962;const Double_t kY4=0.379;
 
@@ -285,10 +305,10 @@ Double_t AliRICHParam::Mathieson(Double_t xMin,Double_t yMin,Double_t xMax,Doubl
   return 4*kX4*(TMath::ATan(ux2)-TMath::ATan(ux1))*kY4*(TMath::ATan(uy2)-TMath::ATan(uy1));
 }  
 //__________________________________________________________________________________________________
-TVector AliRICHParam::Loc2Area(TVector2 x2)
+TVector AliRICHParam::Loc2Area(const TVector2 &x2)
 {
-// Calculates the area of disintegration for a given point. It's assumed here that this points lays on anode wire.
-// Area is a rectangulare set of pads defined by its left-down and right-up coners.
+//Calculates the area of disintegration for a given point. It's assumed here that this points lays on anode wire.
+//Area is a rectangulare set of pads defined by its left-down and right-up coners.
   TVector area(4);
   TVector pad=Loc2Pad(x2); 
   area[0]=area[2]=pad[0]; area[1]=area[3]=pad[1];//area is just a pad fired  
@@ -301,7 +321,9 @@ TVector AliRICHParam::Loc2Area(TVector2 x2)
 //__________________________________________________________________________________________________
 Bool_t AliRICHParam::IsOverTh(Int_t ,TVector ,Double_t q)
 {
-// Checks if the current q is over threshold and FEE will save this value to data concentrator.
+//Checks if the current q is over threshold and FEE will save this value to data concentrator.
   return (q>NsigmaTh()*(SigmaThMean()+(1.-2*gRandom->Rndm())*SigmaThSpread()));
 }
 #endif //AliRICHParam_h
+
+
index 08cab2e..9982e2e 100644 (file)
@@ -13,18 +13,19 @@ ClassImp(AliRICHTracker)
 //__________________________________________________________________________________________________
 Int_t AliRICHTracker::PropagateBack(AliESD *pESD)
 {
-  //invoked by AliRecontruction for RICH
-  //if ESD doesn't contain tracks, try to reconstruct with particle from STACK 
-  //(this case is just to forsee a standalone RICH simulation
-  TNtupleD *hn=0;
+// Interface callback methode invoked by AliRecontruction during tracking after TOF
+// It steers to different way to provide the final reconstructed information sutable for analisys:
+// 1. AliESD  - reconstructed tracks are used     
+// 2. RICH private ntuple for debug- stack particles used instead of reconstructed tracks     
   AliDebug(1,"Start pattern recognition");
-  if(pESD->GetNumberOfTracks()) RecWithESD(pESD);
+  if(pESD->GetNumberOfTracks())
+    RecWithESD(pESD);
   else
-    RecWithStack(hn);
+    RecWithStack(0);
   AliDebug(1,"Stop pattern recognition");
 
   return 0; // error code: 0=no error;
-} //pure virtual from AliTracker
+}//PropagateBack()
 //__________________________________________________________________________________________________
 void AliRICHTracker::RecWithESD(AliESD *pESD)
 {
@@ -57,15 +58,17 @@ void AliRICHTracker::RecWithESD(AliESD *pESD)
     for(Int_t iClusN=0;iClusN<pRich->Clusters(iChamber)->GetEntries();iClusN++){//clusters loop for intersected chamber
       AliRICHcluster *pClus=(AliRICHcluster*)pRich->Clusters(iChamber)->UncheckedAt(iClusN);//get pointer to current cluster
       Double_t distCurrent=pClus->DistTo(helix.PosPc());//ditance between current cluster and helix intersection with PC
-      if(distCurrent<distMip){distMip=distCurrent;iMipId=iClusN;}//find cluster nearest to the track 
-      
+      if(distCurrent<distMip){
+        distMip=distCurrent;
+        iMipId=iClusN;
+      }//find cluster nearest to the track       
       AliDebug(1,Form("Ploc (%f,%f,%f) dist= %f",helix.Ploc().Mag(),helix.Ploc().Theta()*TMath::RadToDeg(),
-                                                                    helix.Ploc().Phi()*TMath::RadToDeg(),pClus->DistTo(helix.PosPc())));
-    }////clusters loop for intersected chamber
+                                       helix.Ploc().Phi()*TMath::RadToDeg(),pClus->DistTo(helix.PosPc())));
+    }//clusters loop for intersected chamber
     
     AliDebug(1,Form("Min distance cluster: %i dist is %f",iMipId,distMip));
     
-    AliRICHRecon recon(&helix,pRich->Clusters(iChamber),iMipId);
+    AliRICHRecon recon(&helix,pRich->Clusters(iChamber),iMipId); //create reconstruction object for helix, list of clusters
     Double_t thetaCerenkov=recon.ThetaCerenkov(); //search for mean Cerenkov angle for this track
     AliDebug(1,Form("FINAL Theta Cerenkov=%f",thetaCerenkov));
     pTrack->SetRICHsignal(thetaCerenkov);
@@ -80,8 +83,7 @@ void AliRICHTracker::RecWithESD(AliESD *pESD)
 //__________________________________________________________________________________________________
 void AliRICHTracker::RecWithStack(TNtupleD *hn)
 {
-  // reconstruction for particles from STACK
-  //
+// Reconstruction for particles from STACK
   AliRICH *pRich=((AliRICH*)gAlice->GetDetector("RICH"));
   
 //  pRich->GetLoader()->GetRunLoader()->LoadHeader();
@@ -101,7 +103,7 @@ void AliRICHTracker::RecWithStack(TNtupleD *hn)
   if(pRich->GetLoader()->LoadRecPoints()) {AliDebug(1,Form("No clusters found in RICH"));return;}
   pRich->GetLoader()->TreeR()->GetEntry(0);
 
-  for(Int_t iTrackN=0;iTrackN<iNtracks;iTrackN++){//ESD tracks loop
+  for(Int_t iTrackN=0;iTrackN<iNtracks;iTrackN++){//stack particles loop
     TParticle *pParticle = pStack->Particle(iTrackN);
     if(!pParticle) {AliDebug(1,Form("Not a valid TParticle pointer. Track skipped"));continue;}
     AliDebug(1,Form(" PDG code : %i",pParticle->GetPdgCode()));
@@ -150,7 +152,7 @@ void AliRICHTracker::RecWithStack(TNtupleD *hn)
       
       AliDebug(1,Form("Ploc (%f,%f,%f) dist= %f",helix.Ploc().Mag(),helix.Ploc().Theta()*TMath::RadToDeg(),
                                                                     helix.Ploc().Phi()*TMath::RadToDeg(),pClus->DistTo(helix.PosPc())));
-    }////clusters loop for intersected chamber
+    }//clusters loop for intersected chamber
     
     AliDebug(1,Form("Min distance cluster: %i dist is %f",iMipId,distMip));
     hnvec[6]=mipX;hnvec[7]=mipY;
@@ -164,32 +166,27 @@ void AliRICHTracker::RecWithStack(TNtupleD *hn)
     hnvec[13]=(Double_t)pParticle->GetPdgCode();
     if(hn) hn->Fill(hnvec);
     AliDebug(1,Form("FINAL Theta Cerenkov=%f",thetaCerenkov));
-//    pTrack->SetRICHsignal(thetaCerenkov);
-
-//    Double_t richPID[5]={0.2,0.2,0.2,0.2,0.2}; //start with equal probs for (e,mu,pi,k,p)
-//    CalcProb(thetaCerenkov,p0.Mag(),richPID);
-    
-  }//ESD tracks loop
+  }//stack particles loop
   
   pRich->GetLoader()->UnloadRecPoints();
-
   AliDebug(1,"Stop.");  
-} //RecWithStack
-
+}//RecWithStack
+//__________________________________________________________________________________________________
 Int_t AliRICHTracker::LoadClusters(TTree *pTree)
 {
 // Load clusters for RICH
   AliDebug(1,"Start.");  pTree->GetEntry(0);  AliDebug(1,"Stop."); return 0;
 }
-
 //__________________________________________________________________________________________________
 void AliRICHTracker::CalcProb(Double_t thetaCer,Double_t pmod, Double_t *richPID)
 {
-// 
+// Calculates probability to be a electron-muon-pion-kaon-proton 
+// from the given Cerenkov angle and momentum assuming no initial particle composition
+// (i.e. apriory probability to be the particle of the given sort is the same for all sorts)  
   Double_t height[AliPID::kSPECIES];Double_t totalHeight=0;
   for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++){
     Double_t mass = AliPID::ParticleMass(iPart);
-    Double_t refIndex=AliRICHParam::IndOfRefC6F14(6.755);
+    Double_t refIndex=AliRICHParam::RefIdxC6F14(6.755);
     Double_t cosThetaTh = TMath::Sqrt(mass*mass+pmod*pmod)/(refIndex*pmod);
     if(cosThetaTh>=1) {break;}
     Double_t thetaTh = TMath::ACos(cosThetaTh);
@@ -200,3 +197,4 @@ void AliRICHTracker::CalcProb(Double_t thetaCer,Double_t pmod, Double_t *richPID
   }
   for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++) richPID[iPart] = height[iPart]/totalHeight;    
 }//CalcProb
+//__________________________________________________________________________________________________
index bced7de..32d292b 100644 (file)
@@ -102,8 +102,7 @@ void AliRICHv0::StepManager()
   }
   Info("Medium","id=%i (%s)",iTmedId,sTmed);
   
-  TLorentzVector x4;
-  gMC->TrackPosition(x4);
+  TLorentzVector x4; gMC->TrackPosition(x4);
   Float_t glo[3],loc[3];
   glo[0]=x4.X();glo[1]=x4.Y();glo[2]=x4.Z();  
   gMC->Gmtod(glo,loc,1);
@@ -119,4 +118,3 @@ void AliRICHv0::StepManager()
   }
   Info(Form("Step %i",iStepN++),"end of current step\n");
 }//StepManager()
-
index bf5ec4b..dcec638 100644 (file)
@@ -9,8 +9,8 @@
 class AliRICHv0 : public AliRICH 
 {    
 public:
-  AliRICHv0():AliRICH()                                              {;}                      //default ctor
-  AliRICHv0(const char *name, const char *title):AliRICH(name,title) {;}                      //named ctor
+                 AliRICHv0():AliRICH()                                              {;}       //default ctor
+                 AliRICHv0(const char *name, const char *title):AliRICH(name,title) {;}       //named ctor
   virtual       ~AliRICHv0()                                         {;}                      //dtor
   virtual void   Init()                                              {;}                      //interface from AliRICH
   virtual Int_t  IsVersion()                                    const{return 0;}              //interface from AliRICH
index 41a1ff6..1457920 100644 (file)
@@ -9,16 +9,15 @@
 class AliRICHv1 : public AliRICH 
 {
 public:
-                 AliRICHv1():AliRICH()                                               {;}
-                 AliRICHv1(const char *name, const char *title):AliRICH(name,title)  {;}
-  virtual       ~AliRICHv1()                                                         {;}
+                 AliRICHv1():AliRICH()                                               {;}          //default ctor
+                 AliRICHv1(const char *name, const char *title):AliRICH(name,title)  {;}          //named ctor
+  virtual       ~AliRICHv1()                                                         {;}          //dtor
   virtual void   Init()                                                              {;}
-  virtual Int_t  IsVersion()                                                    const{return 1;}
-  
-  virtual void   StepManager();                                         //full slow step manager
-          Bool_t IsLostByFresnel();                                     //checks if the photon lost on Fresnel reflection  
-          void   GenerateFeedbacks(Int_t iChamber,Float_t eloss=0);     //generates feedback photons; eloss=0 for photon
-private:
+  virtual Int_t  IsVersion()                                                    const{return 1;}  
+  virtual void   StepManager();                                                                   //full slow step manager
+          Bool_t IsLostByFresnel();                                                               //checks if the photon lost on Fresnel reflection  
+          void   GenerateFeedbacks(Int_t iChamber,Float_t eloss=0);                               //generates feedback photons; eloss=0 for photon
+protected:
   ClassDef(AliRICHv1,1)                                                 //RICH full version for simulation
 };
                
index 6d3c467..ede7201 100644 (file)
@@ -4,30 +4,30 @@ void Opticals()
   gROOT->Reset();
 #endif   
   int i;
-  const Int_t kNbins=30;
+  const Int_t kNbins=30;//number of photon energy points
   
-  Float_t aPckov[kNbins];  for(i=0;i<kNbins;i++) aPckov[i]=(0.1*i+5.5)*1e-9;//Photons energy bins 5.5 eV - 8.5 eV step 0.1 eV
+  Float_t aPckov[kNbins];  for(i=0;i<kNbins;i++) aPckov[i]=1e-9*AliRICHParam::PhotonEnergy(i);
     
   
-  Float_t aIdxC6F14[kNbins];  
-  Float_t aIdxCH4[kNbins];
-  Float_t aIdxGrid[kNbins];
-  Float_t aIdxSiO2[kNbins];
-  Float_t aIdxOpSiO2[kNbins];
+  Float_t aIdxC6F14[kNbins];  //Freon ref index
+  Float_t aIdxCH4[kNbins];    //Methane ref index
+  Float_t aIdxSiO2[kNbins];   //Quartz ref index
+  Float_t aIdxMetal[kNbins];  //Metal ref index always 0
+  Float_t aIdxGel[kNbins];    //Aerogel ref index
         
   for (i=0;i<kNbins;i++){
-    aIdxC6F14[i] = (Float_t)AliRICHParam::IndOfRefC6F14(aPckov[i]*1e9);
-    aIdxSiO2[i]  = (Float_t)AliRICHParam::IndOfRefSiO2(aPckov[i]*1e9);
-    aIdxCH4[i]   = (Float_t)AliRICHParam::IndOfRefCH4();
-    aIdxGrid[i]    =1;
-    aIdxOpSiO2[i]  =1;
+    aIdxC6F14[i] = AliRICHParam::RefIdxC6F14(aPckov[i]*1e9);
+    aIdxSiO2[i]  = AliRICHParam::RefIdxSiO2(aPckov[i]*1e9);
+    aIdxCH4[i]   = AliRICHParam::RefIdxCH4(5.5);
+    aIdxGel[i]   = AliRICHParam::RefIdxGel(5.5);
+    aIdxMetal[i]  =0;//metal
   } 
     
-  Float_t aAbsC6F14[kNbins]={//New values from A.DiMauro 28.10.03 total 31
+  Float_t aAbsC6F14[kNbins]={//New values from A.DiMauro 28.10.03 total 30
     32701.4219, 17996.1141, 10039.7281, 1799.1230, 1799.1231, 1799.1231, 1241.4091, 179.0987, 179.0986, 179.0987,
       179.0987,   118.9800,    39.5058,   23.7244,   11.1283,    7.1573,    3.6249,   2.1236,   0.7362,   0.5348,
         0.3387,     0.3074,     0.3050,    0.0001,    0.0001,    0.0001,    0.0001,   0.0001,   0.0001,   0.0001};    
-  Float_t aAbsSiO2[kNbins]={//New values from A.DiMauro 28.10.03 total 31
+  Float_t aAbsSiO2[kNbins]={//New values from A.DiMauro 28.10.03 total 30
        34.4338, 30.5424, 30.2584, 31.4928, 31.7868, 17.8397, 9.3410, 6.4492, 6.1128, 5.8128,
         5.5589,  5.2877,  5.0162,  4.7999,  4.5734,  4.2135, 3.7471, 2.6033, 1.5223, 0.9658,
         0.4242,  0.2500,  0.1426,  0.0863,  0.0793,  0.0724, 0.0655, 0.0587, 0.0001, 0.0001};
@@ -36,11 +36,16 @@ void Opticals()
   Float_t aAbsCH4[kNbins];
   Float_t aAbsGrid[kNbins];
   Float_t aAbsCsI[kNbins];
+  Float_t aAbsGel[kNbins];
+  Float_t aAbsRef[kNbins];//abs prob on reflector
+  
     for (i=0;i<kNbins;i++){
-      aAbsCH4[i]         =AbsoCH4(aPckov[i]*1e9); 
+      aAbsCH4[i]         =AliRICHParam::AbsCH4(AliRICHParam::PhotonEnergy(i)); 
       aAbsOpSiO2[i]      =1e-5; 
       aAbsCsI[i]         =1e-4; 
-      aAbsGrid[i]        =1e-4; 
+      aAbsGrid[i]        =1e-4;       
+      aAbsGel[i]         =AliRICHParam::AbsGel(AliRICHParam::PhotonEnergy(i)); 
+      aAbsRef[i]         =0.01;//1% reflector abs prob 
     }
     
   Float_t aQeCsI[kNbins] = {//New values from A.DiMauro 28.10.03 total 31
@@ -52,7 +57,7 @@ void Opticals()
   
   Float_t aQeAll[kNbins];
   for(i=0;i<kNbins;i++){
-    aQeCsI[i]/= (1.0-Fresnel(aPckov[i]*1e9,1.0,0)); //FRESNEL LOSS CORRECTION
+    aQeCsI[i]/= (1.0-Fresnel(AliRICHParam::PhotonEnergy(i),1.0,0)); //FRESNEL LOSS CORRECTION
     aQeAll[i]=1; //QE for all other materials except for PC must be 1.
   }
        
@@ -73,84 +78,98 @@ void Opticals()
     0.1859, 0.1876, 0.1892, 0.1909, 0.2075, 0.2158, 0.0001, 0.0001, 0.0001, 0.0001 };      
 
 //Now plot all the thigs  
-//Freon, Quartz, Opaque ,Methane,CsI,Grid   
-  const Int_t kC6F14M= 24; const Int_t kC6F14C=kRed;  
-  const Int_t kCH4M= 25;   const Int_t kCH4Color=kGreen;  
-  const Int_t kSiO2M= 26;  const Int_t kSiO2C=kBlue;  
-  const Int_t kCsIMarker = 2;   const Int_t kCsIColor =kMagenta;  
+//Markers and colors for different curves:
+    
+  const Double_t kWidth=0.25,kHeight=0.2;  
+  const Int_t kC6F14Marker=24 , kC6F14Color=kRed;  
+  const Int_t kCH4Marker  =25 , kCH4Color  =kGreen;  
+  const Int_t kSiO2M      =26 , kSiO2Color =kBlue;  
+  const Int_t kCsIMarker  = 2 , kCsIColor  =kMagenta;  
+  const Int_t kGelMarker  =27 , kGelColor  =46;  
+  const Int_t kRefMarker  =28 , kRefColor  =47;  
   
   TCanvas *pC=new TCanvas("c1","RICH optics to check",1100,900);
 
   pC->Divide(3,2);           
-           
+//Ref index           
   pC->cd(1);                 
-  TGraph *pIdxC6F14G=new TGraph(kNbins,aPckov,aIdxC6F14);pIdxC6F14G->SetMarkerStyle(kC6F14M); pIdxC6F14G->SetMarkerColor(kC6F14C);
-  TGraph *pIdxSiO2G =new TGraph(kNbins,aPckov,aIdxSiO2); pIdxSiO2G->SetMarkerStyle(kSiO2M);   pIdxSiO2G->SetMarkerColor(kSiO2C);  
-  TGraph *pIdxCH4G  =new TGraph(kNbins,aPckov,aIdxCH4);  pIdxCH4G->SetMarkerStyle(kCH4M);     pIdxCH4G->SetMarkerColor(kCH4Color);
-  TMultiGraph *pIdxMG=new TMultiGraph();  TLegend *pIdxLe=new TLegend(0.6,0.2,0.85,0.4);
-  pIdxMG->Add(pIdxC6F14G);     pIdxLe->AddEntry(pIdxC6F14G,   "C6F14","p");            
-  pIdxMG->Add(pIdxSiO2G);      pIdxLe->AddEntry(pIdxSiO2G,    "SiO2","p");          
-  pIdxMG->Add(pIdxCH4G);       pIdxLe->AddEntry(pIdxCH4G,     "CH4","p");          
-  pIdxMG->Draw("APL");pIdxMG->GetXaxis()->SetTitle("ref index versus energy, GeV");pIdxMG->Draw("APL");
+  TGraph *pIdxC6F14=new TGraph(kNbins,aPckov,aIdxC6F14);pIdxC6F14->SetMarkerStyle(kC6F14Marker);pIdxC6F14->SetMarkerColor(kC6F14Color);
+  TGraph *pIdxSiO2 =new TGraph(kNbins,aPckov,aIdxSiO2); pIdxSiO2 ->SetMarkerStyle(kSiO2M)      ;pIdxSiO2 ->SetMarkerColor(kSiO2Color);  
+  TGraph *pIdxCH4  =new TGraph(kNbins,aPckov,aIdxCH4);  pIdxCH4  ->SetMarkerStyle(kCH4Marker)       ;pIdxCH4  ->SetMarkerColor(kCH4Color);
+  TGraph *pIdxGel  =new TGraph(kNbins,aPckov,aIdxGel);  pIdxGel  ->SetMarkerStyle(kGelMarker)  ;pIdxGel  ->SetMarkerColor(kGelColor);
+  
+  TMultiGraph *pIdxMG=new TMultiGraph();  TLegend *pIdxLe=new TLegend(0.5,0.21,0.5+kWidth,0.21+kHeight);
+  pIdxMG->Add(pIdxC6F14); pIdxLe->AddEntry(pIdxC6F14,"C6F14"  ,"p");            
+  pIdxMG->Add(pIdxSiO2) ; pIdxLe->AddEntry(pIdxSiO2 ,"SiO2"   ,"p");          
+  pIdxMG->Add(pIdxCH4)  ; pIdxLe->AddEntry(pIdxCH4  ,"CH4"    ,"p");          
+  pIdxMG->Add(pIdxGel)  ; pIdxLe->AddEntry(pIdxGel  ,"Aerogel","p");          
+  pIdxMG->Draw("APL");
+  pIdxMG->SetTitle("Refractive index");  pIdxMG->GetXaxis()->SetTitle("energy, GeV");
+  pIdxMG->Draw("APL");
   pIdxLe->Draw();      
-
-
+//Absorbtion
   pC->cd(2);
-  TGraph *pAbsC6F14G=new TGraph(kNbins,aPckov,aAbsC6F14);pAbsC6F14G->SetMarkerStyle(kC6F14M); pAbsC6F14G->SetMarkerColor(kC6F14C);
-  TGraph *pAbsSiO2G =new TGraph(kNbins,aPckov,aAbsSiO2); pAbsSiO2G->SetMarkerStyle(kSiO2M);   pAbsSiO2G->SetMarkerColor(kSiO2C);
-  TMultiGraph *pAbsMG=new TMultiGraph();  TLegend *pAbsLegend=new TLegend(0.6,0.3,0.85,0.5);
-  pAbsMG->Add(pAbsC6F14G);      pAbsLegend->AddEntry(pAbsC6F14G,  "C6F14","p"); 
-  pAbsMG->Add(pAbsSiO2G);       pAbsLegend->AddEntry(pAbsSiO2G,  "SiO2", "p"); 
-  pAbsMG->Draw("APL"); pAbsMG->GetXaxis()->SetTitle("absorbtion length, cm versus energy, GeV");  pAbsMG->Draw("APL");
-  pAbsLegend->Draw();   
-
-  pC->cd(3);      
-  TGraph *pAbsCH4Gr=new TGraph(kNbins,aPckov,aAbsCH4);
-  pAbsCH4Gr->SetMarkerStyle(kCH4M); pAbsCH4Gr->SetMarkerColor(kCH4Color);
-  pAbsCH4Gr->Draw("APL");
-  pAbsCH4Gr->GetXaxis()->SetTitle("energy, GeV");
-  pAbsCH4Gr->SetTitle("CH4 absorption length, cm");
-  pAbsCH4Gr->Draw("APL");
+  gPad->SetLogy();
+  TGraph *pAbsC6F14=new TGraph(kNbins,aPckov,aAbsC6F14);pAbsC6F14->SetMarkerStyle(kC6F14Marker); pAbsC6F14->SetMarkerColor(kC6F14Color);
+  TGraph *pAbsSiO2 =new TGraph(kNbins,aPckov,aAbsSiO2) ;pAbsSiO2 ->SetMarkerStyle(kSiO2M)      ; pAbsSiO2 ->SetMarkerColor(kSiO2Color);
+  TGraph *pAbsCH4  =new TGraph(kNbins,aPckov,aAbsCH4)  ;pAbsCH4  ->SetMarkerStyle(kCH4Marker)  ; pAbsCH4  ->SetMarkerColor(kCH4Color);
+  TGraph *pAbsGel  =new TGraph(kNbins,aPckov,aAbsGel)  ;pAbsGel  ->SetMarkerStyle(kGelMarker)  ; pAbsGel  ->SetMarkerColor(kGelColor);
+  TGraph *pAbsRef  =new TGraph(kNbins,aPckov,aAbsRef)  ;pAbsRef  ->SetMarkerStyle(kRefMarker)  ; pAbsRef  ->SetMarkerColor(kRefColor);
   
+  TMultiGraph *pAbsMG=new TMultiGraph();  TLegend *pAbsLe=new TLegend(0.2,0.15,0.2+kWidth,0.15+kHeight);
+  pAbsMG->Add(pAbsC6F14);      pAbsLe->AddEntry(pAbsC6F14,  "C6F14"    ,"p"); 
+  pAbsMG->Add(pAbsSiO2);       pAbsLe->AddEntry(pAbsSiO2 ,  "SiO2"     ,"p"); 
+  pAbsMG->Add(pAbsCH4);        pAbsLe->AddEntry(pAbsCH4  ,  "CH4"      ,"p"); 
+  pAbsMG->Add(pAbsGel);        pAbsLe->AddEntry(pAbsGel  ,  "Aerogel"  ,"p"); 
+  pAbsMG->Add(pAbsRef);        pAbsLe->AddEntry(pAbsRef  ,  "Reflector","p"); 
+  pAbsMG->Draw("APL"); 
+  pAbsMG->SetTitle("Absorbtion length,cm"); pAbsMG->GetXaxis()->SetTitle("energy, GeV");  
+  pAbsMG->Draw("APL");
+  pAbsLe->Draw();   
+//QE  
   pC->cd(4);
-  TGraph *pQeCsIG=new TGraph(kNbins,aPckov,aQeCsI);
-  pQeCsIG->SetMarkerStyle(kCsIMarker); pQeCsIG->SetMarkerColor(kCsIColor);
-  pQeCsIG->Draw("APL");
-  pQeCsIG->GetXaxis()->SetTitle("energy, GeV");
-  pQeCsIG->SetTitle("QE");
-  pQeCsIG->Draw("APL");
-  
+  TGraph *pQeCsI=new TGraph(kNbins,aPckov,aQeCsI); pQeCsI->SetMarkerStyle(kCsIMarker); pQeCsI->SetMarkerColor(kCsIColor);
+  pQeCsI->Draw("APL");
+  pQeCsI->SetTitle("QE");  pQeCsI->GetXaxis()->SetTitle("energy, GeV");
+  pQeCsI->Draw("APL");  
 //transmission  
   pC->cd(5);
   Float_t mm =0.1;Float_t cm=1.; 
   Float_t aTrC6F14[kNbins],aTrSiO2[kNbins],aTrCH4[kNbins];
   Float_t aTotTr[kNbins];
-  for(Int_t i=0;i<kNbins;i++){
+  for(Int_t i=0;i<kNbins;i++){//calculate probability for photon to survive during transversing a volume of material with absorption length  
     aTrC6F14[i]=TMath::Exp(-15*mm/(aAbsC6F14[i]+0.0001));
     aTrSiO2[i] =TMath::Exp(-5*mm/(aAbsSiO2[i] +0.0001));
     aTrCH4[i]  =TMath::Exp(-8*cm/(aAbsCH4[i]  +0.0001));    
     aTotTr[i]    =aTrC6F14[i]*aTrSiO2[i]*aTrCH4[i]*aQeCsI[i];
   }
-  TGraph *pTrC6F14G=new TGraph(kNbins,aPckov,aTrC6F14);pTrC6F14G->SetMarkerStyle(kC6F14M);pTrC6F14G->SetMarkerColor(kC6F14C);  
-  TGraph *pTrSiO2G=new TGraph(kNbins,aPckov,aTrSiO2);pTrSiO2G->SetMarkerStyle(kSiO2M); pTrSiO2G->SetMarkerColor(kSiO2C);  
-  TGraph *pTrCH4G=new TGraph(kNbins,aPckov,aTrCH4);pTrCH4G->SetMarkerStyle(kCH4M); pTrCH4G->SetMarkerColor(kCH4Color);   
-  TGraph *pTotTrG=new TGraph(kNbins,aPckov,aTotTr);pTotTrG->SetMarkerStyle(30);pTotTrG->SetMarkerColor(kYellow);
-  TMultiGraph *pTrMG=new TMultiGraph();
-  TLegend *pTrLegend=new TLegend(0.2,0.4,0.35,0.6);
-  pTrMG->Add(pQeCsIG);       pTrLegend->AddEntry(pQeCsIG,   "CsI QE", "p");  
-  pTrMG->Add(pTrC6F14G);     pTrLegend->AddEntry(pTrC6F14G, "C6F14", "p");  
-  pTrMG->Add(pTrSiO2G);      pTrLegend->AddEntry(pTrSiO2G,  "SiO2",  "p");          
-  pTrMG->Add(pTrCH4G);       pTrLegend->AddEntry(pTrCH4G,   "CH4",   "p");          
-  pTrMG->Add(pTotTrG);       pTrLegend->AddEntry(pTotTrG,   "total", "p");          
-  pTrMG->Draw("APL");pTrMG->GetXaxis()->SetTitle("transmission versus energy, GeV");pTrMG->Draw("APL");
-  pTrLegend->Draw();      
-
+  TGraph *pTrC6F14=new TGraph(kNbins,aPckov,aTrC6F14);pTrC6F14->SetMarkerStyle(kC6F14Marker);pTrC6F14->SetMarkerColor(kC6F14Color);  
+  TGraph *pTrSiO2 =new TGraph(kNbins,aPckov,aTrSiO2) ;pTrSiO2 ->SetMarkerStyle(kSiO2M)      ;pTrSiO2 ->SetMarkerColor(kSiO2Color);  
+  TGraph *pTrCH4  =new TGraph(kNbins,aPckov,aTrCH4)  ;pTrCH4  ->SetMarkerStyle(kCH4Marker)  ;pTrCH4  ->SetMarkerColor(kCH4Color);   
+  TGraph *pTotTr  =new TGraph(kNbins,aPckov,aTotTr)  ;pTotTr  ->SetMarkerStyle(30)          ;pTotTr  ->SetMarkerColor(kYellow);
+  
+  TMultiGraph *pTrMG=new TMultiGraph();  TLegend *pTrLe=new TLegend(0.2,0.4,0.2+kWidth,0.4+kHeight);
+  pTrMG->Add(pQeCsI);       pTrLe->AddEntry(pQeCsI,   "CsI QE", "p");  
+  pTrMG->Add(pTrC6F14);     pTrLe->AddEntry(pTrC6F14, "C6F14", "p");  
+  pTrMG->Add(pTrSiO2);      pTrLe->AddEntry(pTrSiO2,  "SiO2",  "p");          
+  pTrMG->Add(pTrCH4);       pTrLe->AddEntry(pTrCH4,   "CH4",   "p");          
+  pTrMG->Add(pTotTr);       pTrLe->AddEntry(pTotTr,   "total", "p");          
+  pTrMG->Draw("APL");
+  pTrMG->SetTitle("Transmission");  pTrMG->GetXaxis()->SetTitle("energy, GeV");
+  pTrMG->Draw("APL");
+  pTrLe->Draw();      
+//comparison of new and old staff
   pC->cd(6);
+  TGraph *pQeCsIold=new TGraph(kNbins,aPckov,aQeCsIold); pQeCsIold->SetMarkerStyle(kC6F14Marker);pQeCsIold->SetMarkerColor(kC6F14Color);  
   
-  TMultiGraph *pCompMG=new TMultiGraph;
-  pCompMG->Add(pQeCsIG);
-  pCompMG->Add(new TGraph(kNbins,aPckov,aQeCsIold));
-  pCompMG->Draw("APL");pCompMG->GetXaxis()->SetTitle("comparison of QE versus energy, GeV");pCompMG->Draw("APL");
+  TMultiGraph *pCompMG=new TMultiGraph;  TLegend *pCompLe=new TLegend(0.2,0.6,0.2+kWidth,0.6+kHeight);
+  pCompMG->Add(pQeCsI);       pCompLe->AddEntry(pQeCsI,    "QE new 30.10.03", "p");  
+  pCompMG->Add(pQeCsIold);    pCompLe->AddEntry(pQeCsIold, "QE old 01.01.02", "p");
+  pCompMG->Draw("APL");
+  pCompMG->SetTitle("Comparison of new and old staff");
+  pCompMG->GetXaxis()->SetTitle("energy, GeV");
+  pCompMG->Draw("APL");
+  pCompLe->Draw();      
     
   return;
   TCanvas *pQeC=new TCanvas("pQeC","CsI QE currently all the same",800,900); 
@@ -172,38 +191,6 @@ void Opticals()
     pQeCsIGr->Draw("APL");
   }
 }//main
-
-//__________________________________________________________________________________________________
-Float_t AbsoCH4(Float_t x)
-{//Evaluate the absorbtion lenght of CH4
-  Float_t sch4[9] = {.12,.16,.23,.38,.86,2.8,7.9,28.,80.};              //MB X 10^22
-  Float_t em[9] = {8.1,8.158,8.212,8.267,8.322,8.378,8.435,8.493,8.55};
-  const Float_t kLoschmidt=2.686763e19;                                      // LOSCHMIDT NUMBER IN CM-3
-  const Float_t kPressure=750.,kTemperature=283.;                                      
-  const Float_t pn=kPressure/760.;
-  const Float_t tn=kTemperature/273.16;
-  const Float_t c0=-1.655279e-1;
-  const Float_t c1=6.307392e-2;
-  const Float_t c2=-8.011441e-3;
-  const Float_t c3=3.392126e-4;
-               
-  Float_t crossSection=0;                        
-  if (x<7.75) 
-    crossSection=.06e-22;
-  else if(x>=7.75 && x<=8.1){                 //------ METHANE CROSS SECTION cm-2 ASTROPH. J. 214, L47 (1978)                                               
-       crossSection=(c0+c1*x+c2*x*x+c3*x*x*x)*1.e-18;
-  }else if (x> 8.1){
-    Int_t j=0;
-    while (x<=em[j] || x>=em[j+1]){
-      j++;
-      Float_t a=(sch4[j+1]-sch4[j])/(em[j+1]-em[j]);
-      crossSection=(sch4[j]+a*(x-em[j]))*1e-22;
-    }
-  }//if
-    
-    Float_t density=kLoschmidt*pn/tn; //CH4 molecular density 1/cm-3
-    return 1./(density*crossSection);
-}//AbsoCH4()
 //__________________________________________________________________________________________________
 Float_t Fresnel(Float_t ene,Float_t pdoti, Bool_t pola)
 {//ENE(EV), PDOTI=COS(INC.ANG.), PDOTR=COS(POL.PLANE ROT.ANG.)
diff --git a/RICH/RICHLinkDef.h b/RICH/RICHLinkDef.h
deleted file mode 100644 (file)
index aaeed71..0000000
+++ /dev/null
@@ -1,23 +0,0 @@
-#ifdef __CINT__
-#pragma link off all globals;
-#pragma link off all classes;
-#pragma link off all functions;
-#pragma link C++ class  AliRICH+;
-#pragma link C++ class  AliRICHv0+;
-#pragma link C++ class  AliRICHv1+;
-#pragma link C++ class  AliRICHParam+;
-#pragma link C++ class  AliRICHhit+;
-#pragma link C++ class  AliRICHdigit+;
-#pragma link C++ class  AliRICHcluster+;
-#pragma link C++ class  AliRICHChamber+;
-#pragma link C++ class  AliRICHMap+;
-#pragma link C++ class  AliRICHClusterFinder+;
-#pragma link C++ class  AliRICHRecon+;
-#pragma link C++ class  AliRICHDigitizer+;
-#pragma link C++ class  AliRICHDisplFast+;
-#pragma link C++ class  AliRICHReconstructor+;
-#pragma link C++ class  AliRICHHelix+;
-#pragma link C++ class  AliGenRadioactive+;
-#pragma link C++ class  AliRICHTracker+;
-#pragma link C++ enum   RadioNuclei_t;
-#endif
diff --git a/RICH/RICHbaseLinkDef.h b/RICH/RICHbaseLinkDef.h
new file mode 100644 (file)
index 0000000..0e3a1c8
--- /dev/null
@@ -0,0 +1,12 @@
+#ifdef __CINT__
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+#pragma link C++ class  AliRICHParam+;
+#pragma link C++ class  AliRICHChamber+;
+#pragma link C++ class  AliRICH+;
+#pragma link C++ class  AliRICHhit+;
+#pragma link C++ class  AliRICHdigit+;
+#pragma link C++ class  AliRICHcluster+;
+#pragma link C++ class  AliRICHDisplFast+;
+#endif
diff --git a/RICH/RICHrecLinkDef.h b/RICH/RICHrecLinkDef.h
new file mode 100644 (file)
index 0000000..0d813a3
--- /dev/null
@@ -0,0 +1,11 @@
+#ifdef __CINT__
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+#pragma link C++ class  AliRICHMap+;
+#pragma link C++ class  AliRICHReconstructor+;
+#pragma link C++ class  AliRICHClusterFinder+;
+#pragma link C++ class  AliRICHTracker+;
+#pragma link C++ class  AliRICHRecon+;
+#pragma link C++ class  AliRICHHelix+;
+#endif
diff --git a/RICH/RICHsimLinkDef.h b/RICH/RICHsimLinkDef.h
new file mode 100644 (file)
index 0000000..7eb0194
--- /dev/null
@@ -0,0 +1,10 @@
+#ifdef __CINT__
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+#pragma link C++ class  AliRICHv0+;
+#pragma link C++ class  AliRICHv1+;
+#pragma link C++ class  AliRICHDigitizer+;
+#pragma link C++ class  AliGenRadioactive+;
+#pragma link C++ enum   RadioNuclei_t;
+#endif
index 61a9536..adbc0ff 100644 (file)
@@ -1,21 +1,51 @@
 Module                 :=RICH
 
-include        lib$(Module).pkg
+include        lib$(Module)base.pkg
+SrcBase        :=$(SRCS)
 
-Sources        :=$(SRCS)
+include        lib$(Module)sim.pkg
+SrcSim         :=$(SRCS)
+
+include        lib$(Module)rec.pkg
+SrcRec         :=$(SRCS)
 
 RootTarget      :=$(shell root-config --arch)
-OutputDir      :=$(MY_LIB)/$(Module)
-SoLib          :=$(MY_LIB)/lib$(Module).so
+DirOut         :=$(LIB_MY)/$(Module)
+LibBase                :=$(LIB_MY)/lib$(Module)base.so
+LibSim         :=$(LIB_MY)/lib$(Module)sim.so
+LibRec         :=$(LIB_MY)/lib$(Module)rec.so
+
+
+HdrBase       := $(SrcBase:.cxx=.h)  $(Module)baseLinkDef.h
+HdrSim        := $(SrcSim:.cxx=.h)   $(Module)simLinkDef.h
+HdrRec        := $(SrcRec:.cxx=.h)   $(Module)recLinkDef.h
+
+DictSrcBase   := $(DirOut)/Dict$(Module)base.cxx
+DictObjBase   := $(DictSrcBase:.cxx=.o)
+
+DictSrcSim    := $(DirOut)/Dict$(Module)sim.cxx
+DictObjSim    := $(DictSrcSim:.cxx=.o)
+
+DictSrcRec    := $(DirOut)/Dict$(Module)rec.cxx
+DictObjRec    := $(DictSrcRec:.cxx=.o)
+
+
+
+ObjBase       := $(patsubst %.cxx,$(DirOut)/%.o,$(SrcBase)) $(DictObjBase)
+ObjSim        := $(patsubst %.cxx,$(DirOut)/%.o,$(SrcSim))  $(DictObjSim)
+ObjRec        := $(patsubst %.cxx,$(DirOut)/%.o,$(SrcRec))  $(DictObjRec)
+
+
+DepFile              := $(DirOut)/$(Module).depend
 
 ifeq ($(RootTarget),linuxicc)
- Compiler      :=icc
- CompilerOptions:=-O0 -fpstkchk -I$(shell root-config --incdir) -I$(ALICE_ROOT)/include
- SoOptions      := -Wl,-soname,lib$(Module).so -shared -O -g   
+ Compiler    :=icc
+ CompilerOpt :=-O0 -fpstkchk -I$(shell root-config --incdir) -I$(ALICE_ROOT)/include
+ LibOpt      :=-O -g -shared -Wl
 else
- Compiler       :=g++
- CompilerOptions:=-g -W -Wall -Werror -fPIC -pipe -fmessage-length=0 -Wno-long-long -pedantic-errors -ansi -I$(shell root-config --incdir) -I$(ALICE_ROOT)/include
- SoOptions      :=-O -g -shared -Wl,-soname,lib$(Module).so  
+ Compiler    :=g++
+ CompilerOpt :=-g -W -Wall -Werror -fPIC -pipe -fmessage-length=0 -Wno-long-long -pedantic-errors -ansi -I$(shell root-config --incdir) -I$(ALICE_ROOT)/include
+ LibOpt      :=-O -g -shared -Wl
 endif
 
 ifdef ALIVERBOSE
@@ -24,57 +54,80 @@ else
        Mute            :=@        
 endif        
 
-Headers          := $(Sources:.cxx=.h)  $(Module)LinkDef.h
-
-DictSource    := $(OutputDir)/G__$(Module).cxx
-DictHeader    := $(DictSource:.cxx=.h)
-DictObject    := $(DictSource:.cxx=.o)
-
-
-Objects         := $(patsubst %.cxx,$(OutputDir)/%.o,$(Sources)) $(DictObject)
-
-
-DepFile                := $(OutputDir)/$(Module).d
-
 ##### TARGETS #####
 
-$(SoLib):       $(Objects)
+all:           $(LibBase) $(LibSim) $(LibRec)
+$(LibBase):    $(ObjBase)
        @echo "Creating $@"
-       $(Mute)$(Compiler) $(SoOptions) $^ -o $@        
+       $(Mute)$(Compiler) $(LibOpt) $^ -o $@   
 
-$(DepFile):    $(Headers)
-       @[ -d $(OutputDir) ] || mkdir -p $(OutputDir)
+$(LibSim):     $(ObjSim)
+       @echo "Creating $@"
+       $(Mute)$(Compiler) $(LibOpt) $^ -o $@   
+        
+$(LibRec):     $(ObjRec)
+       @echo "Creating $@"
+       $(Mute)$(Compiler) $(LibOpt) $^ -o $@   
+        
+$(DepFile):    $(HdrBase) $(HdrSim) $(HdrRec)
+       @[ -d $(DirOut) ] || mkdir -p $(DirOut)
        @touch $(DepFile)
        @echo "Generating dependency $@"
-       $(Mute)rmkdepend -f$(DepFile) -p$(OutputDir)/ -- $(CxxFlags) -- $(Sources) 2>/dev/null
+       $(Mute)rmkdepend -f$(DepFile) -p$(DirOut)/ -- $(CompilerOpt) -- $(SrcBase) $(SrcSim) $(SrcRec) 2>/dev/null
 
-$(DictSource):  $(Headers)
+$(DictSrcBase):  $(HdrBase)
+       @echo "Generating $@"
+       $(Mute)rootcint -f $@ -c $(filter -I%,$(CompilerOpt)) $^
+        
+$(DictSrcSim):  $(HdrSim)
+       @echo "Generating $@"
+       $(Mute)rootcint -f $@ -c $(filter -I%,$(CompilerOpt)) $^
+        
+$(DictSrcRec):  $(HdrRec)
        @echo "Generating $@"
-       $(Mute)rootcint -f $@ -c $(filter -I%,$(CompilerOptions)) $^
+       $(Mute)rootcint -f $@ -c $(filter -I%,$(CompilerOpt)) $^
 
-$(DictObject) : $(DictSource)
+$(DictObjBase) : $(DictSrcBase)
        @echo "Compiling $^"        
-       $(Mute)$(Compiler) $(CompilerOptions) -I. -c $^ -o $@
+       $(Mute)$(Compiler) $(CompilerOpt) -I. -c $^ -o $@
+        
+$(DictObjSim) : $(DictSrcSim)
+       @echo "Compiling $^"        
+       $(Mute)$(Compiler) $(CompilerOpt) -I. -c $^ -o $@
+        
+$(DictObjRec) : $(DictSrcRec)
+       @echo "Compiling $^"        
+       $(Mute)$(Compiler) $(CompilerOpt) -I. -c $^ -o $@
 
 show:
-       @echo "Headers: $(Headers)"
-       @echo "Sources: $(Sources)"
-       @echo "Depend : $(DepFile)"
-       @echo "Objects: $(Objects)"
-       @echo "Library: $(SoLib)"
-
-spec:  $(Sources)
+       @echo    "Base Headers: $(HdrBase)"
+       @echo    "Base Sources: $(SrcBase)"
+       @echo    "Base Objects: $(ObjBase)"
+       @echo    "Base    Dict: $(DictSrcBase)"
+       @echo -e "Base Library: $(LibBase)\n"
+        
+       @echo    "Sim Headers: $(HdrSim)"
+       @echo    "Sim Sources: $(SrcSim)"
+       @echo    "Sim Objects: $(ObjSim)"
+       @echo -e "Sim Library: $(LibSim)\n"
+
+       @echo    "Rec Headers: $(HdrRec)"
+       @echo    "Rec Sources: $(SrcRec)"
+       @echo    "Rec Objects: $(ObjRec)"
+       @echo -e "Rec Library: $(LibRec)\n"
+
+spec:  $(SrcBase)
        @echo $^
        @echo $@
                 
 clean:
        @echo "Cleaning..."
-       $(Mute)rm -rf $(Objects) $(DictSource) $(DictHeader) $(SoLib) $(DepFile)                  
+       $(Mute)rm -rf $(DirOut) $(LibBase) $(LibSim) $(LibRec)
 
 ############################ cxx rule #########################################
-$(OutputDir)/%.o : %.cxx
+$(DirOut)/%.o : %.cxx
        @echo $*.cxx
-       $(Mute)$(Compiler) $(CompilerOptions) -c $*.cxx -o $(OutputDir)/$*.o
+       $(Mute)$(Compiler) $(CompilerOpt) -c $*.cxx -o $(DirOut)/$*.o
 ############################ Dependencies #####################################
 
 -include $(DepFile) 
diff --git a/RICH/libRICH.pkg b/RICH/libRICH.pkg
deleted file mode 100644 (file)
index 4b1a095..0000000
+++ /dev/null
@@ -1,9 +0,0 @@
-SRCS   =  AliRICH.cxx AliRICHv0.cxx AliRICHv1.cxx\
-        AliRICHParam.cxx \
-        AliRICHClusterFinder.cxx AliRICHMap.cxx\
-        AliRICHChamber.cxx AliRICHRecon.cxx\
-         AliRICHDisplFast.cxx\
-        AliRICHDigitizer.cxx AliGenRadioactive.cxx AliRICHHelix.cxx AliRICHTracker.cxx AliRICHReconstructor.cxx
-
-HDRS =  $(SRCS:.cxx=.h)
-DHDR= RICHLinkDef.h
diff --git a/RICH/libRICHbase.pkg b/RICH/libRICHbase.pkg
new file mode 100644 (file)
index 0000000..8e358b4
--- /dev/null
@@ -0,0 +1,4 @@
+SRCS:= AliRICHParam.cxx AliRICHChamber.cxx AliRICH.cxx AliRICHDisplFast.cxx
+
+HDRS:= $(SRCS:.cxx=.h)
+DHDR:= RICHbaseLinkDef.h
diff --git a/RICH/libRICHrec.pkg b/RICH/libRICHrec.pkg
new file mode 100644 (file)
index 0000000..0c26b99
--- /dev/null
@@ -0,0 +1,4 @@
+SRCS:=  AliRICHReconstructor.cxx AliRICHClusterFinder.cxx AliRICHTracker.cxx AliRICHRecon.cxx AliRICHHelix.cxx AliRICHMap.cxx 
+
+HDRS:= $(SRCS:.cxx=.h)
+DHDR:= RICHrecLinkDef.h
diff --git a/RICH/libRICHsim.pkg b/RICH/libRICHsim.pkg
new file mode 100644 (file)
index 0000000..01d18de
--- /dev/null
@@ -0,0 +1,5 @@
+SRCS:= AliRICHv0.cxx AliRICHv1.cxx AliRICHDigitizer.cxx \
+       AliGenRadioactive.cxx
+
+HDRS:= $(SRCS:.cxx=.h)
+DHDR:= RICHsimLinkDef.h
index e8a4be4..48ef166 100644 (file)
@@ -10,13 +10,13 @@ void ps(Int_t event=0)  {r->PrintSDigits(event);} //utility print sdigits
 void pd(Int_t event=0)  {r->PrintDigits(event);}  //utility print digits
 void pc(Int_t event=0)  {r->PrintClusters(event);}//utility print clusters
 void pt(Int_t event=0)  {r->PrintTracks(event);}  //utility print tracks
-Int_t ne(Int_t event=0) {r->Nparticles(kElectron,event);}   //utility number of electrons
-Int_t npi0(Int_t event=0) {r->Nparticles(kPi0,event);}   //utility number of electrons
-Int_t npip(Int_t event=0) {r->Nparticles(kPiPlus,event);}   //utility number of electrons
-Int_t npim(Int_t event=0) {r->Nparticles(kPiMinus,event);}   //utility number of electrons
-Int_t nk0(Int_t event=0) {r->Nparticles(kK0,event);}   //utility number of electrons
-Int_t nkp(Int_t event=0) {r->Nparticles(kKPlus,event);}   //utility number of electrons
-Int_t nkm(Int_t event=0) {r->Nparticles(kKMinus,event);}   //utility number of electrons
+Int_t ne(Int_t event=0)   {AliRICH::Nparticles(kElectron,event,al);} //utility number of electrons
+Int_t npi0(Int_t event=0) {AliRICH::Nparticles(kPi0,event,al);}      //utility number of electrons
+Int_t npip(Int_t event=0) {AliRICH::Nparticles(kPiPlus,event,al);}   //utility number of electrons
+Int_t npim(Int_t event=0) {AliRICH::Nparticles(kPiMinus,event,al);}  //utility number of electrons
+Int_t nk0(Int_t event=0)  {AliRICH::Nparticles(kK0,event,al);}       //utility number of electrons
+Int_t nkp(Int_t event=0)  {AliRICH::Nparticles(kKPlus,event,al);}    //utility number of electrons
+Int_t nkm(Int_t event=0)  {AliRICH::Nparticles(kKMinus,event,al);}   //utility number of electrons
 //__________________________________________________________________________________________________
 void pp(int tid)
 {
@@ -61,7 +61,7 @@ Int_t prim(Int_t tid)
 //__________________________________________________________________________________________________
 void Show()
 {  
-//  CreateHists();
+  CreateHists();
   Info("","\n\n\n");
 //load all trees  
   al->LoadHeader(); 
@@ -104,7 +104,6 @@ void Show()
     }//TreeH loop
     Info("Show-HIT","Evt %i->   %i particles %i primaries  %i entries in TreeH %i hits",
                      iEventN,   iNparticles,    iNprims,      iNentries,         iHitsCounter);
-    FillContribs(pPart->GetPdgCode(),pHit->C(),kTRUE);
     
     if(isSdigits){
       rl->TreeS()->GetEntry(0);