/*
$Log$
+ Revision 1.9 2000/05/31 08:19:38 jbarbosa
+ Fixed bug in StepManager
+
Revision 1.8 2000/05/26 17:30:08 jbarbosa
Cerenkov angle now stored within cerenkov data structure.
-////////////////////////////////////////////////////////
-// Manager and hits classes for set:RICH version 0 //
-/////////////////////////////////////////////////////////
+/////////////////////////////////////////////////////////////
+// Manager and hits classes for set: RICH default version //
+/////////////////////////////////////////////////////////////
#include <TTUBE.h>
#include <TNode.h>
#include <TRandom.h>
#include "AliRICHv0.h"
+#include "AliRICHSegResV0.h"
#include "AliRun.h"
#include "AliMC.h"
#include "iostream.h"
AliRICHv0::AliRICHv0(const char *name, const char *title)
: AliRICH(name,title)
{
+ //
+// Version 0
+// Default Segmentation
+ AliRICHSegmentationV0* SegmentationV0 = new AliRICHSegmentationV0;
+//
+// Segmentation parameters
+ SegmentationV0->SetPadSize(0.84,0.80);
+ SegmentationV0->SetDAnod(0.84/2);
+//
+// Geometry parameters
+ AliRICHGeometry* GeometryV0 = new AliRICHGeometryV0;
+ GeometryV0->SetGapThickness(8);
+ GeometryV0->SetProximityGapThickness(.4);
+ GeometryV0->SetQuartzLength(131);
+ GeometryV0->SetQuartzWidth(126.2);
+ GeometryV0->SetQuartzThickness(.5);
+ GeometryV0->SetOuterFreonLength(131);
+ GeometryV0->SetOuterFreonWidth(40.3);
+ GeometryV0->SetInnerFreonLength(131);
+ GeometryV0->SetInnerFreonWidth(40.3);
+ GeometryV0->SetFreonThickness(1);
+//
+// Response parameters
+ AliRICHResponseV0* Rresponse0 = new AliRICHResponseV0;
+ Rresponse0->SetSigmaIntegration(5.);
+ Rresponse0->SetChargeSlope(40.);
+ Rresponse0->SetChargeSpread(0.18, 0.18);
+ Rresponse0->SetMaxAdc(1024);
+ Rresponse0->SetAlphaFeedback(0.05);
+ Rresponse0->SetEIonisation(26.e-9);
+ Rresponse0->SetSqrtKx3(0.77459667);
+ Rresponse0->SetKx2(0.962);
+ Rresponse0->SetKx4(0.379);
+ Rresponse0->SetSqrtKy3(0.77459667);
+ Rresponse0->SetKy2(0.962);
+ Rresponse0->SetKy4(0.379);
+ Rresponse0->SetPitch(0.25);
+//
+//
+// AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH");
+
fCkov_number=0;
fFreon_prod=0;
-
+
fChambers = new TObjArray(7);
for (Int_t i=0; i<7; i++) {
-
- (*fChambers)[i] = new AliRICHChamber();
-
+
+ (*fChambers)[i] = new AliRICHChamber();
+
+ }
+
+ for (Int_t i=0; i<7; i++) {
+ SetGeometryModel(i,GeometryV0);
+ SetSegmentationModel(i, SegmentationV0);
+ SetResponseModel(i, Rresponse0);
+ SetNsec(i,1);
}
}
Int_t nspacers = (Int_t)(TMath::Abs(geometry->GetInnerFreonLength()/14.4));
//printf("\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n Spacers:%d\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n",nspacers);
- printf("Nspacers: %d", nspacers);
+ //printf("Nspacers: %d", nspacers);
//for (i = 1; i <= 9; ++i) {
//zs = (5 - i) * 14.4; //Original settings
gMC->Gspos("META", 1, "SRIC", 0., 1.276, 0., 0, "ONLY");
gMC->Gspos("CSI ", 1, "SRIC", 0., 1.276 + geometry->GetGapThickness()/2 + .25, 0., 0, "ONLY");
- printf("Position of the gap: %f to %f\n", 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 - .2, 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 + .2);
+ //printf("Position of the gap: %f to %f\n", 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 - .2, 1.276 + geometry->GetGapThickness()/2 - geometry->GetProximityGapThickness()/2 + .2);
// Place RICH inside ALICE apparatus
void AliRICHv0::Init()
{
- printf("\n\n\n Start Init for version 0 - CPC chamber type \n\n\n");
-
+
+ printf("*********************************** RICH_INIT ***********************************\n");
+ printf("* *\n");
+ printf("* AliRICHv0 Default version started *\n");
+ printf("* *\n");
+
+
+ AliRICHSegmentation* segmentation;
+ AliRICHGeometry* geometry;
+ AliRICHResponse* response;
+
+
//
// Initialize Tracking Chambers
//
Float_t pos7[3]={-161.399,443.3999,-165.3};
Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],new TRotMatrix("rot999","rot999",90,20,108.1999,110,18.2,110));
- printf("\n\n\n Finished Init for version 0 - CPC chamber type\n\n\n");
+ segmentation=Chamber(0).GetSegmentationModel(0);
+ geometry=Chamber(0).GetGeometryModel();
+ response=Chamber(0).GetResponseModel();
+
+
+ printf("* Pads : %3dx%3d *\n",segmentation->Npx(),segmentation->Npy());
+ printf("* Pad size : %5.2f x%5.2f mm2 *\n",segmentation->Dpx(),segmentation->Dpy());
+ printf("* Gap Thickness : %5.1f mm *\n",geometry->GetGapThickness());
+ printf("* Radiator Width : %5.1f mm *\n",geometry->GetQuartzWidth());
+ printf("* Radiator Length : %5.1f mm *\n",geometry->GetQuartzLength());
+ printf("* Freon Thickness : %5.1f mm *\n",geometry->GetFreonThickness());
+ printf("* Charge Slope : %5.1f ADC *\n",response->ChargeSlope());
+ printf("* Feedback Prob. : %5.2f %% *\n",response->AlphaFeedback());
+ printf("* *\n");
+ printf("* Success! *\n");
+ printf("* *\n");
+ printf("*********************************************************************************\n");
+
}
//___________________________________________
void AliRICHv0::StepManager()
{
- Int_t copy, id;
- static Int_t idvol;
- static Int_t vol[2];
- Int_t ipart;
- static Float_t hits[18];
- static Float_t Ckov_data[19];
- TLorentzVector Position;
- TLorentzVector Momentum;
- Float_t pos[3];
- Float_t mom[4];
- Float_t Localpos[3];
- Float_t Localmom[4];
- Float_t Localtheta,Localphi;
- Float_t theta,phi;
- Float_t destep, step;
- Float_t ranf[2];
- Int_t NPads;
- Float_t coscerenkov;
- static Float_t eloss, xhit, yhit, tlength;
- const Float_t big=1.e10;
-
- TClonesArray &lhits = *fHits;
- TGeant3 *geant3 = (TGeant3*) gMC;
- TParticle *current = (TParticle*)(*gAlice->Particles())[gAlice->CurrentTrack()];
-
- //if (current->Energy()>1)
- //{
-
- // Only gas gap inside chamber
- // Tag chambers and record hits when track enters
-
- idvol=-1;
- id=gMC->CurrentVolID(copy);
- Float_t cherenkov_loss=0;
- //gAlice->KeepTrack(gAlice->CurrentTrack());
-
- gMC->TrackPosition(Position);
- pos[0]=Position(0);
- pos[1]=Position(1);
- pos[2]=Position(2);
- Ckov_data[1] = pos[0]; // X-position for hit
- Ckov_data[2] = pos[1]; // Y-position for hit
- Ckov_data[3] = pos[2]; // Z-position for hit
- //Ckov_data[11] = gAlice->CurrentTrack();
-
- AliRICH *RICH = (AliRICH *) gAlice->GetDetector("RICH");
-
- /********************Store production parameters for Cerenkov photons************************/
-//is it a Cerenkov photon?
- if (gMC->TrackPid() == 50000050) {
-
- //if (gMC->VolId("GAP ")==gMC->CurrentVolID(copy))
- //{
- Float_t Ckov_energy = current->Energy();
- //energy interval for tracking
- if (Ckov_energy > 5.6e-09 && Ckov_energy < 7.8e-09 )
- //if (Ckov_energy > 0)
- {
- if (gMC->IsTrackEntering()){ //is track entering?
- if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
- { //is it in freo?
- if (geant3->Gctrak()->nstep<1){ //is it the first step?
- Int_t mother = current->GetFirstMother();
-
- //printf("Second Mother:%d\n",current->GetSecondMother());
-
- Ckov_data[10] = mother;
- Ckov_data[11] = gAlice->CurrentTrack();
- Ckov_data[12] = 1; //Media where photon was produced 1->Freon, 2->Quarz
- fCkov_number++;
- fFreon_prod=1;
- //printf("Index: %d\n",fCkov_number);
- } //first step question
- } //freo question
-
- if (geant3->Gctrak()->nstep<1){ //is it first step?
- if (gMC->VolId("QUAR")==gMC->CurrentVolID(copy)) //is it in quarz?
- {
- Ckov_data[12] = 2;
- } //quarz question
- } //first step question
-
- //printf("Before %d\n",fFreon_prod);
- } //track entering question
-
- if (Ckov_data[12] == 1) //was it produced in Freon?
- //if (fFreon_prod == 1)
- {
- if (gMC->IsTrackEntering()){ //is track entering?
- //printf("Got in");
- if (gMC->VolId("META")==gMC->CurrentVolID(copy)) //is it in gap?
- {
- //printf("Got in\n");
- gMC->TrackMomentum(Momentum);
- mom[0]=Momentum(0);
- mom[1]=Momentum(1);
- mom[2]=Momentum(2);
- mom[3]=Momentum(3);
- // Z-position for hit
-
-
- /**************** Photons lost in second grid have to be calculated by hand************/
-
- Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
- Float_t t = (1. - .025 / cophi) * (1. - .05 / cophi);
- gMC->Rndm(ranf, 1);
- //printf("grid calculation:%f\n",t);
- if (ranf[0] > t) {
- geant3->StopTrack();
- Ckov_data[13] = 5;
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- //printf("Lost one in grid\n");
- }
- /**********************************************************************************/
- } //gap
-
- if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy)) //is it in csi?
- {
- gMC->TrackMomentum(Momentum);
- mom[0]=Momentum(0);
- mom[1]=Momentum(1);
- mom[2]=Momentum(2);
- mom[3]=Momentum(3);
-
- /********* Photons lost by Fresnel reflection have to be calculated by hand********/
- /***********************Cerenkov phtons (always polarised)*************************/
-
- Float_t cophi = TMath::Cos(TMath::ATan2(mom[0], mom[1]));
- Float_t t = Fresnel(Ckov_energy*1e9,cophi,1);
- gMC->Rndm(ranf, 1);
- if (ranf[0] < t) {
- geant3->StopTrack();
- Ckov_data[13] = 6;
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- //printf("Lost by Fresnel\n");
- }
- /**********************************************************************************/
- }
- } //track entering?
-
-
- /********************Evaluation of losses************************/
- /******************still in the old fashion**********************/
-
- Int_t i1 = geant3->Gctrak()->nmec; //number of physics mechanisms acting on the particle
- for (Int_t i = 0; i < i1; ++i) {
- // Reflection loss
- if (geant3->Gctrak()->lmec[i] == 106) { //was it reflected
- Ckov_data[13]=10;
- if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
- Ckov_data[13]=1;
- if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
- Ckov_data[13]=2;
- geant3->StopTrack();
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- } //reflection question
-
-
- // Absorption loss
- else if (geant3->Gctrak()->lmec[i] == 101) { //was it absorbed?
- Ckov_data[13]=20;
- if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
- Ckov_data[13]=11;
- if (gMC->CurrentVolID(copy) == gMC->VolId("QUAR"))
- Ckov_data[13]=12;
- if (gMC->CurrentVolID(copy) == gMC->VolId("META"))
- Ckov_data[13]=13;
- if (gMC->CurrentVolID(copy) == gMC->VolId("GAP "))
- Ckov_data[13]=13;
-
- if (gMC->CurrentVolID(copy) == gMC->VolId("SRIC"))
- Ckov_data[13]=15;
-
- // CsI inefficiency
- if (gMC->CurrentVolID(copy) == gMC->VolId("CSI ")) {
- Ckov_data[13]=16;
- }
- geant3->StopTrack();
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- //printf("Added cerenkov %d\n",fCkov_number);
- } //absorption question
-
-
- // Photon goes out of tracking scope
- else if (geant3->Gctrak()->lmec[i] == 30) { //is it below energy treshold?
- Ckov_data[13]=21;
- geant3->StopTrack();
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- } // energy treshold question
- } //number of mechanisms cycle
- /**********************End of evaluation************************/
- } //freon production question
- } //energy interval question
- //}//inside the proximity gap question
- } //cerenkov photon question
-
- /**************************************End of Production Parameters Storing*********************/
-
-
- /*******************************Treat photons that hit the CsI (Ckovs and Feedbacks)************/
-
- if (gMC->TrackPid() == 50000050 || gMC->TrackPid() == 50000051) {
- //printf("Cerenkov\n");
- if (gMC->VolId("CSI ")==gMC->CurrentVolID(copy))
- {
-
- if (gMC->Edep() > 0.){
- gMC->TrackPosition(Position);
- gMC->TrackMomentum(Momentum);
- pos[0]=Position(0);
- pos[1]=Position(1);
- pos[2]=Position(2);
- mom[0]=Momentum(0);
- mom[1]=Momentum(1);
- mom[2]=Momentum(2);
- mom[3]=Momentum(3);
- Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
- Double_t rt = TMath::Sqrt(tc);
- theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
- phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
- gMC->Gmtod(pos,Localpos,1);
- gMC->Gmtod(mom,Localmom,2);
-
- gMC->CurrentVolOffID(2,copy);
- vol[0]=copy;
- idvol=vol[0]-1;
-
- //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
- //->Sector(Localpos[0], Localpos[2]);
- //printf("Sector:%d\n",sector);
-
- /*if (gMC->TrackPid() == 50000051){
- fFeedbacks++;
- printf("Feedbacks:%d\n",fFeedbacks);
- }*/
-
- ((AliRICHChamber*) (*fChambers)[idvol])
- ->SigGenInit(Localpos[0], Localpos[2], Localpos[1]);
- if(idvol<7) {
- Ckov_data[0] = gMC->TrackPid(); // particle type
- Ckov_data[1] = pos[0]; // X-position for hit
- Ckov_data[2] = pos[1]; // Y-position for hit
- Ckov_data[3] = pos[2]; // Z-position for hit
- Ckov_data[4] = theta; // theta angle of incidence
- Ckov_data[5] = phi; // phi angle of incidence
- Ckov_data[8] = (Float_t) fNPadHits; // first padhit
- Ckov_data[9] = -1; // last pad hit
- Ckov_data[13] = 4; // photon was detected
- Ckov_data[14] = mom[0];
- Ckov_data[15] = mom[1];
- Ckov_data[16] = mom[2];
-
- destep = gMC->Edep();
- gMC->SetMaxStep(big);
- cherenkov_loss += destep;
- Ckov_data[7]=cherenkov_loss;
-
- NPads = MakePadHits(Localpos[0],Localpos[2],cherenkov_loss,idvol,cerenkov);
- if (fNPadHits > (Int_t)Ckov_data[8]) {
- Ckov_data[8]= Ckov_data[8]+1;
- Ckov_data[9]= (Float_t) fNPadHits;
- }
-
- Ckov_data[17] = NPads;
- //printf("Npads:%d",NPads);
-
- //TClonesArray *Hits = RICH->Hits();
- AliRICHHit *mipHit = (AliRICHHit*) (fHits->UncheckedAt(0));
- if (mipHit)
- {
- mom[0] = current->Px();
- mom[1] = current->Py();
- mom[2] = current->Pz();
- Float_t Mip_px = mipHit->fMomX;
- Float_t Mip_py = mipHit->fMomY;
- Float_t Mip_pz = mipHit->fMomZ;
-
- Float_t r = mom[0]*mom[0] + mom[1]*mom[1] + mom[2]*mom[2];
- Float_t rt = TMath::Sqrt(r);
- Float_t Mip_r = Mip_px*Mip_px + Mip_py*Mip_py + Mip_pz*Mip_pz;
- Float_t Mip_rt = TMath::Sqrt(Mip_r);
- if ((rt*Mip_rt) > 0)
- {
- coscerenkov = (mom[0]*Mip_px + mom[1]*Mip_py + mom[2]*Mip_pz)/(rt*Mip_rt);
- }
- else
- {
- coscerenkov = 0;
- }
- Float_t cherenkov = TMath::ACos(coscerenkov);
- Ckov_data[18]=cherenkov;
- }
- //if (sector != -1)
- //{
- AddHit(gAlice->CurrentTrack(),vol,Ckov_data);
- AddCerenkov(gAlice->CurrentTrack(),vol,Ckov_data);
- //}
- }
- }
- }
- }
-
- /***********************************************End of photon hits*********************************************/
-
-
- /**********************************************Charged particles treatment*************************************/
-
- else if (gMC->TrackCharge())
- //else if (1 == 1)
- {
-//If MIP
- /*if (gMC->IsTrackEntering())
- {
- hits[13]=20;//is track entering?
- }*/
- if (gMC->VolId("FRE1")==gMC->CurrentVolID(copy) || gMC->VolId("FRE2")==gMC->CurrentVolID(copy))
- {
- fFreon_prod=1;
- }
-
- if (gMC->VolId("GAP ")== gMC->CurrentVolID(copy)) {
-// Get current particle id (ipart), track position (pos) and momentum (mom)
-
- gMC->CurrentVolOffID(3,copy);
- vol[0]=copy;
- idvol=vol[0]-1;
-
- //Int_t sector=((AliRICHChamber*) (*fChambers)[idvol])
- //->Sector(Localpos[0], Localpos[2]);
- //printf("Sector:%d\n",sector);
-
- gMC->TrackPosition(Position);
- gMC->TrackMomentum(Momentum);
- pos[0]=Position(0);
- pos[1]=Position(1);
- pos[2]=Position(2);
- mom[0]=Momentum(0);
- mom[1]=Momentum(1);
- mom[2]=Momentum(2);
- mom[3]=Momentum(3);
- gMC->Gmtod(pos,Localpos,1);
- gMC->Gmtod(mom,Localmom,2);
-
- ipart = gMC->TrackPid();
- //
- // momentum loss and steplength in last step
- destep = gMC->Edep();
- step = gMC->TrackStep();
-
- //
- // record hits when track enters ...
- if( gMC->IsTrackEntering()) {
-// gMC->SetMaxStep(fMaxStepGas);
- Double_t tc = mom[0]*mom[0]+mom[1]*mom[1];
- Double_t rt = TMath::Sqrt(tc);
- theta = Float_t(TMath::ATan2(rt,Double_t(mom[2])))*kRaddeg;
- phi = Float_t(TMath::ATan2(Double_t(mom[1]),Double_t(mom[0])))*kRaddeg;
-
-
- Double_t Localtc = Localmom[0]*Localmom[0]+Localmom[2]*Localmom[2];
- Double_t Localrt = TMath::Sqrt(Localtc);
- Localtheta = Float_t(TMath::ATan2(Localrt,Double_t(Localmom[1])))*kRaddeg;
- Localphi = Float_t(TMath::ATan2(Double_t(Localmom[2]),Double_t(Localmom[0])))*kRaddeg;
-
- hits[0] = Float_t(ipart); // particle type
- hits[1] = Localpos[0]; // X-position for hit
- hits[2] = Localpos[1]; // Y-position for hit
- hits[3] = Localpos[2]; // Z-position for hit
- hits[4] = Localtheta; // theta angle of incidence
- hits[5] = Localphi; // phi angle of incidence
- hits[8] = (Float_t) fNPadHits; // first padhit
- hits[9] = -1; // last pad hit
- hits[13] = fFreon_prod; // did id hit the freon?
- hits[14] = mom[0];
- hits[15] = mom[1];
- hits[16] = mom[2];
-
- tlength = 0;
- eloss = 0;
- fFreon_prod = 0;
-
- Chamber(idvol).LocaltoGlobal(Localpos,hits+1);
-
-
- //To make chamber coordinates x-y had to pass LocalPos[0], LocalPos[2]
- xhit = Localpos[0];
- yhit = Localpos[2];
- // Only if not trigger chamber
- if(idvol<7) {
- //
- // Initialize hit position (cursor) in the segmentation model
- ((AliRICHChamber*) (*fChambers)[idvol])
- ->SigGenInit(Localpos[0], Localpos[2], Localpos[1]);
- }
- }
-
- //
- // Calculate the charge induced on a pad (disintegration) in case
- //
- // Mip left chamber ...
- if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
- gMC->SetMaxStep(big);
- eloss += destep;
- tlength += step;
-
-
- // Only if not trigger chamber
- if(idvol<7) {
- if (eloss > 0)
- {
- if(gMC->TrackPid() == kNeutron)
- printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
- NPads = MakePadHits(xhit,yhit,eloss,idvol,mip);
- hits[17] = NPads;
- //printf("Npads:%d",NPads);
- }
- }
-
- hits[6]=tlength;
- hits[7]=eloss;
- if (fNPadHits > (Int_t)hits[8]) {
- hits[8]= hits[8]+1;
- hits[9]= (Float_t) fNPadHits;
- }
-
- //if(sector !=-1)
- new(lhits[fNhits++]) AliRICHHit(fIshunt,gAlice->CurrentTrack(),vol,hits);
- eloss = 0;
- //
- // Check additional signal generation conditions
- // defined by the segmentation
- // model (boundary crossing conditions)
- } else if
- (((AliRICHChamber*) (*fChambers)[idvol])
- ->SigGenCond(Localpos[0], Localpos[2], Localpos[1]))
- {
- ((AliRICHChamber*) (*fChambers)[idvol])
- ->SigGenInit(Localpos[0], Localpos[2], Localpos[1]);
- if (eloss > 0)
- {
- if(gMC->TrackPid() == kNeutron)
- printf("\n\n\n\n\n Neutron Making Pad Hit!!! \n\n\n\n");
- NPads = MakePadHits(xhit,yhit,eloss,idvol,mip);
- hits[17] = NPads;
- //printf("Npads:%d",NPads);
- }
- xhit = Localpos[0];
- yhit = Localpos[2];
- eloss = destep;
- tlength += step ;
- //
- // nothing special happened, add up energy loss
- } else {
- eloss += destep;
- tlength += step ;
- }
- }
- }
- /*************************************************End of MIP treatment**************************************/
- //}
+ //Dummy step manager
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff