CGEM BOSS 6.6.5.h
BESIII Offline Software System
Loading...
Searching...
No Matches
BesTofDigitizerEcV4_dbs Class Reference

#include <BesTofDigitizerEcV4_dbs.hh>

+ Inheritance diagram for BesTofDigitizerEcV4_dbs:

Classes

struct  adc_info
 

Public Member Functions

 BesTofDigitizerEcV4_dbs ()
 
 ~BesTofDigitizerEcV4_dbs ()
 
virtual void Digitize (ScintSingle *, BesTofDigitsCollection *)
 
G4double Smear_gaussian (G4double, G4double, G4double)
 
G4int Get_gapnumber (G4int partID, G4double z_position)
 
G4double Distance_for_avalanche (G4int gap, G4double z_position, G4int partID)
 
void simulate_avalanche (G4int actual_strip, vector< adc_info > *avalanche_info, double *avalanche_threshold_time, double *induced_charge)
 
double adcsignal_get_sigma (double alpha, double eta, double stepwidth)
 
long int adcsignal_get_n_electron (double s, double stepwidth, double alpha, double eta)
 
long int adcsignal_simulate_step (int n_elec, double stepwidth, double alpha, double eta, bool saturated)
 
double chargetotime (double)
 
double Get_NINO_pulselength_resolution (double induced_charge_fc)
 
double Get_NINO_leadingedge_resolution (double induced_charge_fc)
 
- Public Member Functions inherited from BesTofDigitizerV
 BesTofDigitizerV ()
 
 ~BesTofDigitizerV ()
 
void Initialize ()
 

Static Public Member Functions

static G4int Calculate_Readoutstrip_number (G4double x, G4double y, G4int partId_f, G4int module_mrpc_f)
 
static G4int Calculate_Readoutstrip_number_continuum (G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
 
static G4int Produce_unique_identifier (G4int module_mrpc_f, G4int readoutstripnumber_f)
 
static G4int Get_module_mrpc_from_unique_identifier (G4int unique_identifier_f)
 
static G4int Get_stripnumber_from_unique_identifier (G4int unique_identifier_f)
 
static G4double Calculate_strip_transition_time_1 (G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
 
static G4double Calculate_strip_transition_time_2 (G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
 
static G4double Calculate_resulting_phi (G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
 
static double GetTransitionTime_extrap_track (G4double x_mm, G4double y_mm, int partId_f, int module_mrpc_f, int my_strip)
 
static double GetTransitionTime_dbsmatch (double delta_t, int my_strip)
 
static double Average_transition_time_dbs (G4int my_module)
 

Additional Inherited Members

- Protected Attributes inherited from BesTofDigitizerV
BesTofDigitsCollectionm_besTofDigitsCollection
 
BesTofHitsCollectionm_THC
 
ITofCaliSvcm_tofCaliSvc
 
ITofSimSvcm_tofSimSvc
 
ITofQElecSvcm_tofQElecSvc
 
G4double m_ADC [2]
 
G4double m_TDC [2]
 
G4int m_trackIndex
 
G4double m_globalTime
 
- Static Protected Attributes inherited from BesTofDigitizerV
static bool m_booked = false
 
static NTuple::Tuple * m_tupleTof1 = 0
 
static NTuple::Item< double > m_partId
 
static NTuple::Item< double > m_scinNb
 
static NTuple::Item< double > m_edep
 
static NTuple::Item< double > m_nHits
 
static NTuple::Item< double > m_time1st0
 
static NTuple::Item< double > m_time1st1
 
static NTuple::Item< double > m_timelast0
 
static NTuple::Item< double > m_timelast1
 
static NTuple::Item< double > m_totalPhot0
 
static NTuple::Item< double > m_totalPhot1
 
static NTuple::Item< double > m_NphAllSteps
 
static NTuple::Item< double > m_max0
 
static NTuple::Item< double > m_max1
 
static NTuple::Item< double > m_tdc0
 
static NTuple::Item< double > m_adc0
 
static NTuple::Item< double > m_tdc1
 
static NTuple::Item< double > m_adc1
 
static NTuple::Tuple * m_tupleTof2 = 0
 
static NTuple::Item< double > m_eTotal
 
static NTuple::Item< double > m_nDigi
 
static NTuple::Item< double > m_partIdMPV
 
static NTuple::Item< double > m_scinNbMPV
 
static NTuple::Item< double > m_edepMPV
 
static NTuple::Item< double > m_nDigiOut
 
static NTuple::Tuple * m_tupleTof3 = 0
 
static NTuple::Item< double > m_forb
 
static NTuple::Item< double > m_timeFlight
 
static NTuple::Item< double > m_ddT
 
static NTuple::Item< double > m_scinSwim
 
static NTuple::Item< double > m_scinTime
 
static NTuple::Item< double > m_transitTime
 
static NTuple::Item< double > m_endTime
 
static NTuple::Item< double > m_edepHit
 

Detailed Description

Definition at line 25 of file BesTofDigitizerEcV4_dbs.hh.

Constructor & Destructor Documentation

◆ BesTofDigitizerEcV4_dbs()

BesTofDigitizerEcV4_dbs::BesTofDigitizerEcV4_dbs ( )

Definition at line 30 of file BesTofDigitizerEcV4_dbs.cc.

31{
32
33
34 m_storedata=false;
35
36 if(m_storedata)
37 {
38 FileName = "output_mrpc_digitizerecv4.root";
39 m_TreeFile=new TFile(TString(FileName),"RECREATE","TreeFile");
40 m_ecv4tree = new TTree("ecv4tree", ""); //Create Tree
41 m_ecv4tree->Branch("m_partId",&m_partID,"m_partId/D");
42 m_ecv4tree->Branch("m_stripidentifier",&m_stripidentifier,"m_stripidentifier/D");
43 m_ecv4tree->Branch("m_trackIndex",&m_trackIndex,"m_trackIndex/D");
44 m_ecv4tree->Branch("m_signal_pc",&m_signal_pc,"m_signal_pc/D");
45 m_ecv4tree->Branch("m_time_threshold",&m_time_threshold,"m_time_threshold/D");
46 m_ecv4tree->Branch("m_time_1sthit",&m_time_1sthit,"m_time_1sthit/D");
47 m_ecv4tree->Branch("m_time_1",&m_time_1,"m_time_1/D");
48 m_ecv4tree->Branch("m_time_2",&m_time_2,"m_time_2/D");
49 m_ecv4tree->Branch("m_trans_time_1",&m_trans_time_1,"m_trans_time_1/D");
50 m_ecv4tree->Branch("m_trans_time_2",&m_trans_time_2,"m_trans_time_2/D");
51 m_ecv4tree->Branch("m_firedstrip",&m_firedstrip,"m_firedstrip/D");
52 m_ecv4tree->Branch("m_numberions",&m_numberions,"m_numberions/D");
53 m_ecv4tree->Branch("m_numberofgaps_with_avalanches",&m_numberofgaps_with_avalanches,"m_numberofgaps_with_avalanches/D");
54 m_ecv4tree->Branch("m_numberofstripsfired",&m_numberofstripsfired,"m_numberofstripsfired/D");
55 m_ecv4tree->Branch("m_multihit",&m_multihit,"m_multihit/D");
56 m_ecv4tree->Branch("m_inefficiency",&m_inefficiency,"m_inefficiency/D");
57 m_ecv4tree->Branch("m_particle_true",&m_particle_true,"m_particle_true/D");
58 m_ecv4tree->Branch("m_particle_has_signal",&m_particle_has_signal,"m_particle_has_signal/D");
59
60 }
61
62}

◆ ~BesTofDigitizerEcV4_dbs()

BesTofDigitizerEcV4_dbs::~BesTofDigitizerEcV4_dbs ( )

Definition at line 65 of file BesTofDigitizerEcV4_dbs.cc.

66{
67
68 if(m_storedata)
69 {
70 m_TreeFile->Write();
71 m_TreeFile->Close();
72 }
73
74}

Member Function Documentation

◆ adcsignal_get_n_electron()

long int BesTofDigitizerEcV4_dbs::adcsignal_get_n_electron ( double s,
double stepwidth,
double alpha,
double eta )

Definition at line 1288 of file BesTofDigitizerEcV4_dbs.cc.

1289{
1290
1291 double nbar = exp((alpha-eta)*stepwidth);
1292 double s_border= eta/alpha*(nbar-1.)/(nbar-eta/alpha);
1293 double k = eta/alpha;
1294
1295
1296 if(s<s_border)
1297 {
1298 return 0;
1299
1300 }
1301 else
1302 {
1303 return (long int) (1. + (1/log(((nbar-1)/(nbar-k))) * log(((nbar-k)*(s-1)/(k-1)/nbar))));
1304 }
1305
1306
1307}//close get_n_electron
EvtComplex exp(const EvtComplex &c)
const double alpha
XmlRpcServer s

Referenced by adcsignal_simulate_step().

◆ adcsignal_get_sigma()

double BesTofDigitizerEcV4_dbs::adcsignal_get_sigma ( double alpha,
double eta,
double stepwidth )

Definition at line 1310 of file BesTofDigitizerEcV4_dbs.cc.

1311{
1312
1313 double k=eta/alpha;
1314 double nbar = exp((alpha-eta)*stepwidth);
1315
1316 return (sqrt( (1.+k)/(1.-k) * nbar * (nbar-1.) )) ;
1317
1318}//close get_sigma

Referenced by adcsignal_simulate_step().

◆ adcsignal_simulate_step()

long int BesTofDigitizerEcV4_dbs::adcsignal_simulate_step ( int n_elec,
double stepwidth,
double alpha,
double eta,
bool saturated )

Definition at line 1244 of file BesTofDigitizerEcV4_dbs.cc.

1245{
1246
1247 if(n_elec<=150 && (saturated == false))
1248 {
1249
1250 long int n_elec_step =0;
1251
1252 for(int j=0;j<n_elec;j++)
1253 {
1254 double s=G4UniformRand();
1255
1256 n_elec_step = n_elec_step + adcsignal_get_n_electron(s, stepwidth, alpha, eta);
1257
1258 }
1259 return n_elec_step;
1260 }
1261
1262 if(n_elec>150 && saturated == false)
1263 {
1264 double sigma=adcsignal_get_sigma(alpha,eta,stepwidth);
1265 double nbar = exp((alpha-eta)*stepwidth);
1266 double my_sigma = G4RandGauss::shoot(0,(sqrt(n_elec)*sigma));
1267
1268 return (long int)(n_elec*nbar+ my_sigma);
1269 }
1270
1271 if(saturated==true)
1272 {
1273
1274 return n_elec;
1275 }//close if(is_saturated==true)
1276
1277}//closesimulate step
long int adcsignal_get_n_electron(double s, double stepwidth, double alpha, double eta)
double adcsignal_get_sigma(double alpha, double eta, double stepwidth)

Referenced by simulate_avalanche().

◆ Average_transition_time_dbs()

G4double BesTofDigitizerEcV4_dbs::Average_transition_time_dbs ( G4int my_module)
static

Definition at line 849 of file BesTofDigitizerEcV4_dbs.cc.

850{
851 G4double signal_velocity_copper_mrpc_strip= 0.8*0.299792458; // unit: percentage * c_0 * mm/ps
852
853 G4double length_to_border_of_the_readoutstrip=0.;
854
855 if(my_module==1 || my_module==2) length_to_border_of_the_readoutstrip=44;
856 if(my_module==3 || my_module==4) length_to_border_of_the_readoutstrip=46;
857 if(my_module==5 || my_module==6) length_to_border_of_the_readoutstrip=49;
858 if(my_module==7 || my_module==8) length_to_border_of_the_readoutstrip=51;
859 if(my_module==9 || my_module==10) length_to_border_of_the_readoutstrip=54;
860 if(my_module==11 || my_module==12) length_to_border_of_the_readoutstrip=56;
861 if(my_module==13 || my_module==14) length_to_border_of_the_readoutstrip=59;
862 if(my_module==15 || my_module==16) length_to_border_of_the_readoutstrip=61;
863 if(my_module==17 || my_module==18) length_to_border_of_the_readoutstrip=64;
864 if(my_module==19 || my_module==20) length_to_border_of_the_readoutstrip=66;
865 if(my_module==21 || my_module==22) length_to_border_of_the_readoutstrip=69;
866 if(my_module==23 || my_module==24) length_to_border_of_the_readoutstrip=71;
867
868 double transition_time = length_to_border_of_the_readoutstrip/signal_velocity_copper_mrpc_strip/1000.;
869 return transition_time;
870}

Referenced by MrpcDBSTrack::setCalibration().

◆ Calculate_Readoutstrip_number()

G4int BesTofDigitizerEcV4_dbs::Calculate_Readoutstrip_number ( G4double x,
G4double y,
G4int partId_f,
G4int module_mrpc_f )
static

Definition at line 400 of file BesTofDigitizerEcV4_dbs.cc.

402{
403
404//Within this functions all numbers are in mm or in degree!!!!
405
406 G4int readoutstripnumber=0;
407
408 G4double result_radius= sqrt(x_mm*x_mm+y_mm*y_mm);
409 G4double phi= atan2(y_mm,x_mm)*180./pi;
410
411
412 if(-180<phi && phi<90){phi=phi+270;} // Enter a coordinate system, which has on its highest point the 0 degree, and increases to the left/right side! The increasing depends on the endcap
413 else {phi=phi-90;}
414
415
416 G4double result_phi=0;
417 if(partId_f!=3 && partId_f!=6) result_phi=0;
418 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
419
420
421
422 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20 - phi + result_phi);
423 else result_phi = (module_mrpc_f-1)*20 - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
424
425 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
426 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
427 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
428
429
430
431 result_radius= result_radius*cos(result_phi*pi/180.0) - 634.56 + 159.5 - 3.82; //634.56 is the origin of coordinates from the module, 159.5 is the distance from the transformation of the readout strips...
432 // 3.82 is the ks_trafo
433
434 double result_x=fabs((result_radius/cos(result_phi*pi/180.0)+634.56-159.5+3.82)*sin(result_phi*pi/180)); //Calculate the resulting position x!
435
436
437
438 if( -12.5 <result_radius && result_radius < 12.5 && result_x <= 44 ) readoutstripnumber=1;
439 if( (-12.5 +(25.0+4.0)*1 ) <result_radius && result_radius < (12.5+(25.0+4.0)*1) && result_x <=46) readoutstripnumber=2;
440 if( (-12.5 +(25.0+4.0)*2 ) <result_radius && result_radius < (12.5+(25.0+4.0)*2) && result_x <=49) readoutstripnumber=3;
441 if( (-12.5 +(25.0+4.0)*3 ) <result_radius && result_radius < (12.5+(25.0+4.0)*3) && result_x <=51) readoutstripnumber=4;
442 if( (-12.5 +(25.0+4.0)*4 ) <result_radius && result_radius < (12.5+(25.0+4.0)*4) && result_x <=54) readoutstripnumber=5;
443 if( (-12.5 +(25.0+4.0)*5 ) <result_radius && result_radius < (12.5+(25.0+4.0)*5) && result_x <=56) readoutstripnumber=6;
444 if( (-12.5 +(25.0+4.0)*6 ) <result_radius && result_radius < (12.5+(25.0+4.0)*6) && result_x <=59) readoutstripnumber=7;
445 if( (-12.5 +(25.0+4.0)*7 ) <result_radius && result_radius < (12.5+(25.0+4.0)*7) && result_x <=61) readoutstripnumber=8;
446 if( (-12.5 +(25.0+4.0)*8 ) <result_radius && result_radius < (12.5+(25.0+4.0)*8) && result_x <=64) readoutstripnumber=9;
447 if( (-12.5 +(25.0+4.0)*9 ) <result_radius && result_radius < (12.5+(25.0+4.0)*9) && result_x <=66) readoutstripnumber=10;
448 if( (-12.5 +(25.0+4.0)*10 ) <result_radius && result_radius < (12.5+(25.0+4.0)*10) && result_x <=69) readoutstripnumber=11;
449 if( (-12.5 +(25.0+4.0)*11 ) <result_radius && result_radius < (12.5+(25.0+4.0)*11) && result_x <=71) readoutstripnumber=12;
450
451
452 return readoutstripnumber;
453
454}
double sin(const BesAngle a)
Definition BesAngle.h:210
double cos(const BesAngle a)
Definition BesAngle.h:213
const float pi
Definition vector3.h:133

Referenced by Digitize().

◆ Calculate_Readoutstrip_number_continuum()

G4int BesTofDigitizerEcV4_dbs::Calculate_Readoutstrip_number_continuum ( G4double x_mm,
G4double y_mm,
G4int partId_f,
G4int module_mrpc_f )
static

Definition at line 461 of file BesTofDigitizerEcV4_dbs.cc.

462{
463
464 //Within this functions all numbers are in mm or in degree!!!!
465
466 G4int readoutstripnumber=0;
467
468 G4double result_radius= sqrt(x_mm*x_mm+y_mm*y_mm);
469 G4double phi= atan2(y_mm,x_mm)*180./pi;
470
471
472 if(-180<phi && phi<90){phi=phi+270;}//See funtion above!
473 else {phi=phi-90;}
474
475
476 G4double result_phi=0;
477 if(partId_f!=3 && partId_f!=6) result_phi=0;
478 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
479 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20 - phi + result_phi);
480 else result_phi = (module_mrpc_f-1)*20 - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
481
482 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
483 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
484 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
485
486
487
488 result_radius= result_radius*cos(result_phi*pi/180.0) - 634.56 + 159.5-3.82; //634.56 is the origin of coordinates from the module, 159.5 is the distance from the transformation of the readout strips...
489
490 //The +-12.0 are chosen arbitrarily. In this region we assume reconstructed tracks as MRPC tracks!
491 if( (-12.5 -12.0) <=result_radius && result_radius < 12.5+2.0 ) readoutstripnumber=1;
492 if( (-12.5 +(25.0+4.0)*1-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*1)+2.0 ) readoutstripnumber=2;
493 if( (-12.5 +(25.0+4.0)*2-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*2) +2.0) readoutstripnumber=3;
494 if( (-12.5 +(25.0+4.0)*3-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*3) +2.0) readoutstripnumber=4;
495 if( (-12.5 +(25.0+4.0)*4-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*4) +2.0) readoutstripnumber=5;
496 if( (-12.5 +(25.0+4.0)*5-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*5) +2.0) readoutstripnumber=6;
497 if( (-12.5 +(25.0+4.0)*6-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*6) +2.0) readoutstripnumber=7;
498 if( (-12.5 +(25.0+4.0)*7-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*7) +2.0) readoutstripnumber=8;
499 if( (-12.5 +(25.0+4.0)*8-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*8) +2.0) readoutstripnumber=9;
500 if( (-12.5 +(25.0+4.0)*9-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*9) +2.0) readoutstripnumber=10;
501 if( (-12.5 +(25.0+4.0)*10-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*10) +2.0) readoutstripnumber=11;
502 if( (-12.5 +(25.0+4.0)*11-2.0 ) <=result_radius && result_radius < (12.5+(25.0+4.0)*11) +12.0) readoutstripnumber=12;
503
504
505 return readoutstripnumber;
506
507}

Referenced by ExtSteppingAction::UserSteppingAction().

◆ Calculate_resulting_phi()

G4double BesTofDigitizerEcV4_dbs::Calculate_resulting_phi ( G4double x_mm,
G4double y_mm,
G4int partId_f,
G4int module_mrpc_f )
static

Definition at line 368 of file BesTofDigitizerEcV4_dbs.cc.

369{
370
371
372
373
374 G4double phi= atan2(y_mm,x_mm)*180./pi;
375
376
377 if(-180<phi && phi<90){phi=phi+270;} // Enter a coordinate system, which has on its highest point the 0 degree, and increases to the left/right side! The increasing depends on the endcap
378 else {phi=phi-90;}
379
380
381 G4double result_phi=0;
382 if(partId_f!=3 && partId_f!=6) result_phi=0;
383 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
384
385 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20 - phi + result_phi);
386 else result_phi = (module_mrpc_f-1)*20 - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
387
388 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
389 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
390 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
391
392 return result_phi;
393 // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
394
395}

Referenced by Digitize().

◆ Calculate_strip_transition_time_1()

G4double BesTofDigitizerEcV4_dbs::Calculate_strip_transition_time_1 ( G4double x_mm,
G4double y_mm,
G4int partId_f,
G4int module_mrpc_f )
static

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Definition at line 516 of file BesTofDigitizerEcV4_dbs.cc.

517{
518 G4double signal_velocity_copper_mrpc_strip= 0.8*0.299792458; // unit: percentage * c_0 * mm/ps
519
520 //Within this functions all numbers are in mm or in degree!!!!
521 G4double length_to_border_of_the_readoutstrip=0;
522 G4double result_radius= sqrt(x_mm*x_mm+y_mm*y_mm); //Radius where the particle did hit the detector
523 G4double phi= atan2(y_mm,x_mm)*180./pi; //Phi angle of the particle hitting the detector
524
525 if(-180<phi && phi<90){phi=phi+270;} // Transformation into a coordinate-system, whose coordinates do increase to the right/left side. Left or right does depend on the endcap!
526 else {phi=phi-90;}
527
528 G4double result_phi=0;
529 if(partId_f!=3 && partId_f!=6) result_phi=0;
530 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
531
532
533 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20 - phi + result_phi);
534 else result_phi = (module_mrpc_f-1)*20 - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
535 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
536
537 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
538 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
539
540 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
541 //The Phi is now the is the phi angle for each module and not for the whole assembly
542 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
543
544
545 result_radius= result_radius*cos(result_phi*pi/180.0) - 634.56 + 159.5 - 3.82; //634.56 is the origin of coordinates from the module, 159.5 is the distance from the transformation of the readout strips...
546 // 3.82 is the ks_trafo
547
548 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
549 //The radius is now the the radius within one module. The point of origin for the
550 //radius is the middle of the lowest readout cell
551 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
552
553
554 double result_x=fabs((result_radius/cos(result_phi*pi/180.0)+634.56-159.5+3.82)*sin(result_phi*pi/180)); //Calculate the resulting position x!
555
556
557 if( -12.5 <result_radius && result_radius < 12.5 && result_x <= 44 )
558 {
559 length_to_border_of_the_readoutstrip=sqrt( (44-result_x)*(44-result_x) + fabs((result_radius-0.)*(result_radius-0.)));
560 }
561 if( (-12.5 +(25.0+4.0)*1 ) <result_radius && result_radius < (12.5+(25.0+4.0)*1) && result_x <=46)
562 {
563 length_to_border_of_the_readoutstrip=sqrt( (46-result_x)*(46-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*1)));
564 }
565 if( (-12.5 +(25.0+4.0)*2 ) <result_radius && result_radius < (12.5+(25.0+4.0)*2) && result_x <=49)
566 {
567 length_to_border_of_the_readoutstrip=sqrt( (49-result_x)*(49-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*2)));
568 }
569 if( (-12.5 +(25.0+4.0)*3 ) <result_radius && result_radius < (12.5+(25.0+4.0)*3) && result_x <=51)
570 {
571 length_to_border_of_the_readoutstrip=sqrt( (51-result_x)*(51-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*3)));
572 }
573 if( (-12.5 +(25.0+4.0)*4 ) <result_radius && result_radius < (12.5+(25.0+4.0)*4) && result_x <=54)
574 {
575 length_to_border_of_the_readoutstrip=sqrt( (54-result_x)*(54-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*4)));
576 }
577 if( (-12.5 +(25.0+4.0)*5 ) <result_radius && result_radius < (12.5+(25.0+4.0)*5) && result_x <=56)
578 {
579 length_to_border_of_the_readoutstrip=sqrt( (56-result_x)*(56-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*5)));
580 }
581 if( (-12.5 +(25.0+4.0)*6 ) <result_radius && result_radius < (12.5+(25.0+4.0)*6) && result_x <=59)
582 {
583 length_to_border_of_the_readoutstrip=sqrt( (59-result_x)*(59-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*6)));
584 }
585 if( (-12.5 +(25.0+4.0)*7 ) <result_radius && result_radius < (12.5+(25.0+4.0)*7) && result_x <=61)
586 {
587 length_to_border_of_the_readoutstrip=sqrt( (61-result_x)*(61-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*7)));
588 }
589 if( (-12.5 +(25.0+4.0)*8 ) <result_radius && result_radius < (12.5+(25.0+4.0)*8) && result_x <=64)
590 {
591 length_to_border_of_the_readoutstrip=sqrt( (64-result_x)*(64-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*8)));
592 }
593 if( (-12.5 +(25.0+4.0)*9 ) <result_radius && result_radius < (12.5+(25.0+4.0)*9) && result_x <=66)
594 {
595 length_to_border_of_the_readoutstrip=sqrt( (66-result_x)*(66-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*9)));
596 }
597 if( (-12.5 +(25.0+4.0)*10 ) <result_radius && result_radius < (12.5+(25.0+4.0)*10) && result_x <=69)
598 {
599 length_to_border_of_the_readoutstrip=sqrt( (69-result_x)*(69-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*10)));
600 }
601 if( (-12.5 +(25.0+4.0)*11 ) <result_radius && result_radius < (12.5+(25.0+4.0)*11) && result_x <=71)
602 {
603 length_to_border_of_the_readoutstrip=sqrt( (71-result_x)*(71-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*11)));
604 }
605 double transition_time = length_to_border_of_the_readoutstrip/signal_velocity_copper_mrpc_strip/1000.; //length is in mm , 3*10^8 m/s = 0.3 mm/ps --> Unit of transitionstime is ns (due to division by 1000.)
606 return transition_time;
607
608}//close calculate_strip_transition_time_1

Referenced by Digitize().

◆ Calculate_strip_transition_time_2()

G4double BesTofDigitizerEcV4_dbs::Calculate_strip_transition_time_2 ( G4double x_mm,
G4double y_mm,
G4int partId_f,
G4int module_mrpc_f )
static

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

Definition at line 613 of file BesTofDigitizerEcV4_dbs.cc.

614{
615 G4double signal_velocity_copper_mrpc_strip= 0.8*0.299792458; // unit: percentage * c_0 * mm/ps
616
617 //Within this functions all numbers are in mm or in degree!!!!
618 G4double length_to_border_of_the_readoutstrip=0;
619 G4double result_radius= sqrt(x_mm*x_mm+y_mm*y_mm); //Radius where the particle did hit the detector
620 G4double phi= atan2(y_mm,x_mm)*180./pi; //Phi angle of the particle hitting the detector
621
622 if(-180<phi && phi<90){phi=phi+270;} // Transformation into a coordinate-system, whose coordinates do increase to the right/left side. Left or right does depend on the endcap!
623 else {phi=phi-90;}
624
625 G4double result_phi=0;
626 if(partId_f!=3 && partId_f!=6) result_phi=0;
627 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
628
629
630 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20 - phi + result_phi);
631 else result_phi = (module_mrpc_f-1)*20 - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
632 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
633
634 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
635 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
636
637 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
638 //The Phi is now the is the phi angle for each module and not for the whole assembly
639 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
640
641
642 result_radius= result_radius*cos(result_phi*pi/180.0) - 634.56 + 159.5 - 3.82; //634.56 is the origin of coordinates from the module, 159.5 is the distance from the transformation of the readout strips...
643 // 3.82 is the ks_trafo
644
645 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
646 //The radius is now the the radius within one module. The point of origin for the
647 //radius is the middle of the lowest readout cell
648 //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
649
650
651 double result_x=fabs((result_radius/cos(result_phi*pi/180.0)+634.56-159.5+3.82)*sin(result_phi*pi/180)); //Calculate the resulting position x!
652
653
654 if( -12.5 <result_radius && result_radius < 12.5 && result_x <= 44 )
655 {
656 length_to_border_of_the_readoutstrip=sqrt( (48+result_x)*(48+result_x) + fabs((result_radius-0.)*(result_radius-0.)));
657 }
658 if( (-12.5 +(25.0+4.0)*1 ) <result_radius && result_radius < (12.5+(25.0+4.0)*1) && result_x <=46)
659 {
660 length_to_border_of_the_readoutstrip=sqrt( (50+result_x)*(50+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*1)));
661 }
662 if( (-12.5 +(25.0+4.0)*2 ) <result_radius && result_radius < (12.5+(25.0+4.0)*2) && result_x <=49)
663 {
664 length_to_border_of_the_readoutstrip=sqrt( (53+result_x)*(53+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*2)));
665 }
666 if( (-12.5 +(25.0+4.0)*3 ) <result_radius && result_radius < (12.5+(25.0+4.0)*3) && result_x <=51)
667 {
668 length_to_border_of_the_readoutstrip=sqrt( (55+result_x)*(55+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*3)));
669 }
670 if( (-12.5 +(25.0+4.0)*4 ) <result_radius && result_radius < (12.5+(25.0+4.0)*4) && result_x <=54)
671 {
672 length_to_border_of_the_readoutstrip=sqrt( (58+result_x)*(58+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*4)));
673 }
674 if( (-12.5 +(25.0+4.0)*5 ) <result_radius && result_radius < (12.5+(25.0+4.0)*5) && result_x <=56)
675 {
676 length_to_border_of_the_readoutstrip=sqrt( (60+result_x)*(60+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*5)));
677 }
678 if( (-12.5 +(25.0+4.0)*6 ) <result_radius && result_radius < (12.5+(25.0+4.0)*6) && result_x <=59)
679 {
680 length_to_border_of_the_readoutstrip=sqrt( (63+result_x)*(63+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*6)));
681 }
682 if( (-12.5 +(25.0+4.0)*7 ) <result_radius && result_radius < (12.5+(25.0+4.0)*7) && result_x <=61)
683 {
684 length_to_border_of_the_readoutstrip=sqrt( (65+result_x)*(65+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*7)));
685 }
686 if( (-12.5 +(25.0+4.0)*8 ) <result_radius && result_radius < (12.5+(25.0+4.0)*8) && result_x <=64)
687 {
688 length_to_border_of_the_readoutstrip=sqrt( (68+result_x)*(68+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*8)));
689 }
690 if( (-12.5 +(25.0+4.0)*9 ) <result_radius && result_radius < (12.5+(25.0+4.0)*9) && result_x <=66)
691 {
692 length_to_border_of_the_readoutstrip=sqrt( (70+result_x)*(70+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*9)));
693 }
694 if( (-12.5 +(25.0+4.0)*10 ) <result_radius && result_radius < (12.5+(25.0+4.0)*10) && result_x <=69)
695 {
696 length_to_border_of_the_readoutstrip=sqrt( (73+result_x)*(73+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*10)));
697 }
698 if( (-12.5 +(25.0+4.0)*11 ) <result_radius && result_radius < (12.5+(25.0+4.0)*11) && result_x <=71)
699 {
700 length_to_border_of_the_readoutstrip=sqrt( (75+result_x)*(75+result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*11)));
701 }
702
703 double transition_time = length_to_border_of_the_readoutstrip/signal_velocity_copper_mrpc_strip/1000.; //length is in mm , 3*10^8 m/s = 0.3 mm/ps --> Unit of transitionstime is ns (due to division by 1000.)
704 return transition_time;
705
706}//close calculate_strip_transition_time_2

