V0A air volume at 325cm
[u/mrichter/AliRoot.git] / FIT / FITsim / AliFITv2.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16
17 /////////////////////////////////////////////////////////////////////
18 //                                                                 //
19 // FIT detector full geometry  version 1                             //
20 //
21 //Begin Html       
22 /*
23 <img src="gif/AliFITv2Class.gif">
24 */
25 //End Html
26 //                                                                  //
27 //                                                                  //
28 //////////////////////////////////////////////////////////////////////
29
30 #include <Riostream.h>
31 #include <stdlib.h>
32
33 #include "TGeoCompositeShape.h"
34 #include "TGeoManager.h"
35 #include "TGeoMatrix.h"
36 #include "TGeoVolume.h"
37 #include "TGeoTube.h"
38 #include "TGeoBBox.h"
39 #include "TGeoNode.h"
40
41
42 #include <TGeoGlobalMagField.h>
43 #include <TGraph.h>
44 #include <TLorentzVector.h>
45 #include <TMath.h>
46 #include <TVirtualMC.h>
47 #include <TString.h>
48
49 #include "AliLog.h"
50 #include "AliMagF.h"
51 #include "AliRun.h"
52
53 #include "AliFITHits.h"
54 #include "AliFITv2.h"
55
56 #include "AliMC.h"
57 #include "AliCDBLocal.h"
58 #include "AliCDBStorage.h"
59 #include "AliCDBManager.h"
60 #include "AliCDBEntry.h"
61 #include "AliTrackReference.h"
62
63 ClassImp(AliFITv2)
64
65
66 //--------------------------------------------------------------------
67 AliFITv2::AliFITv2():  AliFIT(),
68                      fIdSens1(0),
69                      fIdSens2(0),
70                      fPMTeff(0x0)
71
72 {
73   //
74   // Standart constructor for T0 Detector version 0
75 }
76 //--------------------------------------------------------------------
77 AliFITv2::AliFITv2(const char *name, const char *title):
78   AliFIT(name,title),
79   fIdSens1(0),
80   fIdSens2(0),
81   fPMTeff(0x0)
82
83 {
84   //
85   // Standart constructor for T0 Detector version 0
86   //
87   fIshunt = 2; 
88   SetPMTeff();
89 }
90 //_____________________________________________________________________________
91
92 AliFITv2::~AliFITv2() 
93 {
94   // desctructor  
95 }
96
97 //-------------------------------------------------------------------------
98 void AliFITv2::CreateGeometry()
99 {
100   //
101   // Create the geometry of FIT Detector version 1 full geometry
102   //
103   // begin Html
104   //
105
106   Int_t *idtmed = fIdtmed->GetArray();
107   Float_t zdetC = 85;
108   Float_t zdetA = 335;
109   
110   Int_t idrotm[999];
111   Double_t x,y,z;
112   Float_t pstartC[3] = {6., 20 ,5};
113   Float_t pstartA[3] = {2.55, 20 ,5};
114   // Float_t pinstart[3] = {3.2,3.2,3.9};
115   Float_t pinstart[3] = {3.2,3.2,4.};
116   Float_t  pmcp[3] = {3.19, 3.19, 2.8}; //MCP
117   Float_t ptop[3] = {1.324, 1.324, 1.};//cherenkov radiator
118   Float_t preg[3] = {1.324, 1.324, 0.05};//photcathode 
119
120   Float_t zV0A = 325.;
121   Float_t pV0Amother[3] = {4.25, 41.25, 0.6};
122   Float_t pV0A[3] = {4.3, 41.2, 0.5};
123
124   AliMatrix(idrotm[901], 90., 0., 90., 90., 180., 0.);
125   
126   //-------------------------------------------------------------------
127   //  T0 volume 
128   //-------------------------------------------------------------------
129   //C side
130
131    Float_t zc[36] = {2,0, 2,0, 2,0, 2,0, 2,2, 2,2, 2,2, 2,2, 0,0,0,0};
132    Float_t za[25] = {2,0, 2,0, 2,2, 0,2, 0,2, 2,0, 2,0, 0,0, 0,0,0,0};
133   
134    
135   TGeoVolumeAssembly * stlinA = new TGeoVolumeAssembly("0STL");  // A side mother 
136   TGeoVolumeAssembly * stlinC = new TGeoVolumeAssembly("0STR");  // C side mother 
137  //T0 interior
138   TVirtualMC::GetMC()->Gsvolu("0INS","BOX",idtmed[kAir],pinstart,3);
139   TGeoVolume *ins = gGeoManager->GetVolume("0INS");
140  // 
141   TGeoTranslation *tr[40];
142   TString nameTr;
143   //C side
144   Float_t xa=-12.8;
145   Float_t ya=-12.8;
146   Int_t itr=0;
147   for (Int_t itrx=0; itrx<5; itrx++) {
148     for (Int_t itry=0; itry<5; itry++) {
149       nameTr = Form("0TR%i",itr+1);
150       z=-pstartA[2]+pinstart[2]/* +za[itr]*/;
151       if(itr !=12){
152         if(TMath::Abs(xa)<10 && TMath::Abs(ya)<10) z= z-2;
153         tr[itr] = new TGeoTranslation(nameTr.Data(),xa,ya, z );
154         tr[itr]->RegisterYourself();
155         stlinA->AddNode(ins,itr,tr[itr]);
156       }
157       printf(" A %i %f %f %f \n",itr, xa, ya, z+zdetA);
158       itr++;
159       ya+=6.4;
160     }
161     ya=-12.8;
162     xa+=6.4;
163   }
164   Float_t xc=-16;
165   Float_t yc=-16;
166   for (Int_t itrx=0; itrx<6; itrx++) {
167     for (Int_t itry=0; itry<6; itry++) {
168       nameTr = Form("0TR%i",itr+1);
169       z=-pstartC[2]+pinstart[2] /*+ zc[itr] */;
170       if (itr!=39 && itr!=40 && itr!=45 &&itr!=46) {
171         if( TMath::Abs(xc)<10 &&  TMath::Abs(yc)<10) z= z+2;
172         tr[itr] = new TGeoTranslation(nameTr.Data(),xc,yc, z );
173         tr[itr]->RegisterYourself();
174         stlinC->AddNode(ins,itr,tr[itr]);
175       }
176       printf(" C %i %f %f %f \n",itr,xc, yc, z+zdetC);
177       itr++;
178       yc+=6.4;
179     }
180     yc=-16;
181     xc+=6.4;
182   }
183     TGeoVolume *alice = gGeoManager->GetVolume("ALIC");
184   alice->AddNode(stlinA,1,new TGeoTranslation(0,0, zdetA ) );
185   //  alice->AddNode(stlinC,1,new TGeoTranslation(0,0, zdetC ) );
186   TGeoRotation * rotC = new TGeoRotation( "rotC",90., 0., 90., 90., 180., 0.);
187   alice->AddNode(stlinC,1, new TGeoCombiTrans(0., 0., -zdetC , rotC) );
188   
189   x=0;
190   y=0;
191    
192   // Entry window (glass)
193   TVirtualMC::GetMC()->Gsvolu("0TOP","BOX",idtmed[kOpGlass],ptop,3); //glass
194   TGeoVolume *top = gGeoManager->GetVolume("0TOP");
195   TVirtualMC::GetMC()->Gsvolu ("0REG", "BOX", idtmed[kOpGlassCathode], preg, 3); 
196   TGeoVolume *cat = gGeoManager->GetVolume("0REG");
197   TVirtualMC::GetMC()->Gsvolu("0MCP","BOX",idtmed[kGlass],pmcp,3); //glass
198   TGeoVolume *mcp = gGeoManager->GetVolume("0MCP");
199  
200   Int_t ntops=0;
201    Float_t xin=0, yin=0;
202    for (Int_t ix=0; ix<2; ix++) {
203      xin = - pinstart[0] + 0.55 + (ix+0.5)*2*ptop[0] ;
204      for (Int_t iy=0; iy<2 ; iy++) {
205        z = - pinstart[2]+ptop[2];
206        yin = - pinstart[1] + 0.55 + (iy+0.5)*2*ptop[1];
207        ntops++;
208        ins->AddNode(top, ntops, new TGeoTranslation(xin,yin,z) );
209        // printf(" 0TOP  full x %f y %f z %f \n", xin, yin, z);
210        z = -pinstart[2] + 2 * ptop[2] + preg[2];
211        ins->AddNode(cat, ntops, new TGeoTranslation(xin,yin,z) );
212
213        // printf(" GEOGEO  %i %i %i %f %f %f %f %f %f \n", ntops, ix, iy,
214        //  xin,yin,x1[ntops],y1[ntops],x1[ntops]+xin,y1[ntops]+yin);
215      }
216    }
217 // MCP
218    z=-pinstart[2] + 2*ptop[2] + 2*preg[2] + pmcp[2];
219   ins->AddNode(mcp, 1 , new TGeoTranslation(0,0,z) );
220
221   //V0A 
222    TVirtualMC::GetMC()->Gsvolu("0V0AM","TUBE",idtmed[kAir],pV0Amother,3);
223    TVirtualMC::GetMC()->Gspos ("0V0AM",1, "ALIC", 0,0,zV0A , 0, "ONLY");
224    TVirtualMC::GetMC()->Gsvolu("0V0A","TUBE",idtmed[kSensAir],pV0A,3);
225    TVirtualMC::GetMC()->Gspos ("0V0A",1, "0V0AM", 0, 0, 0, 0, "ONLY");
226
227
228  
229 }    
230 //------------------------------------------------------------------------
231 void AliFITv2::AddAlignableVolumes() const
232 {
233   //
234   // Create entries for alignable volumes associating the symbolic volume
235   // name with the corresponding volume path. Needs to be syncronized with
236   // eventual changes in the geometry.
237   //
238   TString volPath;
239   TString symName, sn;
240   TString vpAalign = "/ALIC_1/0STL_1";
241   TString vpCalign = "/ALIC_1/0STR_1";
242   for (Int_t imod=0; imod<2; imod++)  {
243     if (imod==0) {volPath  = vpCalign; symName="/ALIC_1/0STL"; }
244     if (imod==1) {volPath  = vpAalign; symName="/ALIC_1/0STR"; }
245     
246     AliDebug(2,"--------------------------------------------");
247     AliDebug(2,Form("volPath=%s\n",volPath.Data()));
248     AliDebug(2,Form("symName=%s\n",symName.Data()));
249     AliDebug(2,"--------------------------------------------");
250     if(!gGeoManager->SetAlignableEntry(symName.Data(),volPath.Data()))
251       AliFatal(Form("Alignable entry %s not created. Volume path %s not valid",
252 symName.Data(),volPath.Data()));
253    }
254 }   
255 //------------------------------------------------------------------------
256 void AliFITv2::CreateMaterials()
257 {
258
259    Int_t isxfld   = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ();
260    Float_t sxmgmx = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max();
261    //   Float_t a,z,d,radl,absl,buf[1];
262    // Int_t nbuf;
263 // AIR
264                                                                                 
265    Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
266    Float_t zAir[4]={6.,7.,8.,18.};
267    Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
268    Float_t dAir = 1.20479E-3;
269    Float_t dAir1 = 1.20479E-11;
270  // Radiator  glass SiO2
271    Float_t aglass[2]={28.0855,15.9994};
272    Float_t zglass[2]={14.,8.};
273    Float_t wglass[2]={1.,2.};
274    Float_t dglass=2.65;
275  // MCP glass SiO2
276    Float_t dglass_mcp=1.3;
277 //*** Definition Of avaible T0 materials ***
278    AliMixture(1, "Vacuum$", aAir, zAir, dAir1,4,wAir);
279    AliMixture(2, "Air$", aAir, zAir, dAir,4,wAir);
280    AliMixture( 4, "MCP glass   $",aglass,zglass,dglass_mcp,-2,wglass);
281    AliMixture( 24, "Radiator Optical glass$",aglass,zglass,dglass,-2,wglass);
282    
283    AliMedium(1, "Air$", 2, 0, isxfld, sxmgmx, 10., .1, 1., .003, .003);
284    AliMedium(3, "Vacuum$", 1, 0, isxfld, sxmgmx, 10., .01, .1, .003, .003);
285    AliMedium(6, "Glass$", 4, 0, isxfld, sxmgmx, 10., .01, .1, .003, .003);
286   
287    AliMedium(16, "OpticalGlass$", 24, 1, isxfld, sxmgmx, 10., .01, .1, .003, .003);
288     AliMedium(19, "OpticalGlassCathode$", 24, 1, isxfld, sxmgmx, 10., .01, .1, .003, .003);
289    AliMedium(22, "SensAir$", 2, 1, isxfld, sxmgmx, 10., .1, 1., .003, .003);
290    
291    AliDebugClass(1,": ++++++++++++++Medium set++++++++++");
292    
293    
294 }
295
296 //-------------------------------------------------------------------
297 void AliFITv2::DefineOpticalProperties()
298 {
299
300
301 // Optical properties definition.
302    Int_t *idtmed = fIdtmed->GetArray();
303 // Definition Cherenkov parameters
304    int i;
305    const Int_t kNbins=31;
306    
307    Float_t rindexSiO2[kNbins],  efficAll[kNbins], rindexAir[kNbins], absorAir[kNbins],rindexCathodeNext[kNbins], absorbCathodeNext[kNbins];
308    Double_t efficMet[kNbins], aReflMet[kNbins];
309            
310    // quartz 20mm
311    Float_t aAbsSiO2[kNbins]={29.0, 28.6, 28.3, 27.7, 27.3, 26.7, 26.4, 
312                              25.9, 25.3, 24.9, 24.5, 23.7, 
313                              23.2, 22.8, 22.4, 21.8, 21.3,
314                              22.8, 22.1, 21.7, 21.2, 20.5, 
315                              19.9, 19.3, 18.7, 18.0, 17.1, 
316                              16.3, 15.3, 14.3, 14.3   };
317           
318    Float_t aPckov[kNbins]  ={3.87, 3.94, 4.02, 4.11, 4.19, 4.29, 4.38,
319                              4.48, 4.58, 4.69, 4.81, 4.93, 
320                              5.05, 5.19, 5.33, 5.48, 5.63,
321                              5.8,  5.97, 6.16, 6.36, 6.57, 
322                              6.8,  7.04, 7.3,  7.58, 7.89, 
323                              8.22, 8.57, 8.97, 9.39 };  
324   Double_t dPckov[kNbins]  ={3.87, 3.94, 4.02, 4.11, 4.19, 4.29, 4.38,
325                              4.48, 4.58, 4.69, 4.81, 4.93,
326                              5.05, 5.19, 5.33, 5.48, 5.63,
327                              5.8,  5.97, 6.16, 6.36, 6.57,
328                              6.8,  7.04, 7.3,  7.58, 7.89,
329                              8.22, 8.57, 8.97, 9.39 };
330
331    /*     
332    Float_t effCathode[kNbins]={0.11, 0.13, 0.15, 0.16, 0.18, 0.19, 0.20,
333                               0.21, 0.22, 0.23, 0.24, 0.26, 
334                               0.27, 0.29, 0.30, 0.29, 0.29, 
335                               0.28, 0.28, 0.27, 0.26, 0.25, 
336                               0.25, 0.23, 0.20, 0.19, 0.17,
337                               0.17, 0.17, 0.2, 0.23};
338    */     
339    //  Float_t aAbsSiO2[kNbins]; //quartz 30mm
340  for(i=0;i<kNbins;i++)
341
342     {
343       aPckov[i]=aPckov[i]*1e-9;//Photons energy bins 4 eV - 8.5 eV step 0.1 eV
344       dPckov[i]=dPckov[i]*1e-9;//Photons energy bins 4 eV - 8.5 eV step 0.1 eV 
345       //      rindexAir[i]=0.0001;
346       rindexAir[i] = 1.;
347       rindexSiO2[i]=1.458; //refractive index for qwarts
348       rindexCathodeNext[i]=0;
349       efficAll[i]=1.;
350       efficMet[i]=0.;
351       aReflMet[i]=1.;
352       //      aAbsSiO2[i]=28.5; //quartz 30mm
353       absorAir[i]=0.3;      
354       absorbCathodeNext[i]=0;
355
356     }
357   
358   TVirtualMC::GetMC()->SetCerenkov (idtmed[kOpGlass], kNbins, aPckov, aAbsSiO2, efficAll, rindexSiO2 );
359    // TVirtualMC::GetMC()->SetCerenkov (idtmed[kOpGlassCathode], kNbins, aPckov, aAbsSiO2, effCathode, rindexSiO2 );
360    TVirtualMC::GetMC()->SetCerenkov (idtmed[kOpGlassCathode], kNbins, aPckov, aAbsSiO2,efficAll , rindexSiO2 );
361    //  TVirtualMC::GetMC()->SetCerenkov (idtmed[kOpAir], kNbins, aPckov,absorAir , efficAll,rindexAir );
362    //   TVirtualMC::GetMC()->SetCerenkov (idtmed[kOpAirNext], kNbins, aPckov,absorbCathodeNext , efficAll, rindexCathodeNext);
363
364    //Define a boarder for radiator optical properties
365    TVirtualMC::GetMC()->DefineOpSurface("surfRd", kUnified /*kGlisur*/,kDielectric_metal,kPolished, 0.);
366    TVirtualMC::GetMC()->SetMaterialProperty("surfRd", "EFFICIENCY", kNbins, dPckov, efficMet);
367    TVirtualMC::GetMC()->SetMaterialProperty("surfRd", "REFLECTIVITY", kNbins, dPckov, aReflMet);
368
369
370 }
371
372 //-------------------------------------------------------------------
373 void AliFITv2::Init()
374 {
375 // Initialises version 0 of the Forward Multiplicity Detector
376 //
377   AliFIT::Init();
378   fIdSens1=TVirtualMC::GetMC()->VolId("0REG");
379   fIdSens2=TVirtualMC::GetMC()->VolId("0V0A");
380
381    AliDebug(1,Form("%s: *** FIT version 1 initialized ***\n",ClassName()));
382 }
383
384 //-------------------------------------------------------------------
385
386 void AliFITv2::StepManager()
387 {
388   //
389   // Called for every step in the T0 Detector
390   //
391   Int_t id,copy,copy1;
392   static Float_t hits[6];
393   static Int_t vol[3];
394   TLorentzVector pos;
395   TLorentzVector mom;
396   
397   //   TClonesArray &lhits = *fHits;
398   
399   if(!TVirtualMC::GetMC()->IsTrackAlive()) return; // particle has disappeared
400   
401   id=TVirtualMC::GetMC()->CurrentVolID(copy);  
402   // Check the sensetive volume
403   if(id==fIdSens1 ) { 
404     if(TVirtualMC::GetMC()->IsTrackEntering()) {
405       TVirtualMC::GetMC()->CurrentVolOffID(1,copy1);
406       vol[1] = copy1;
407       vol[0]=copy;
408       TVirtualMC::GetMC()->TrackPosition(pos);
409       hits[0] = pos[0];
410       hits[1] = pos[1];
411       hits[2] = pos[2];
412       if(pos[2]<0) vol[2] = 0;
413       else vol[2] = 1 ;
414       //      printf(" volumes pmt %i mcp %i side %i x %f y %f z %f\n",  vol[0], vol[1],  vol[2], hits[0], hits[1], hits[2] );
415       
416       Float_t etot=TVirtualMC::GetMC()->Etot();
417       hits[3]=etot;
418       Int_t iPart= TVirtualMC::GetMC()->TrackPid();
419       Int_t partID=TVirtualMC::GetMC()->IdFromPDG(iPart);
420       hits[4]=partID;
421       Float_t ttime=TVirtualMC::GetMC()->TrackTime();
422       hits[5]=ttime*1e12;
423       if (TVirtualMC::GetMC()->TrackPid() == 50000050)   // If particles is photon then ...
424         {
425           //      if(RegisterPhotoE(vol[1]-1,hits[3])) {
426           if(RegisterPhotoE(hits[3])) {
427             AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber(),vol,hits);
428             // Create a track reference at the exit of photocatode
429           }         
430         }
431       
432       //charge particle 
433       if ( TVirtualMC::GetMC()->TrackCharge() )
434         AddTrackReference(gAlice->GetMCApp()->GetCurrentTrackNumber(), AliTrackReference::kFIT);
435       
436      }// trck entering          
437   } //sensitive
438   //V0A
439   if(id==fIdSens2 ) { 
440     if ( TVirtualMC::GetMC()->TrackCharge()  ) {
441       if(TVirtualMC::GetMC()->IsTrackEntering()) {
442         TVirtualMC::GetMC()->TrackPosition(pos);
443         hits[0] = pos[0];
444         hits[1] = pos[1];
445         hits[2] = pos[2];
446         vol[0]=0;
447         vol[1]=0;
448         vol[2]=2;
449         
450         Float_t etot=TVirtualMC::GetMC()->Etot();
451         hits[3]=etot;
452         Int_t iPart= TVirtualMC::GetMC()->TrackPid();
453         Int_t partID=TVirtualMC::GetMC()->IdFromPDG(iPart);
454         hits[4]=partID;
455         Float_t ttime=TVirtualMC::GetMC()->TrackTime();
456         hits[5]=ttime*1e12;
457         AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber(),vol,hits);
458         //      printf(" volumes pmt %i mcp %i vol %i x %f y %f z %f particle %i all \n",  vol[0], vol[1],  vol[2], hits[0], hits[1], hits[2], hits[4]);
459       }
460     }    
461   }
462
463 }
464
465
466 //------------------------------------------------------------------------
467 Bool_t AliFITv2::RegisterPhotoE(Double_t energy)
468 {
469   
470   
471   //  Float_t hc=197.326960*1.e6; //mev*nm
472   Double_t hc=1.973*1.e-6; //gev*nm
473   Float_t lambda=hc/energy;
474   Float_t eff = fPMTeff->Eval(lambda);
475   Double_t  p = gRandom->Rndm();
476   
477   if (p > eff)
478     return kFALSE;
479   
480   return kTRUE;
481 }
482
483 //----------------------------------------------------------------------------
484
485 void AliFITv2::SetPMTeff()
486 {
487   Float_t lambda[50];
488   Float_t eff[50 ] = {0,        0,       0.23619,  0.202909, 0.177913, 
489                     0.175667, 0.17856, 0.190769, 0.206667, 0.230286,
490                     0.252276, 0.256267,0.26,     0.27125,  0.281818,
491                     0.288118, 0.294057,0.296222, 0.301622, 0.290421, 
492                     0.276615, 0.2666,  0.248,    0.23619,  0.227814, 
493                     0.219818, 0.206667,0.194087, 0.184681, 0.167917, 
494                     0.154367, 0.1364,  0.109412, 0.0834615,0.0725283, 
495                     0.0642963,0.05861, 0.0465,   0.0413333,0.032069, 
496                     0.0252203,0.02066, 0.016262, 0.012,    0.00590476,
497                     0.003875, 0.00190, 0,        0,        0          } ;
498   for (Int_t i=0; i<50; i++) lambda[i]=200+10*i; 
499
500   fPMTeff = new TGraph(50,lambda,eff);
501  
502 }
503
504
505
506
507