/[PAMELA software]/PamVMC/nd/src/PamVMCNDDig.cxx
ViewVC logotype

Annotation of /PamVMC/nd/src/PamVMCNDDig.cxx

Parent Directory Parent Directory | Revision Log Revision Log

Revision 1.5 - (hide annotations) (download)
Fri Jun 12 18:39:25 2009 UTC (15 years, 6 months ago) by pam-rm2
Branch: MAIN
CVS Tags: v1r0, HEAD
Changes since 1.1: +1 -1 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 "PamVMCNDDig.h"
2     #include <TRandom.h>
3     #include <TMath.h>
4    
5    
6     #define DEBUG 0
7    
8     using TMath::Log;
9     using TMath::Exp;
10     using TMath::Abs;
11    
12     ClassImp(PamVMCNDDig)
13    
14     void PamVMCNDDig::Digitize(){
15    
16     cout<<"Starting ND Digitization...."<<endl;
17    
18     Int_t NdN = 0;
19    
20     const Double_t Showercut = 10.e-6; // 10 us
21     const Double_t Window = 200.e-6; //200 us
22     const Double_t Tau = 17.e-6; // 17 us. Time when amplitude of signal reduce
23     // in e times
24     const Double_t SigmaTau = 4.e-6; //Sigma of fall time (obtained from
25     //Lebedev's oscillograms)
26     const Double_t Tresh = 130; //According Lebedev th = 130 mV
27     const Double_t keVtomV = Tresh/200.; //We suppose that 130 mV = 200 keV... everything linear
28    
29    
30     Double_t deadtime[36];
31     Double_t Time,U, Uback;
32     Int_t pos;
33    
34     for (Int_t i = 0; i<36; i++) deadtime[i] = Showercut;
35    
36     TClonesArray * hc = (TClonesArray*)fhitscolmap.GetValue("NDTI");
37    
38    
39     PamVMCDetectorHit * hit = 0;
40     PamVMCDetectorHit * hit1 = 0;
41     PamVMCDetectorHit *hit0 =0;
42    
43    
44     if(hc){
45    
46     Int_t HitsNum = hc->GetEntriesFast();
47    
48     //sorting by time and put results into vector fpHits
49     if (HitsNum > 1){
50     for(Int_t i=0; i<hc->GetEntriesFast(); i++){
51     Time = 0.;
52     for(Int_t j=0; j<hc->GetEntriesFast(); j++){
53     hit = (PamVMCDetectorHit*)hc->At(j);
54     if((hit->GetTOF()>Time) && (!IsInsideVector(hit))){
55     hit1 = hit;
56     Time = hit->GetTOF();
57     }
58     }
59    
60     fpHits.push_back(hit1);
61     }
62     }
63    
64     //now we separate hits by tubes
65    
66     for(Int_t i = fpHits.size()-1; i>=0; i--){
67     hit = (PamVMCDetectorHit*)fpHits.at(i);
68     pos = hit->GetPOS()-1;
69     tubeHits[pos].push_back(hit);
70     }
71     fpHits.clear();// we don't need it
72    
73     //From this moment all hits sorted by time and present in tubes
74    
75     for(Int_t i = 0; i<NTUBES; i++){
76     if( tubeHits[i].size() > 0 ){
77     for(UInt_t j = 0; j<tubeHits[i].size(); j++){
78     Uback = 0.;
79     hit =(PamVMCDetectorHit*)tubeHits[i].at(j);
80     U = keVtomV*(hit->GetEREL())*1.e6; //initial amplitude of this signal;
81     Time = hit->GetTOF();
82     for(UInt_t k = 1; k<=j; k++){
83     hit0 = (PamVMCDetectorHit*)tubeHits[i].at(k-1);
84     //hit1 = (PamVMCDetectorHit*)tubeHits[i].at(j);
85     if (tauHits[i].at(k-1)> 0)
86     Uback += keVtomV*(hit0->GetEREL())*1.e6*Exp((hit0->GetTOF()-/*hit1->GetTOF()*/hit->GetTOF())/Abs(tauHits[i].at(k-1)));
87     }
88    
89 pam-rm2 1.5 tauHits[i].push_back(frandom->Gaus(Tau,SigmaTau));
90 nikolas 1.1 if (DEBUG){
91     cout<<"++++++HIT+++++++"<<endl;
92     cout<<"TubeNo="<<hit->GetPOS()<<endl
93     <<"Particle="<<((TParticlePDG*)(TDatabasePDG::Instance()->GetParticle(hit->GetPDG())))->GetName()<<endl
94     <<" Erel="<<hit->GetEREL()*1.0e6<<" (keV)"<<endl
95     <<" Time="<<hit->GetTOF()*1.e6<<" (us)"<<endl
96     <<" Patn="<<hit->GetPATH()*10.<<" (mm)"<<endl
97     <<" Tau="<<tauHits[i].at(j)*1.e6<<" (us)"<<endl
98     <<" U="<<U<<" (mV)"<<endl
99     <<" Uback="<<Uback<<" (mV)"<<endl;
100     if ((Uback < Tresh) && (U+Uback > Tresh) && (Time<Showercut+Window) && (Time>Showercut))
101     cout<<" Particle detected as NEUTRON"<<endl;
102     }
103     if ((Uback < Tresh) && (U+Uback > Tresh) && (Time<Showercut+Window) && (Time>Showercut)) NdN++;
104    
105     }
106     }
107     }
108    
109     if(DEBUG)
110     cout<<"==== Total nuber of detected neutrons for for event ==="<<NdN<<" ==="<<endl;
111    
112     }
113    
114    
115    
116     //Cleaning of everything...
117     for (Int_t i = 0; i<NTUBES; i++){
118     tauHits[i].clear();
119     tubeHits[i].clear();
120     }
121    
122     fData.clear();
123    
124     if (DEBUG) cout<<"REGISTERED: "<< NdN <<" neutrons"<<endl;
125     for(Int_t i=0;i<3;i++){
126     fData.push_back(0x0000);
127     fData.push_back(0x010F);
128     }
129     fData.at(0) = ((0xFF00 & (256*NdN)));
130     }
  ViewVC Help
Powered by ViewVC 1.1.23