Referenced by Digitize().

◆ chargetotime()

double BesTofDigitizerEcV4_dbs::chargetotime ( double charge)

Definition at line 1325 of file BesTofDigitizerEcV4_dbs.cc.

1326{
1327
1328
1329 double time=0.;
1330 if(charge>=13.) time=(-4.50565/log(charge*0.0812208)+16.6653)*ns;
1331
1332 return time;
1333
1334}//close chargetotime
Double_t time
#define ns(x)
Definition xmltok.c:1504

Referenced by Digitize().

◆ Digitize()

void BesTofDigitizerEcV4_dbs::Digitize ( ScintSingle * single_module,
BesTofDigitsCollection * DC )
virtual

Reimplemented from BesTofDigitizerV.

Definition at line 76 of file BesTofDigitizerEcV4_dbs.cc.

77{
78
79 trackIndex=0;
80 partId=-999;
81 module_mrpc=-999;
82 time=-999.;
83
84
86 G4DigiManager* digiManager = G4DigiManager::GetDMpointer();
87 G4int THCID = digiManager->GetHitsCollectionID("BesTofHitsCollection");
88 m_THC = (BesTofHitsCollection*) (digiManager->GetHitsCollection(THCID));
89 if (m_THC)
90 {
91
92
93
94 G4int n_hit=single_module->GetHitIndexes_mrpc()->size();
95
96
97 vector<adc_info> avalanche_info;//Will store the basic information for signal creation etc.
98 avalanche_info.clear();
99
100 partId = single_module->GetPartId();
101 module_mrpc = single_module->GetModule_mrpc();
102
103
104
105 //Process the hits
106 int deadarea=0;//DeadArea Counter: How often can't we associate an readoutstrip?
107
108 for(G4int i=0;i<n_hit;i++)
109 {
110
111
112 BesTofHit* hit = (*m_THC)[( *(single_module->GetHitIndexes_mrpc()) )[i]];
113
114
115 G4ThreeVector help = hit->GetPos();
116 G4int readoutstripnumber = Calculate_Readoutstrip_number(help.x()/mm,help.y()/mm,partId,module_mrpc);
117
118 if( readoutstripnumber==0 ) {deadarea++; continue;}
119
120
121 adc_info struct_adc_info;
122 struct_adc_info.gap =0; struct_adc_info.readoutstrip =0; struct_adc_info.ions =0; struct_adc_info.zdistance =0.; struct_adc_info.trackindex=0;struct_adc_info.time=0.;
123 struct_adc_info.pdg_code=0;//Intialize all values
124
125 struct_adc_info.gap=Get_gapnumber(partId, help.z());
126 struct_adc_info.readoutstrip=readoutstripnumber;
127 struct_adc_info.ions=hit->GetIons();
128 struct_adc_info.zdistance=Distance_for_avalanche(struct_adc_info.gap,help.z(), partId);
129 struct_adc_info.time = hit->GetTime();
130 struct_adc_info.trackindex=hit->GetG4Index();
131 struct_adc_info.pdg_code=hit->GetPDGcode();
132 struct_adc_info.transtime_1=Calculate_strip_transition_time_1(help.x()/mm,help.y()/mm,partId,module_mrpc);
133 struct_adc_info.transtime_2=Calculate_strip_transition_time_2(help.x()/mm,help.y()/mm,partId,module_mrpc);
134 struct_adc_info.x=help.x()/mm;
135 struct_adc_info.y=help.y()/mm;
136 struct_adc_info.partId=partId;
137 struct_adc_info.module=module_mrpc;
138
139 if(struct_adc_info.ions>0 && struct_adc_info.readoutstrip!=0 && struct_adc_info.time!=0.) avalanche_info.push_back(struct_adc_info);
140
141 }//Close for-loop over hits
142
143
144
145 G4int avalanche_info_size = avalanche_info.size();
146
147
148 //Determine which readoutstrips has been fired and push the depeding on the strip number into one vector
149 vector <int> fired_readoutstrips;
150
151 fired_readoutstrips.clear();
152 for(G4int i=0; i<avalanche_info_size; i++)
153 {
154
155 bool dlh=true;
156 for(G4int j=0;j<fired_readoutstrips.size();j++)
157 {
158 if(fired_readoutstrips[j]==avalanche_info[i].readoutstrip){dlh=false; break;}
159 }
160 if(dlh==true)
161 {
162 fired_readoutstrips.push_back(avalanche_info[i].readoutstrip);
163 }
164 }
165
166
167
168 //Determine the adcinfo for each fired strips with the smallest time -> These information will be stored within the digis
169 vector <int> mybest_adcinfo_inreadoutstrips(fired_readoutstrips.size());
170 for(G4int i=0;i< fired_readoutstrips.size();i++)
171 {
172 double time_help=10000.;
173 for(G4int j=0; j<avalanche_info_size; j++)
174 {
175 if((time_help>avalanche_info[j].time) && avalanche_info[j].readoutstrip==fired_readoutstrips[i])
176 {
177 time_help = avalanche_info[j].time;
178 mybest_adcinfo_inreadoutstrips[i]=j;
179 }//close if
180 }//close for j
181 }//close for i
182
183
184
185
186 //Calculate the signal strength and determine the time information and store these information into the digi collection
187 double stripidentifier_1=0; //Each strip will be read out by 2 NINO chips, each on one side!
188 double stripidentifier_2=0;
189
190 for(G4int i=0; i<fired_readoutstrips.size();i++)
191 {
192
193
194
195 //Determine the total number of ions!
196 G4int sumions=0;
197 for(G4int j=0; j<avalanche_info_size; j++)
198 {
199 if(avalanche_info[j].readoutstrip==fired_readoutstrips[i])
200 {
201 sumions = sumions+ avalanche_info[j].ions;
202 }
203
204 }
205
206
207 //Determine the total number of gaps with avalanches
208 G4int numberofgaps_with_avalanches=0;
209 for(G4int j=1;j<=12;j++)
210 {
211 for(G4int k=0; k<avalanche_info_size; k++) { if(avalanche_info[k].gap ==j && avalanche_info[k].ions >0){numberofgaps_with_avalanches++; break;} }
212 }
213
214
215
216
217 G4double signal_pc=0; G4double time_of_threshold_exceedance=0;
218 simulate_avalanche(fired_readoutstrips[i], &avalanche_info,&time_of_threshold_exceedance,&signal_pc);
219 signal_pc=signal_pc;//1e-12*A*s;
220 time_of_threshold_exceedance=time_of_threshold_exceedance*s;//Stores the time when the threshold is reached.
221
222
223 if(signal_pc<0.015) continue; //We only store signals wich are over the threshold
224
225
226 trackIndex = avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].trackindex;
227
228
229
230
231 // We add the threshold crossing time for the avalanche. This will be corrected later within the reconstruction.
232 double resulting_phi = Calculate_resulting_phi(avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].x,avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].y
233 ,avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].partId,avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].module);
234
235 if(resulting_phi<0)//The ninos one the left side of a module (odd numbers) will see the signal earlier!
236 //Time 1 is the smaller transition time.
237 {
238 time = avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].time + time_of_threshold_exceedance + avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_1;
239 time_digi2=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].time + time_of_threshold_exceedance + avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_2;
240 }
241 else//The ninos one the right side of a module (even numbers) will see the signal earlier!
242 //Time 1 is the smaller transition time.
243 {
244 time = avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].time + time_of_threshold_exceedance + avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_2;
245 time_digi2=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].time + time_of_threshold_exceedance + avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_1;
246 }
247
248
249 //time smear #Correct
250 //Uncertainities which are not considered in the digitization process. Following values are assumed:
251 //Avalanche Process : Already considered with time of threshold!
252 //FEE (NINO) : Depending on the charge, will be calculated in the next step! = NINO_RES
253 //Hit position : 10 ps
254 //Readout : 27 ps
255 //Electronics/cable : 20 ps
256 //-----------------------------------------
257 // Sum =sqrt(10**2 + 27**2 +20**2 NINO_RES**2 ) time smear
258
259 double femtocolomb= 1000.;
260
261 double NINO_RES=0;
262 NINO_RES = Get_NINO_leadingedge_resolution(signal_pc*femtocolomb);
263 double smear_time = sqrt(NINO_RES*NINO_RES + 27.0*27.0 + 10.0*10.0 + 20.0*20.0); //This is in ps
264
265 time = Smear_gaussian(time,0,smear_time*0.001*ns);
266 time_digi2=Smear_gaussian(time_digi2,0,smear_time*0.001*ns);
267
268
269 stripidentifier_1= Produce_unique_identifier(module_mrpc,(avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].readoutstrip)*2-1);
270 stripidentifier_2= Produce_unique_identifier(module_mrpc,(avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].readoutstrip)*2);
271
272
273 //Resolution of Pulse length:
274 // 27 ps TDC
275 // Timeresolution of pulse length
276 double NINO_pulselength_resolution = Get_NINO_pulselength_resolution( signal_pc*femtocolomb );
277 double smear_time_2 = sqrt(NINO_pulselength_resolution*NINO_pulselength_resolution + 27.0*27.0 );
278
279 double time2 =chargetotime(signal_pc*femtocolomb)+time; //The function returns the time in ns
280 time2 = Smear_gaussian(time2,0,smear_time_2*0.001*ns);
281
282
283 //so far I do assume,that both Nino chips see the same signals!
284 double time2_digi2= chargetotime(signal_pc*femtocolomb)+time_digi2;
285 time2_digi2 = Smear_gaussian(time2_digi2,0,smear_time_2*0.001*ns);
286
287
288 if(m_storedata)
289 {
290 m_partID=partId;
291 m_stripidentifier=stripidentifier_1;
292 m_trackIndex=trackIndex;
293 m_signal_pc=signal_pc;
294 m_time_threshold= time_of_threshold_exceedance;
295 m_time_1sthit=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].time;
296 m_time_1=time;
297 m_time_2=time2;
298 m_trans_time_1=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_1;
299 m_trans_time_2=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].transtime_2;
300 m_firedstrip=fired_readoutstrips[i];
301 m_numberions=sumions;
302 m_numberofgaps_with_avalanches=numberofgaps_with_avalanches;
303 m_numberofstripsfired = fired_readoutstrips.size();
304 m_x=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].x;
305 m_y=avalanche_info[(mybest_adcinfo_inreadoutstrips[i])].y;
306 m_ecv4tree->Fill();
307 }
308
309
310
311
312
313 BesTofDigi *digi = new BesTofDigi;
314 digi->SetTrackIndex(trackIndex);
315 digi->SetPartId(partId);
316 digi->SetModule_mrpc(stripidentifier_1);
317 digi->SetTime1_mrpc(time);
318 digi->SetTime2_mrpc(time2);
319 m_besTofDigitsCollection->insert(digi);
320
321
322
323
324 BesTofDigi *digi2 = new BesTofDigi;
325 digi2->SetTrackIndex(trackIndex);
326 digi2->SetPartId(partId);
327 digi2->SetModule_mrpc(stripidentifier_2);
328 digi2->SetTime1_mrpc(time_digi2);
329 digi2->SetTime2_mrpc(time2_digi2);
330 m_besTofDigitsCollection->insert(digi2);
331
332
333
334 }//close i
335
336
337
338
339 }//close if(m_THC)
340
341} //Digitize()
G4THitsCollection< BesTofHit > BesTofHitsCollection
Definition BesTofHit.hh:108
std::string help()
void SetTime1_mrpc(G4double time_mrpc)
Definition BesTofDigi.hh:47
void SetPartId(G4int partId)
Definition BesTofDigi.hh:39
void SetTrackIndex(G4int index)
Definition BesTofDigi.hh:38
void SetTime2_mrpc(G4double time_mrpc2)
Definition BesTofDigi.hh:48
void SetModule_mrpc(G4int module_mrpc)
Definition BesTofDigi.hh:46
static G4int Calculate_Readoutstrip_number(G4double x, G4double y, G4int partId_f, G4int module_mrpc_f)
static G4double Calculate_resulting_phi(G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
static G4double Calculate_strip_transition_time_2(G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
G4int Get_gapnumber(G4int partID, G4double z_position)
static G4double Calculate_strip_transition_time_1(G4double x_mm, G4double y_mm, G4int partId_f, G4int module_mrpc_f)
void simulate_avalanche(G4int actual_strip, vector< adc_info > *avalanche_info, double *avalanche_threshold_time, double *induced_charge)
double Get_NINO_pulselength_resolution(double induced_charge_fc)
static G4int Produce_unique_identifier(G4int module_mrpc_f, G4int readoutstripnumber_f)
G4double Distance_for_avalanche(G4int gap, G4double z_position, G4int partID)
double Get_NINO_leadingedge_resolution(double induced_charge_fc)
G4double Smear_gaussian(G4double, G4double, G4double)
BesTofDigitsCollection * m_besTofDigitsCollection
G4int GetIons()
Definition BesTofHit.hh:80
G4int GetPDGcode()
Definition BesTofHit.hh:73
G4double GetTime()
Definition BesTofHit.hh:66
G4ThreeVector GetPos()
Definition BesTofHit.hh:65
G4int GetG4Index()
Definition BesTofHit.hh:59
vector< G4int > * GetHitIndexes_mrpc()
G4int GetModule_mrpc()
G4int GetPartId()

◆ Distance_for_avalanche()

G4double BesTofDigitizerEcV4_dbs::Distance_for_avalanche ( G4int gap,
G4double z_position,
G4int partID )

Definition at line 962 of file BesTofDigitizerEcV4_dbs.cc.

963{
964
965 //All the values are in mm
966 G4double pos_west_1=-1343.5*mm; // = partID 5
967 G4double pos_west_2=-1368.5*mm; // = partID 6
968 G4double pos_east_1=+1343.5*mm; // = partID 4
969 G4double pos_east_2=+1368.5*mm; // = partID 3
970
971
972 //We transform our z-coordinate into the system of one MRPC module. PoV= Interactionpoint:
973 //The gaps 1-6 will have negative values, the gaps 7-12 will have positive values!
974 G4double z_zero=0;
975 if(partID==3)
976 {z_zero= z_position-pos_east_2;}
977 else if(partID==4)
978 {z_zero= z_position-pos_east_1;}
979 else if(partID==5)
980 {z_zero = pos_west_1 - z_position ;}//z_position is negative in this case!!!!
981 else if(partID==6)
982 {z_zero = pos_west_2 - z_position ;}//z_position is negative in this case!!!!
983 else
984 {G4cout << "TofDigitzer: Error in part ID" <<G4endl; return 0;}
985
986 //cout << "z_zero "<< setprecision(100)<< z_zero << G4endl;
987
988 /*
989 if(gap==1) G4cout << (-z_zero-4.60*mm) << " "<< gap << G4endl ;
990 if(gap==2) G4cout <<(-z_zero-3.98*mm)<< " "<< gap << G4endl;
991 if(gap==3) G4cout<< (-z_zero-3.36*mm)<< " "<< gap << G4endl;
992 if(gap==4) G4cout<< (-z_zero-2.74*mm)<< " "<< gap << G4endl;
993 if(gap==5) G4cout<< (-z_zero-2.12*mm)<< " "<< gap << G4endl;
994 if(gap==6) G4cout<< (-z_zero-1.5*mm)<< " "<< gap << G4endl;
995 if(gap==7) G4cout<< (1.72-z_zero*mm)<< " "<< gap << G4endl;
996 if(gap==8) G4cout<< (2.34-z_zero*mm)<< " "<< gap << G4endl;
997 if(gap==9) G4cout<< (2.96-z_zero*mm)<< " "<< gap << G4endl;
998 if(gap==10)G4cout<< (3.58-z_zero*mm)<< " "<< gap << G4endl;
999 if(gap==11)G4cout<< (4.20-z_zero*mm)<< " "<< gap << G4endl;
1000 if(gap==12)G4cout<< (4.82-z_zero*mm)<< " "<< gap << G4endl;
1001 */
1002
1003
1004
1005
1006
1007 if(gap==1) return (-z_zero-4.60*mm);
1008 if(gap==2) return (-z_zero-3.98*mm);
1009 if(gap==3) return (-z_zero-3.36*mm);
1010 if(gap==4) return (-z_zero-2.74*mm);
1011 if(gap==5) return (-z_zero-2.12*mm);
1012 if(gap==6) return (-z_zero-1.5*mm);
1013 if(gap==7) return (1.72-z_zero*mm);
1014 if(gap==8) return (2.34-z_zero*mm);
1015 if(gap==9) return (2.96-z_zero*mm);
1016 if(gap==10)return (3.58-z_zero*mm);
1017 if(gap==11)return (4.20-z_zero*mm);
1018 if(gap==12)return (4.82-z_zero*mm);
1019
1020
1021}

Referenced by Digitize().

◆ Get_gapnumber()

G4int BesTofDigitizerEcV4_dbs::Get_gapnumber ( G4int partID,
G4double z_position )

Definition at line 914 of file BesTofDigitizerEcV4_dbs.cc.

915{
916
917
918
919 G4double pos_west_1=-1343.5*mm; // = partID 5
920 G4double pos_west_2=-1368.5*mm; // = partID 6
921 G4double pos_east_1=+1343.5*mm; // = partID 4
922 G4double pos_east_2=+1368.5*mm; // = partID 3
923
924
925 //We transform our z-coordinate into the system of one MRPC module. PoV= Interactionpoint:
926 //The gaps 1-6 will have negative values, the gaps 7-12 will have positive values!
927 G4double z_zero=0;
928 if(partID==3)
929 {z_zero= z_position-pos_east_2;}
930 else if(partID==4)
931 {z_zero= z_position-pos_east_1;}
932 else if(partID==5)
933 {z_zero = pos_west_1 - z_position ;}//z_position is negative in this case!!!!
934 else if(partID==6)
935 {z_zero = pos_west_2 - z_position ;}//z_position is negative in this case!!!!
936 else
937 {G4cout << "TofDigitzer: Error in part ID" <<G4endl; return 0;}
938
939 //std::cout << "z_pos " << z_position << std::endl;
940 //std::cout << "z_zero " << std::setprecision(100) <<z_zero << std::endl;
941
942
943 //There is a porblem with the numbers. The correct ones are the rounded ones. However in this case too often 0 is returned!
944 if((-4.5999*mm) >= z_zero && z_zero>(-4.8201*mm) ) {return 1;}
945 else if((-3.9799*mm) >= z_zero && z_zero>(-4.2001*mm) ){return 2;}
946 else if((-3.3599*mm) >= z_zero && z_zero>(-3.5801*mm) ){return 3;}
947 else if((-2.7399*mm) >= z_zero && z_zero>(-2.9601*mm) ){return 4;}
948 else if((-2.1199*mm) >= z_zero && z_zero>(-2.3401*mm) ){return 5;}
949 else if((-1.4999*mm) >= z_zero && z_zero>(-1.7201*mm)) {return 6;}
950 else if((1.4999*mm) <= z_zero && z_zero<=(1.7201*mm) ) {return 7;}
951 else if((2.1199*mm) <= z_zero && z_zero<=(2.3401*mm) ) {return 8;}
952 else if((2.7399*mm) <= z_zero && z_zero<=(2.9601*mm) ) {return 9;}
953 else if((3.3599*mm) <= z_zero && z_zero<=(3.5801*mm) ) {return 10;}
954 else if((3.9799*mm) <= z_zero && z_zero<=(4.2001*mm) ) {return 11;}
955 else if((4.5999*mm) <= z_zero && z_zero<=(4.8201*mm) ){return 12;}
956 else {G4cout << "TofDigitzer: Error in Hit-Position! " <<G4endl; return 0; }
957
958}

Referenced by Digitize().

◆ Get_module_mrpc_from_unique_identifier()

G4int BesTofDigitizerEcV4_dbs::Get_module_mrpc_from_unique_identifier ( G4int unique_identifier_f)
static

Definition at line 889 of file BesTofDigitizerEcV4_dbs.cc.

890{
891
892 G4int returnvalue=0;
893
894 returnvalue= unique_identifier_f/25;
895
896
897 return returnvalue;
898}

Referenced by MrpcDBSRec::execute(), MrpcDBSTrack::findTofDataMRPC(), MrpcDBSTrack::get_mrpc_module_id_neighbour(), MrpcDBSTrack::setCalibration(), MrpcDBSTrack::setTofData(), and MrpcDBSTrack::two_digis_from_same_strip().

◆ Get_NINO_leadingedge_resolution()

double BesTofDigitizerEcV4_dbs::Get_NINO_leadingedge_resolution ( double induced_charge_fc)

Definition at line 1353 of file BesTofDigitizerEcV4_dbs.cc.

1354{
1355
1356 double result =0;
1357
1358 if(induced_charge_fc > 50.)
1359 {
1360 result = 67.6737*exp(-induced_charge_fc/50.9995-0.27755)+9.06223;
1361 }
1362 else
1363 {
1364 result = 176.322-2.98345*induced_charge_fc;
1365 }
1366
1367 return result;
1368
1369}

Referenced by Digitize().

◆ Get_NINO_pulselength_resolution()

double BesTofDigitizerEcV4_dbs::Get_NINO_pulselength_resolution ( double induced_charge_fc)

Definition at line 1342 of file BesTofDigitizerEcV4_dbs.cc.

1343{
1344
1345 double result = (64.3326*exp(-induced_charge_fc/25.4638 + 0.944184)+19.4846 );
1346
1347 return result;
1348
1349}

Referenced by Digitize().

◆ Get_stripnumber_from_unique_identifier()

G4int BesTofDigitizerEcV4_dbs::Get_stripnumber_from_unique_identifier ( G4int unique_identifier_f)
static

Definition at line 901 of file BesTofDigitizerEcV4_dbs.cc.

902{
903
904 G4int returnvalue=0;
905
906 returnvalue= unique_identifier_f%25;
907
908 return returnvalue;
909}

Referenced by MrpcDBSRec::execute(), MrpcDBSTrack::findTofDataMRPC(), MrpcDBSTrack::get_mrpc_module_id_neighbour(), MrpcDBSTrack::setCalibration(), MrpcDBSTrack::setTofData(), and MrpcDBSTrack::two_digis_from_same_strip().

◆ GetTransitionTime_dbsmatch()

double BesTofDigitizerEcV4_dbs::GetTransitionTime_dbsmatch ( double delta_t,
int my_strip )
static

Definition at line 719 of file BesTofDigitizerEcV4_dbs.cc.

720{
721 double signal_velocity_copper_mrpc_strip= 0.8*0.299792458; // unit: percentage * c_0 * mm/ps
722 double strip_length=0;
723
724 if(my_strip==1 || my_strip==2 ) strip_length=88;
725 else if(my_strip==3 || my_strip==4 ) strip_length=92;
726 else if(my_strip==5 || my_strip==6 ) strip_length=98;
727 else if(my_strip==7 || my_strip==8 ) strip_length=102;
728 else if(my_strip==9 || my_strip==10 ) strip_length=108;
729 else if(my_strip==11 || my_strip==12 ) strip_length=112;
730 else if(my_strip==13 || my_strip==14 ) strip_length=118;
731 else if(my_strip==15 || my_strip==16 ) strip_length=122;
732 else if(my_strip==17 || my_strip==18 ) strip_length=128;
733 else if(my_strip==19|| my_strip==20 ) strip_length=132;
734 else if(my_strip==21|| my_strip==22 ) strip_length=138;
735 else if(my_strip==23 || my_strip==24 ) strip_length=142;
736
737
738
739 double trans_time =(strip_length/signal_velocity_copper_mrpc_strip/1000+delta_t)/2.;
740
741
742 return trans_time;
743
744}

Referenced by MrpcDBSTrack::setCalibration().

◆ GetTransitionTime_extrap_track()

double BesTofDigitizerEcV4_dbs::GetTransitionTime_extrap_track ( G4double x_mm,
G4double y_mm,
int partId_f,
int module_mrpc_f,
int my_strip )
static

Definition at line 752 of file BesTofDigitizerEcV4_dbs.cc.

753{
754
755 G4double signal_velocity_copper_mrpc_strip= 0.8*0.299792458; // unit: percentage * c_0 * mm/ps
756
757 //Within this functions all numbers are in mm or in degree!!!!
758 G4double length_to_border_of_the_readoutstrip=0;
759 G4double result_radius= sqrt(x_mm*x_mm+y_mm*y_mm); //Radius where the particle did hit the detector
760 G4double phi= atan2(y_mm,x_mm)*180./pi; //Phi angle of the particle hitting the detector
761
762 if(-180<phi && phi<90){phi=phi+270;} // Transformation into a coordinate-system, whose coordinates do increase to the right/left side. Left or right does depend on the endcap!
763 else {phi=phi-90;}
764
765 G4double result_phi=0;
766 if(partId_f!=3 && partId_f!=6) result_phi=0.;
767 else result_phi=9.97312; //This is the angle to which the first and second layer are shifted in degree!!!
768
769 if(partId_f==3 || partId_f==4) result_phi = -((19-module_mrpc_f-1)*20. - phi + result_phi);
770 else result_phi = (module_mrpc_f-1)*20. - phi + result_phi; // If > 0 Readoutstripnumber -> 2,4,6,8,10,12,14.16,18,20,22,24
771 //The 18th and the 1st module are located between 350 and 10 degree. This has to be taken into account
772
773 if(partId_f==4 && module_mrpc_f==18 && result_phi>350) result_phi=result_phi-360;
774 if(partId_f==5 && module_mrpc_f==1 && result_phi<-350) result_phi=result_phi+360;
775
776
777 result_radius= result_radius*cos(result_phi*pi/180.0) - 634.56 + 159.5 - 3.82;
778
779 double result_x=((result_radius/cos(result_phi*pi/180.0)+634.56-159.5+3.82)*sin(result_phi*pi/180));
780
781
782
783 if( (my_strip%2)==0 )
784 {
785 if(my_strip==2)
786 length_to_border_of_the_readoutstrip=sqrt( (44-result_x)*(44-result_x) + fabs((result_radius-0.)*(result_radius-0.)));
787 else if(my_strip==4)
788 length_to_border_of_the_readoutstrip=sqrt( (46-result_x)*(46-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*1)));
789 else if(my_strip==6)
790 length_to_border_of_the_readoutstrip=sqrt( (49-result_x)*(49-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*2)));
791 else if(my_strip==8)
792 length_to_border_of_the_readoutstrip=sqrt( (51-result_x)*(51-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*3)));
793 else if(my_strip==10)
794 length_to_border_of_the_readoutstrip=sqrt( (54-result_x)*(54-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*4)));
795 else if(my_strip==12)
796 length_to_border_of_the_readoutstrip=sqrt( (56-result_x)*(56-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*5)));
797 else if(my_strip==14)
798 length_to_border_of_the_readoutstrip=sqrt( (59-result_x)*(59-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*6)));
799 else if(my_strip==16)
800 length_to_border_of_the_readoutstrip=sqrt( (61-result_x)*(61-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*7)));
801 else if(my_strip==18)
802 length_to_border_of_the_readoutstrip=sqrt( (64-result_x)*(64-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*8)));
803 else if(my_strip==20)
804 length_to_border_of_the_readoutstrip=sqrt( (66-result_x)*(66-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*9)));
805 else if(my_strip==22)
806 length_to_border_of_the_readoutstrip=sqrt( (69-result_x)*(69-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*10)));
807 else if(my_strip==24)
808 length_to_border_of_the_readoutstrip=sqrt( (71-result_x)*(71-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*11)));
809 }
810 else if((my_strip%2)!=0)
811 {
812 if(my_strip==1)
813 length_to_border_of_the_readoutstrip=sqrt( (-44-result_x)*(-44-result_x) + fabs((result_radius-0.)*(result_radius-0.)));
814 else if(my_strip==3)
815 length_to_border_of_the_readoutstrip=sqrt( (-46-result_x)*(-46-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*1)));
816 else if(my_strip==5)
817 length_to_border_of_the_readoutstrip=sqrt( (-49-result_x)*(-49-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*2)));
818 else if(my_strip==7)
819 length_to_border_of_the_readoutstrip=sqrt( (-51-result_x)*(-51-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*3)));
820 else if(my_strip==9)
821 length_to_border_of_the_readoutstrip=sqrt( (-54-result_x)*(-54-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*4)));
822 else if(my_strip==11)
823 length_to_border_of_the_readoutstrip=sqrt( (-56-result_x)*(-56-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*5)));
824 else if(my_strip==13)
825 length_to_border_of_the_readoutstrip=sqrt( (-59-result_x)*(-59-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*6)));
826 else if(my_strip==15)
827 length_to_border_of_the_readoutstrip=sqrt( (-61-result_x)*(-61-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*7)));
828 else if(my_strip==17)
829 length_to_border_of_the_readoutstrip=sqrt( (-64-result_x)*(-64-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*8)));
830 else if(my_strip==19)
831 length_to_border_of_the_readoutstrip=sqrt( (-66-result_x)*(-66-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*9)));
832 else if(my_strip==21)
833 length_to_border_of_the_readoutstrip=sqrt( (-69-result_x)*(-69-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*10)));
834 else if(my_strip==23)
835 length_to_border_of_the_readoutstrip=sqrt( (-71-result_x)*(-71-result_x) + fabs((result_radius-(25.0+4.0)*1.)*(result_radius-(25.0+4.0)*11)));
836 }
837
838 double transition_time = length_to_border_of_the_readoutstrip/signal_velocity_copper_mrpc_strip/1000.; //length is in mm , 3*10^8 m/s = 0.3 mm/ps --> Unit of transitionstime is ns (due to division by 1000.)
839 return transition_time;
840
841
842
843}

