BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
BesEmcDigitizer Class Reference

#include <BesEmcDigitizer.hh>

Inheritance diagram for BesEmcDigitizer:

Public Member Functions

 BesEmcDigitizer (G4String modName)
 ~BesEmcDigitizer ()
virtual void Digitize ()
virtual void GroupHits (BesEmcHitsCollection *)
virtual void AddNoise5x5 (G4double coherentNoise)
virtual void AddNoiseAll (G4double coherentNoise)
 BesEmcDigitizer (G4String modName)
 ~BesEmcDigitizer ()
virtual void Digitize ()
virtual void GroupHits (BesEmcHitsCollection *)
virtual void AddNoise5x5 (G4double coherentNoise)
virtual void AddNoiseAll (G4double coherentNoise)
 BesEmcDigitizer (G4String modName)
 ~BesEmcDigitizer ()
virtual void Digitize ()
virtual void GroupHits (BesEmcHitsCollection *)
virtual void AddNoise5x5 (G4double coherentNoise)
virtual void AddNoiseAll (G4double coherentNoise)

Detailed Description

Constructor & Destructor Documentation

◆ BesEmcDigitizer() [1/3]

BesEmcDigitizer::BesEmcDigitizer ( G4String modName)

Definition at line 29 of file BesEmcDigitizer.cc.

30 : G4VDigitizerModule( modName ), m_emcCalibConstSvc( 0 ) {
31 collectionName.push_back( "BesEmcDigitsCollection" );
32 m_besEmcDigitsCollection = 0;
33
34 // retrieve G4Svc
35 ISvcLocator* svcLocator = Gaudi::svcLocator();
36 IG4Svc* iG4Svc;
37 StatusCode sc = svcLocator->service( "G4Svc", iG4Svc );
38 // m_G4Svc=dynamic_cast<G4Svc *>(iG4Svc);
39 m_G4Svc = iG4Svc;
40
41 // get Emc Ntuple from G4Svc
42 if ( m_G4Svc->EmcRootFlag() )
43 {
44 m_tupleEmc1 = m_G4Svc->GetTupleEmc1();
45 sc = m_tupleEmc1->addItem( "partId", m_partId );
46 sc = m_tupleEmc1->addItem( "nTheta", m_nTheta );
47 sc = m_tupleEmc1->addItem( "nPhi", m_nPhi );
48 sc = m_tupleEmc1->addItem( "edep", m_eDep );
49 sc = m_tupleEmc1->addItem( "nHits", m_nHits );
50 sc = m_tupleEmc1->addItem( "adc", m_adc );
51 sc = m_tupleEmc1->addItem( "tdc", m_tdc );
52
53 m_tupleEmc2 = m_G4Svc->GetTupleEmc2();
54 sc = m_tupleEmc2->addItem( "etot", m_eTot );
55 sc = m_tupleEmc2->addItem( "nDigi", m_nDigi );
56 }
57
58 // Get EmcCalibConstSvc.
59 sc = svcLocator->service( "EmcCalibConstSvc", m_emcCalibConstSvc );
60 if ( sc != StatusCode::SUCCESS )
61 { G4cout << "BesEmcDigitizer Error: Can't get EmcCalibConstSvc." << G4endl; }
62}

◆ ~BesEmcDigitizer() [1/3]

BesEmcDigitizer::~BesEmcDigitizer ( )

Definition at line 64 of file BesEmcDigitizer.cc.

64{}

◆ BesEmcDigitizer() [2/3]

BesEmcDigitizer::BesEmcDigitizer ( G4String modName)

◆ ~BesEmcDigitizer() [2/3]

BesEmcDigitizer::~BesEmcDigitizer ( )

◆ BesEmcDigitizer() [3/3]

BesEmcDigitizer::BesEmcDigitizer ( G4String modName)

◆ ~BesEmcDigitizer() [3/3]

BesEmcDigitizer::~BesEmcDigitizer ( )

Member Function Documentation

◆ AddNoise5x5() [1/3]

void BesEmcDigitizer::AddNoise5x5 ( G4double coherentNoise)
virtual

Definition at line 277 of file BesEmcDigitizer.cc.

277 {
278 BesEmcParameter& emcPara = BesEmcParameter::GetInstance();
279 vector<BesEmcDigi*>* vecDC = m_besEmcDigitsCollection->GetVector();
280 G4int nDigi = m_besEmcDigitsCollection->entries();
281 G4int partMax, thetaMax, phiMax;
282 partMax = thetaMax = phiMax = -99;
283 G4double eMax = 0;
284
285 for ( G4int i = 0; i < nDigi; i++ )
286 {
287 BesEmcDigi* digi = ( *vecDC )[i];
288 double eDigi = digi->GetEnergy();
289 if ( eDigi > eMax )
290 {
291 eMax = eDigi;
292 partMax = digi->GetPartId();
293 thetaMax = digi->GetThetaNb();
294 phiMax = digi->GetPhiNb();
295 }
296 }
297
298 if ( partMax == 1 )
299 { // barrel
300 G4int thetan, thetap, phin, phip;
301 thetan = thetaMax - 2;
302 thetap = thetaMax + 2;
303 phin = phiMax - 2;
304 phip = phiMax + 2;
305
306 if ( thetaMax == 0 )
307 { // #0
308 thetan = thetaMax;
309 }
310 else if ( thetaMax == 1 )
311 { // #1
312 thetan = thetaMax - 1;
313 }
314 else if ( thetaMax == emcPara.GetBSCNbTheta() * 2 - 1 )
315 { // #43
316 thetap = thetaMax;
317 }
318 else if ( thetaMax == emcPara.GetBSCNbTheta() * 2 - 2 )
319 { // #42
320 thetap = thetaMax + 1;
321 }
322
323 if ( phiMax == 0 ) { phin = emcPara.GetBSCNbPhi() - 2; }
324 else if ( phiMax == 1 ) { phin = emcPara.GetBSCNbPhi() - 2; }
325 else if ( phiMax == emcPara.GetBSCNbPhi() - 1 )
326 { // #119
327 phip = 1;
328 }
329 else if ( phiMax == emcPara.GetBSCNbPhi() - 2 )
330 { // #118
331 phip = 0;
332 }
333
334 for ( G4int theta = thetan; theta <= thetap; theta++ )
335 {
336 for ( G4int phi = phin; phi <= phip; phi++ )
337 {
338 G4bool flag = true;
339
340 if ( nDigi > 0 )
341 {
342 for ( G4int i = 0; i < nDigi; i++ )
343 {
344 BesEmcDigi* digi = ( *vecDC )[i];
345 if ( partMax == digi->GetPartId() && theta == digi->GetThetaNb() &&
346 phi == digi->GetPhiNb() )
347 {
348 flag = false;
349 break;
350 }
351 }
352 }
353
354 if ( flag )
355 {
356 BesEmcDigi* digi = new BesEmcDigi;
357 BesEmcWaveform* wave = digi->GetWaveform();
358 digi->SetPartId( partMax );
359 digi->SetThetaNb( theta );
360 digi->SetPhiNb( phi );
361 digi->SetEnergy( 0 );
362 digi->SetTime( m_G4Svc->EmcTime() );
363 digi->SetTrackIndex( -9 );
364
365 wave->updateWaveform( digi );
366 wave->addElecNoise( m_G4Svc->EmcIncoherentNoise(), coherentNoise );
367 wave->digitize();
368
369 G4long bin = 0; // time
370 m_energy = wave->max( bin );
371
372 if ( m_G4Svc->EmcLightOutput() )
373 { m_energy *= emcPara.GetLightOutput( partMax, theta, phi ); }
374
375 digi->SetEnergy( m_energy );
376 digi->SetTime( bin );
377 digi->SetWaveform( wave );
378
379 // Correction of electronics bias
380 G4double ecorr;
381 if ( m_energy < 625. )
382 {
383 ecorr = -0.1285 * log( m_energy / 6805. ); // noise=0.5
384 }
385 else { ecorr = -2.418 + 9.513e-4 * m_energy; }
386
387 if ( m_energy - ecorr > m_G4Svc->EmcNoiseThreshold() )
388 { m_besEmcDigitsCollection->insert( digi ); }
389 else { delete digi; }
390 }
391 } // phi
392 } // theta
393
394 } // part==1
395}
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per bin
Definition FoamA.h:85
static BesEmcParameter & GetInstance()
void updateWaveform(BesEmcHit *)
void addElecNoise(G4double, G4double)
G4double max(G4long &binOfMax) const

Referenced by Digitize().

◆ AddNoise5x5() [2/3]

virtual void BesEmcDigitizer::AddNoise5x5 ( G4double coherentNoise)
virtual

◆ AddNoise5x5() [3/3]

virtual void BesEmcDigitizer::AddNoise5x5 ( G4double coherentNoise)
virtual

◆ AddNoiseAll() [1/3]

void BesEmcDigitizer::AddNoiseAll ( G4double coherentNoise)
virtual

Definition at line 397 of file BesEmcDigitizer.cc.

397 {
398 BesEmcParameter& emcPara = BesEmcParameter::GetInstance();
399 vector<BesEmcDigi*>* vecDC = m_besEmcDigitsCollection->GetVector();
400 G4int nDigi = m_besEmcDigitsCollection->entries();
401 // G4cout<<"nDigi="<<nDigi<<G4endl;
402
403 for ( G4int part = 0; part < 3; part++ )
404 {
405
406 G4int thetaNb;
407 if ( part == 1 )
408 { // barrel
409 thetaNb = emcPara.GetBSCNbTheta() * 2;
410 }
411 else
412 { // endcap
413 thetaNb = 6;
414 }
415
416 for ( G4int theta = 0; theta < thetaNb; theta++ )
417 {
418
419 G4int phiNb;
420 if ( part == 1 ) { phiNb = emcPara.GetBSCNbPhi(); }
421 else { phiNb = emcPara.GetCryInOneLayer( theta ); }
422
423 for ( G4int phi = 0; phi < phiNb; phi++ )
424 {
425
426 G4bool flag = true;
427
428 if ( nDigi > 0 )
429 {
430 // G4cout<<"nDigi="<<nDigi<<"\t";
431
432 for ( G4int i = 0; i < nDigi; i++ )
433 {
434 BesEmcDigi* digi = ( *vecDC )[i];
435 if ( part == digi->GetPartId() && theta == digi->GetThetaNb() &&
436 phi == digi->GetPhiNb() )
437 {
438 // cout<<theta<<"\t"<<phi<<endl;
439 flag = false;
440 break;
441 }
442 }
443 }
444
445 if ( flag )
446 {
447 BesEmcDigi* digi = new BesEmcDigi;
448 digi->SetTrackIndex( -9 );
449 digi->SetPartId( part );
450 digi->SetThetaNb( theta );
451 digi->SetPhiNb( phi );
452
453 bool fastSimulation = true;
454 if ( fastSimulation )
455 {
456
457 m_energy = RandGauss::shoot() * m_G4Svc->EmcNoiseSigma();
458 m_energy += m_G4Svc->EmcNoiseMean();
459 digi->SetTime( (G4int)( G4UniformRand() * 60 ) );
460 }
461 else
462 {
463
464 BesEmcWaveform* wave = digi->GetWaveform();
465 digi->SetTime( m_G4Svc->EmcTime() );
466
467 wave->updateWaveform( digi );
468 wave->addElecNoise( m_G4Svc->EmcIncoherentNoise(), coherentNoise );
469 wave->digitize();
470
471 G4long bin = 0; // time
472 m_energy = wave->max( bin );
473 digi->SetTime( bin );
474 digi->SetWaveform( wave );
475 }
476
477 if ( m_G4Svc->EmcLightOutput() )
478 { m_energy *= emcPara.GetLightOutput( part, theta, phi ); }
479 digi->SetEnergy( m_energy );
480
481 // Correction of electronics bias
482 G4double ecorr;
483 if ( m_energy < 625. )
484 {
485 ecorr = -0.1285 * log( m_energy / 6805. ); // noise=0.5
486 }
487 else { ecorr = -2.418 + 9.513e-4 * m_energy; }
488
489 if ( m_energy - ecorr > m_G4Svc->EmcNoiseThreshold() )
490 { m_besEmcDigitsCollection->insert( digi ); }
491 else { delete digi; }
492 }
493
494 } // phi
495 } // theta
496 } // part
497}

◆ AddNoiseAll() [2/3]

virtual void BesEmcDigitizer::AddNoiseAll ( G4double coherentNoise)
virtual

◆ AddNoiseAll() [3/3]

virtual void BesEmcDigitizer::AddNoiseAll ( G4double coherentNoise)
virtual

◆ Digitize() [1/3]

void BesEmcDigitizer::Digitize ( )
virtual

Definition at line 68 of file BesEmcDigitizer.cc.

68 {
69 Initialize();
70
71 m_besEmcDigitsCollection =
72 new BesEmcDigitsCollection( "BesEmcDigitizer", "BesEmcDigitsCollection" );
73 G4DigiManager* DigiMan = G4DigiManager::GetDMpointer();
74
75 // hits collection ID
76 G4int EHCID;
77 EHCID = DigiMan->GetHitsCollectionID( "BesEmcHitsCollection" );
78
79 // hits collection
80 BesEmcHitsCollection* EHC = 0;
81 EHC = (BesEmcHitsCollection*)( DigiMan->GetHitsCollection( EHCID ) );
82
83 if ( EHC )
84 {
85 // BesEmcParameter& emcPara=BesEmcParameter::GetInstance();
86 m_crystalGroup = new vector<CrystalSingle*>;
87 GroupHits( EHC );
88 G4int size = m_crystalGroup->size();
89 CrystalSingle* cryst;
90 G4int partId, nTheta, nPhi, nHits;
91 G4double eTot = 0, eDigi;
92 BesEmcHit* hit;
93
94 G4double coherentNoise = RandGauss::shoot() * m_G4Svc->EmcCoherentNoise();
95
96 for ( G4int i = 0; i < size; i++ )
97 {
98 cryst = ( *m_crystalGroup )[i]; // all hits in a crystal
99 partId = cryst->GetPartId();
100 nTheta = cryst->GetNTheta();
101 nPhi = cryst->GetNPhi();
102 nHits = cryst->GetHitIndexes()->size();
103 eDigi = cryst->GetEdep();
104 eTot += eDigi;
105
106 BesEmcDigi* digi = new BesEmcDigi;
107 BesEmcWaveform* wave = digi->GetWaveform();
108 G4long bin = 0; // time
109
110 const int indexSize = 200;
111 G4double e[indexSize]; // energy of the same index
112 for ( G4int i = 0; i < indexSize; i++ ) e[i] = 0;
113 G4int index = 0;
114 G4double energy = 0;
115
116 for ( G4int j = 0; j < nHits; j++ )
117 {
118 hit = ( *EHC )[( *( cryst->GetHitIndexes() ) )[j]];
119 energy = hit->GetEdepCrystal();
120 index = hit->GetTrackIndex();
121 if ( index < indexSize && index >= 0 ) e[index] += energy;
122 else G4cout << "Track index overload!" << G4endl;
123 }
124
125 G4double maxi = e[0]; // find the index which gives the most energy in one crystal
126 for ( G4int i = 1; i < indexSize; i++ )
127 {
128 if ( e[i] > maxi )
129 {
130 maxi = e[i];
131 index = i;
132 }
133 }
134
135 if ( eDigi > 0 )
136 {
137 digi->SetPartId( partId );
138 digi->SetThetaNb( nTheta );
139 digi->SetPhiNb( nPhi );
140 digi->SetEnergy( eDigi );
141 digi->SetTime( m_G4Svc->EmcTime() );
142 digi->SetTrackIndex( index );
143
144 wave->updateWaveform( digi );
145 if ( m_G4Svc->EmcNoiseLevel() > 0 )
146 wave->addElecNoise( m_G4Svc->EmcIncoherentNoise(), coherentNoise );
147
148 // to avoid error caused by precision, get energy before digitization
149 m_energy = wave->max( bin );
150 // temp code, subtract pedstal
151 m_energy -= 0.46 * MeV;
152 wave->digitize();
153 wave->max( bin );
154
155 if ( m_G4Svc->EmcLightOutput() )
156 {
157 G4int index = m_emcCalibConstSvc->getIndex( partId, nTheta, nPhi );
158 G4double adc2e = m_emcCalibConstSvc->getDigiCalibConst( index );
159
160 G4double CrystalDeadEcut = m_emcCalibConstSvc->getCrystalDeadEcut( index );
161
162 if ( m_G4Svc->EmcElecSaturation() == 1 )
163 {
164 G4double emaxData = m_emcCalibConstSvc->getCrystalEmaxData( index );
165
166 if ( emaxData > 0 ) { adc2e = emaxData / 2.5; }
167 }
168
169 if ( adc2e <= 1e-5 ) // dead channel
170 { m_energy = 0; }
171 else if ( m_G4Svc->EmcElecSatuDead() == 1 && CrystalDeadEcut > 0 &&
172 m_energy / 1000.0 > CrystalDeadEcut ) //"special" dead channel because of
173 // electronics saturation (GeV)
174 {
175 // cout<<"ixtal="<<index<<" CrystalDeadEcut="<<CrystalDeadEcut<<"
176 // m_energy="<<m_energy<<endl;
177 m_energy = 0;
178 }
179 else
180 {
181
182 m_energy /= adc2e;
183
184 // m_energy /= emcPara.GetLightOutput(partId,nTheta,nPhi);
185 }
186 }
187
188 // fill Emc Ntuple1
189 if ( m_G4Svc->EmcRootFlag() )
190 {
191 m_partId = partId;
192 m_nTheta = nTheta;
193 m_nPhi = nPhi;
194 m_eDep = eDigi;
195 m_nHits = nHits;
196 m_adc = m_energy;
197 m_tdc = bin;
198 m_tupleEmc1->write();
199 }
200
201 digi->SetEnergy( m_energy );
202 digi->SetTime( bin );
203 digi->SetWaveform( wave );
204 m_besEmcDigitsCollection->insert( digi );
205 }
206 }
207
208 // add to noise to no signal crystals
209 if ( m_G4Svc->EmcNoiseLevel() == 2 ) AddNoise5x5( coherentNoise );
210 else if ( m_G4Svc->EmcNoiseLevel() == 3 )
211 ;
212 // AddNoiseAll(coherentNoise);
213
214 // fill Emc Ntuple2
215 if ( m_G4Svc->EmcRootFlag() )
216 {
217 m_eTot = eTot;
218 m_nDigi = size;
219 m_tupleEmc2->write();
220 }
221
222 StoreDigiCollection( m_besEmcDigitsCollection );
223
224 for ( size_t i = 0; i < m_crystalGroup->size(); i++ ) { delete ( *m_crystalGroup )[i]; }
225 m_crystalGroup->clear();
226 delete m_crystalGroup;
227 }
228}
G4TDigiCollection< BesEmcDigi > BesEmcDigitsCollection
G4THitsCollection< BesEmcHit > BesEmcHitsCollection
************Class m_ypar INTEGER m_KeyWgt INTEGER m_KeyIHVP INTEGER m_KeyGPS INTEGER m_IsBeamPolarized INTEGER m_EvtGenInterface DOUBLE PRECISION m_Emin DOUBLE PRECISION m_sphot DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_q2 DOUBLE PRECISION m_PolBeam2 DOUBLE PRECISION m_xErrPb *COMMON c_KK2f $ !CMS energy average $ !Spin Polarization vector first beam $ !Spin Polarization vector second beam $ !Beam energy spread[GeV] $ !minimum hadronization energy[GeV] $ !input READ never touch them !$ !debug facility $ !maximum weight $ !inverse alfaQED $ !minimum real photon energy
Definition KK2f.h:50
virtual void AddNoise5x5(G4double coherentNoise)
virtual void GroupHits(BesEmcHitsCollection *)

◆ Digitize() [2/3]

virtual void BesEmcDigitizer::Digitize ( )
virtual

◆ Digitize() [3/3]

virtual void BesEmcDigitizer::Digitize ( )
virtual

◆ GroupHits() [1/3]

void BesEmcDigitizer::GroupHits ( BesEmcHitsCollection * m_EHC)
virtual

Definition at line 230 of file BesEmcDigitizer.cc.

230 {
231 G4int partId, nTheta, nPhi, size, flag;
232 G4double edep;
233 BesEmcHit* hit;
234 G4int nHits = m_EHC->entries();
235
236 // group the hits which are in the same crystal
237 for ( G4int i = 0; i < nHits; i++ )
238 {
239 hit = ( *m_EHC )[i];
240 partId = hit->GetPartId();
241 nTheta = hit->GetNumThetaCrystal();
242 nPhi = hit->GetNumPhiCrystal();
243 edep = hit->GetEdepCrystal();
244 size = m_crystalGroup->size();
245 flag = 0;
246
247 if ( size > 0 )
248 {
249 CrystalSingle* oldCryst;
250 for ( G4int j = 0; j < size; j++ )
251 {
252 oldCryst = ( *m_crystalGroup )[j];
253 if ( ( oldCryst->GetNTheta() == nTheta ) && ( oldCryst->GetNPhi() == nPhi ) &&
254 ( oldCryst->GetPartId() == partId ) )
255 {
256 oldCryst->GetHitIndexes()->push_back( i );
257 oldCryst->AddEdep( edep );
258 flag = 1;
259 break;
260 }
261 }
262 }
263
264 if ( flag == 0 )
265 {
266 CrystalSingle* newCryst = new CrystalSingle;
267 newCryst->SetPartId( partId );
268 newCryst->SetNTheta( nTheta );
269 newCryst->SetNPhi( nPhi );
270 newCryst->SetEdep( edep );
271 newCryst->GetHitIndexes()->push_back( i );
272 m_crystalGroup->push_back( newCryst );
273 }
274 }
275}

Referenced by Digitize().

◆ GroupHits() [2/3]

virtual void BesEmcDigitizer::GroupHits ( BesEmcHitsCollection * )
virtual

◆ GroupHits() [3/3]

virtual void BesEmcDigitizer::GroupHits ( BesEmcHitsCollection * )
virtual

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