/[PAMELA software]/PamVMC/s4/src/PamVMCS4Dig.cxx
ViewVC logotype

Annotation of /PamVMC/s4/src/PamVMCS4Dig.cxx

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.5 - (hide annotations) (download)
Fri Jun 12 18:39:27 2009 UTC (15 years, 5 months ago) by pam-rm2
Branch: MAIN
CVS Tags: v1r0, HEAD
Changes since 1.1: +0 -0 lines 
- Introduced user-defined names of output files and random seeds number.
Users can do it use options of PamVMCApplication constructor:
PamVMCApplication(const char* name,  const char *title, const char*
filename="pamtest", Int_t seed=0).
The Random object that I use is TRandom3 object which has astronomical
large period (in case of default initialization 0). All random generators
in the code use this object by calling of gRandom singleton which keeps
it.

- Corrected TOF digitization routine. No problems with TDC hits due to
hadronic interactions anymore.

- Some small changes was done to compile code under Root 5.23. +
geant4_vmc v. 2.6 without any warnings

- Some classes of PamG4RunConfiguartion was changed for geant4_vmc v.
2.6.Some obsolete classes was deleted as soon as developers implemented
regions.

- Navigation was changed from "geomRootToGeant4" to "geomRoot", because on
VMC web page written that as soon as Geant4 has no option ONLY/MANY
translation of overlapped geometry to Geant4 through VGM could be wrong.
I'd like to stay with Root navigation:
http://root.cern.ch/root/vmc/Geant4VMC.html. This should be default
option.

- New Tracker digitization routine written by Sergio was implemented

- PamVMC again became compatible with geant4_vmc v.2.5 and ROOT 5.20.
 The problem was that ROOT developers introduced in TVirtualMC class a new
method SetMagField and new base class:TVirtualMagField from which
user-defined classes shoukd be derived
1 nikolas 1.1 #include "PamVMCS4Dig.h"
2     #include <TRandom.h>
3    
4     ClassImp(PamVMCS4Dig)
5    
6     void PamVMCS4Dig::Digitize(){
7    
8     cout<<"Starting S4 Digitization...."<<endl;
9    
10     //This is simple approximation, reflection + atennuation..
11     //I hope in future we will use for generation and propagation of
12     //optical photons abilitiees of montecarlo...
13    
14    
15    
16    
17     Int_t DEBUG=1;
18     // creato: S. Borisov, INFN Roma2 e MEPHI, Sett 2007
19     Int_t i,j,t,NdF,pmt,S4;
20     Float_t E0,X,Y,Z,x,y,z,V[3],q,w,l;
21    
22     const Float_t Xs[2] = {-24.1, 24.1};
23     const Float_t Ys[2] = {-24.1, 24.1};
24     const Float_t Zs[2] = {-0.5, 0.5};
25     const Float_t Yp[6] = {-20., -17., -1., 2., 17., 20.}; //PMT positions
26    
27     const Float_t E1 = 1e-6; //Energy which we spend to produce 1 photon
28     const Float_t p = 0.1;
29     const Int_t l0 = 500;
30    
31     const Int_t S4p = 32; //pedestal
32    
33     UShort_t S4v = 0; // ADC count
34    
35     Int_t NdFT = 0;
36    
37     Float_t Ert = 0.;
38    
39     TClonesArray * hc = (TClonesArray*)fhitscolmap.GetValue("S4");
40     PamVMCDetectorHit * hit = 0;
41    
42     if(hc){
43     for(i=0;i<hc->GetEntriesFast();i++){
44    
45     hit = (PamVMCDetectorHit*)hc->At(i);
46    
47     Ert += hit->GetEREL();
48     NdF=Int_t(hit->GetEREL()/E1);
49     NdFT=0;
50     X=hit->GetXIN();
51     Y=hit->GetYIN();
52     Z=(Float_t)(random())/(Float_t)(0x7fffffff)-0.5;
53     for(j=0;j<NdF;j++){
54     q=(Float_t)random()/(Float_t)0x7fffffff;
55     w=(Float_t)random()/(Float_t)0x7fffffff;
56     V[0]=p*cos(6.28318*q);
57     V[1]=p*sin(6.28318*q);
58     V[2]=p*(2.*w-1.);
59     pmt=0;
60     x=X;
61     y=Y;
62     z=Z;
63     while(pmt==0 && (x>Xs[0] && x<Xs[1])&&(y>Ys[0] && y<Ys[1])&&(z>Zs[0] && z<Zs[1])){
64     l=0;
65     while(pmt==0 && (x>Xs[0] && x<Xs[1])&&(y>Ys[0] && y<Ys[1])&&(z>Zs[0] && z<Zs[1])){
66     x+=V[0];
67     y+=V[1];
68     z+=V[2];
69     l+=p;
70     }
71     if((x<Xs[0]+p || x>Xs[1]-p)&&(y>Ys[0]+p && y<Ys[1]-p)&&(z>Zs[0]+p && z<Zs[1]-p)){
72     for(t=0;t<3;t++){
73     if(y>=Yp[2*t] && y<Yp[2*t+1]){
74     if(pmt==0)NdFT++;
75     pmt=1;
76     break;
77     }
78     }
79     if(pmt==1)break;
80     V[0]=-V[0];
81     }
82     q=(Float_t)random()/(Float_t)0x7fffffff;
83     w=1-exp(-l/l0);
84     if(q<w)break;
85     q=(Float_t)random()/(Float_t)0x7fffffff;
86     w=0.5;
87     if(q<w)break;
88     if((x>Xs[0]+p && x<Xs[1]-p)&&(y<Ys[0]+p || y>Ys[1]-p)&&(z>Zs[0]+p && z<Zs[1]-p))V[1]=-V[1];
89     if((x>Xs[0]+p && x<Xs[1]-p)&&(y>Ys[0]+p && y<Ys[1]-p)&&(z<Zs[0]+p || z>Zs[1]-p))V[2]=-V[2];
90     x+=V[0];
91     y+=V[1];
92     z+=V[2];
93     l=0;
94     }
95     }
96     }
97     Ert=Ert/0.002;
98     q=(Float_t)(random())/(Float_t)0x7fffffff;
99     w=0.7;
100     E0=4064./7.;
101     if(Ert<1) S4=0;
102     else S4=(Int_t)(4064.*(1.-exp(-(Ert-1.)/E0)));
103     i=S4/4;
104     if(S4%4==0)
105     S4v=S4+S4p;
106     else if(S4%4==1){
107     if(q<w) S4v=S4-1+S4p;
108     else S4v=S4+1+S4p;
109     } else if(S4%4==2) S4v=S4+S4p;
110     else if(S4%4==3){
111     if(q<w) S4v=S4+1+S4p;
112     else S4v=S4-1+S4p;
113     }
114     if (DEBUG)
115     cout<<"Ert_S4 = " << Ert << " --- S4v = " << S4v << endl;
116    
117     fData.clear();
118    
119     fData.push_back(S4v);
120     fData.push_back(0xd800);
121     fData.push_back(0x0300);
122     }
123     }
  ViewVC Help
Powered by ViewVC 1.1.23