Referenced by MrpcDBSRec::execute(), and MrpcDBSTrack::setCalibration().

◆ Produce_unique_identifier()

G4int BesTofDigitizerEcV4_dbs::Produce_unique_identifier ( G4int module_mrpc_f,
G4int readoutstripnumber_f )
static

Definition at line 876 of file BesTofDigitizerEcV4_dbs.cc.

877{
878 G4int returnvalue=0;
879
880
881 returnvalue= module_mrpc_f*25 + readoutstripnumber_f;
882
883 return returnvalue;
884}

Referenced by Digitize(), MrpcDBSTrack::setTofData(), and ExtSteppingAction::UserSteppingAction().

◆ simulate_avalanche()

void BesTofDigitizerEcV4_dbs::simulate_avalanche ( G4int actual_strip,
vector< adc_info > * avalanche_info,
double * avalanche_threshold_time,
double * induced_charge )

Attachment coefficient: /m

Definition at line 1033 of file BesTofDigitizerEcV4_dbs.cc.

1034{
1035
1036 //Originally this function was a stand alone program. That's why the style of programming may be a little strange:
1037 //It depends on a method described in the Dissertation of Christian Lippmann and is known in there as 1-D Model.
1038 //Simulation with this model have been done for a chamber of the same type as implemented: Simulation study on the operation of a multi-gap resistive plate chamber - Meas. Sci. Technol. 17 (2006) 123-127
1039 //
1040 // This function use SI values as input. However it returns charge in picoColomb. The time is returned in seconds!
1041 //
1042 //General sequence of the programm:
1043 //1. Initialize all the values: Distribut the avalanche along the gap and set the number of primary electrons.
1044 //2. Invoke stepping action to simulate the avalanche evolution through the gap.
1045 //3. Construct the signal
1046
1047
1048
1049 int number_of_gaps=12;
1050 double alpha=144800.; // Townsend coefficient: /m
1051 double eta=5013; /// Attachment coefficient: /m
1052
1053 double thick_gap = 0.00022; // m //gap thickness
1054 const int steps =200;//Into how many steps is the gap divided?
1055 double user_stepwidth=thick_gap/(double)steps;
1056 double v_drift = 210.9e3; // m/s
1057
1058 double E_weight = 1./(6.*0.22e-3+(5.*0.4e-3+2*0.55e-3)/(8.)+2.*0.07e-3/(3.1)); //epsilon_glass= 8 / epsilon_mylar=3.1 /d_mylar=0.07e-3 , d_g= 0.22e-3 ,
1059 double physics_e= 1.602176462e-19;
1060 double threshold=93622.7;//Threshold in "number of electrons" when the signal will be detected: value_in_fc / physics_e
1061
1062 //double alpha=76750; double eta=8152; double v_drift=163.1e3; //85 kV/cm
1063 //double alpha=82060; double eta=7728; double v_drift=167.3e3; //87.5 kV/cm
1064 //double alpha=87560; double eta=7358; double v_drift=171.1e3;//90 kv/cm
1065 //double alpha=99020; double eta=6725; double v_drift=179.7e3;//95 kv/cm
1066 //double alpha=115200; double eta=6033; double v_drift=191.2e3;//102 kv/cm
1067 //double alpha=124700; double eta=5683; double v_drift=197.8e3;//106 kv/cm
1068 //double alpha=134600; double eta=5374; double v_drift=204.3e3;//110 kv/cm
1069 //double alpha=144500; double eta=5028; double v_drift=210.9e3;//114kv/cm
1070 //double alpha=154700; double eta=4823; double v_drift=217.5e3;//118kv/cm
1071
1072
1073
1074
1075//Intialize
1076 vector<int> primary_clusters_in_gap(number_of_gaps);
1077 for(int i=0;i<number_of_gaps;i++)primary_clusters_in_gap[i]=0;
1078
1079 vector< vector< vector<double> > > step_avalanche(number_of_gaps);
1080 vector< vector< vector<long int> > > elec_avalche(number_of_gaps);
1081
1082 //How many cluster do we have for the different steps:
1083 for(int i=0;i<number_of_gaps;i++)
1084 {
1085 int counter_avalanches=0;
1086 for(int j=0;j<avalanche_info->size();j++)
1087 {
1088 if(actual_strip==(*avalanche_info)[j].readoutstrip && ((i+1)==(*avalanche_info)[j].gap) )
1089 {counter_avalanches++;}
1090 }
1091
1092 primary_clusters_in_gap[i]=counter_avalanches;
1093 step_avalanche[i].resize(primary_clusters_in_gap[i]);
1094 elec_avalche[i].resize(primary_clusters_in_gap[i]);
1095 }
1096
1097
1098 for(int j=0;j<number_of_gaps;j++)
1099 {
1100 for(int i=0;i<primary_clusters_in_gap[j];i++)
1101 {
1102 step_avalanche[j][i].resize(steps+2);
1103 elec_avalche[j][i].resize(steps+2);
1104 }
1105 }
1106
1107 for(int h=0;h<number_of_gaps;h++) {for(int i=0;i<primary_clusters_in_gap[h];i++){for(int j=0;j<steps+2;j++){step_avalanche[h][i][j]=0.;elec_avalche[h][i][j]=0;}}} //Set all arrays equal to zero
1108
1109 vector<int> avalanches_in_gap(number_of_gaps);
1110
1111
1112
1113
1114//Distribute the cluster along the gap and set the number of ions
1115 for(int h=0;h<number_of_gaps;h++)
1116 {
1117
1118
1119 int i=0;
1120 for(int j=0;j<avalanche_info->size();j++)
1121 {
1122 if(actual_strip==(*avalanche_info)[j].readoutstrip && ((h+1)==(*avalanche_info)[j].gap) )
1123 {
1124 step_avalanche[h][i][0]=thick_gap-((*avalanche_info)[j].zdistance/m);
1125 elec_avalche[h][i][0]=(*avalanche_info)[j].ions;
1126 i++;
1127 }
1128 }
1129
1130 avalanches_in_gap[h]=primary_clusters_in_gap[h];
1131
1132 }//close for h (number_of_gaps)
1133
1134 /*
1135 for(int h=0;h<number_of_gaps;h++)
1136 {
1137 for(int j=0;j<primary_clusters_in_gap[h];j++)
1138 {
1139 G4cout << "Gap " << h+1 << " Avalanche " << j+1 << " Pos in gap " << step_avalanche[h][j][0] << G4endl;
1140 }
1141 }
1142 */
1143
1144
1145 for(int h=0;h<number_of_gaps;h++)
1146 {
1147 // G4cout << "Avalacnhes in gap " << h+1 << " are " << avalanches_in_gap[h] << G4endl;
1148 long int tot_elec_in_gap=0;
1149 if(avalanches_in_gap[h]!=0)
1150 {
1151
1152 int i=0;
1153
1154
1155 vector<bool> is_saturated(avalanches_in_gap[h]);
1156 for(int e=0;e<avalanches_in_gap[h];e++) is_saturated[e]=false;
1157
1158
1159
1160 do //Steppingaction: Moves the avalanche through the gap!
1161 {
1162
1163
1164 tot_elec_in_gap=0;
1165 for(int j=0;j<primary_clusters_in_gap[h];j++)
1166 {
1167
1168
1169
1170 elec_avalche[h][j][i+1]=adcsignal_simulate_step(elec_avalche[h][j][i],user_stepwidth,alpha,eta,is_saturated[j]);
1171 step_avalanche[h][j][i+1]= step_avalanche[h][j][i] + user_stepwidth;
1172 tot_elec_in_gap=tot_elec_in_gap + elec_avalche[h][j][i+1];
1173 //if((h+1)==12) G4cout<< "Gap "<<h+1<<" cluster "<<j+1<<" = "<< elec_avalche[h][j][i+1]<<" Pos "<<step_avalanche[h][j][i+1]<< " tot " <<tot_elec_in_gap<<G4endl;
1174
1175 if(elec_avalche[h][j][i+1]>1.5e7){is_saturated[j]=true;}else{is_saturated[j]=false;}
1176 }
1177
1178
1179 i++;
1180
1181
1182 for(int j=0;j<primary_clusters_in_gap[h];j++)
1183 {
1184
1185 if((step_avalanche[h][j][i])>=thick_gap)
1186 {
1187 primary_clusters_in_gap[h]=(primary_clusters_in_gap[h]-1);
1188 }
1189 }
1190
1191
1192 }
1193 while(step_avalanche[h][0][i]<thick_gap); //If the first avalanche leaves the gap, quit the stepping action
1194
1195
1196 }//Close if(avalanches_in_gap!=0)
1197 }//close for over gaps
1198
1199
1200
1201
1202//Produce signal
1203 double signal=0;//stores the total sum of the induced signal!
1204 bool over_threshold = false;
1205
1206 for(int i=0;i<steps;i++)
1207 {
1208 for(int h=0;h<number_of_gaps;h++)
1209 {
1210 for(int j=0;j<avalanches_in_gap[h];j++)
1211 {
1212 signal= signal + elec_avalche[h][j][i];
1213 }
1214 }
1215
1216 if(over_threshold == false && signal>threshold) {over_threshold=true; *avalanche_threshold_time=(thick_gap/steps/v_drift*(i+1));}
1217
1218 }//close the for loop over steps.
1219
1220
1221
1222
1223 //Return the value in pC:
1224
1225 //G4cout <<" *avalanche_threshold_time " << *avalanche_threshold_time << G4endl;
1226 //G4cout << "Induced Charge "<< signal*E_weight*v_drift*physics_e*1e12 << G4endl;
1227 *induced_charge=signal*E_weight*v_drift*physics_e*1e12 * thick_gap/steps/v_drift;
1228
1229
1230
1231}//close simulate_avalanche
long int adcsignal_simulate_step(int n_elec, double stepwidth, double alpha, double eta, bool saturated)

Referenced by Digitize().

◆ Smear_gaussian()

G4double BesTofDigitizerEcV4_dbs::Smear_gaussian ( G4double variable,
G4double gauss_mean,
G4double gauss_sigma )

Definition at line 357 of file BesTofDigitizerEcV4_dbs.cc.

358{
359
360 G4double gaussrand = G4RandGauss::shoot(gauss_mean,gauss_sigma);//returns x-value
361 return (variable+gaussrand);
362
363}

Referenced by Digitize().


The documentation for this class was generated from the following files: