| Geant4 Cross Reference |
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // G4IonTable class implementation
27 //
28 // Author: H.Kurashige, 27 June 1998
29 // --------------------------------------------------------------------
30
31 #include "G4IonTable.hh"
32
33 #include "G4AutoDelete.hh"
34 #include "G4HyperNucleiProperties.hh"
35 #include "G4Ions.hh"
36 #include "G4IsotopeProperty.hh"
37 #include "G4MuonicAtom.hh"
38 #include "G4MuonicAtomHelper.hh"
39 #include "G4NucleiProperties.hh"
40 #include "G4NuclideTable.hh"
41 #include "G4ParticleTable.hh"
42 #include "G4PhysicalConstants.hh"
43 #include "G4StateManager.hh"
44 #include "G4SystemOfUnits.hh"
45 #include "G4Threading.hh"
46 #include "G4UImanager.hh"
47 #include "G4VIsotopeTable.hh"
48 #include "G4ios.hh"
49
50 #include <algorithm>
51 #include <iomanip>
52 #include <iostream>
53 #include <sstream>
54 #include <vector>
55
56 // It is very important for multithreaded Geant4 to keep only one copy of the
57 // particle table pointer and the ion table pointer. However, we try to let
58 // each worker thread hold its own copy of the particle dictionary and the
59 // ion list. This implementation is equivalent to make the ion table thread
60 // private. The two shadow ponters are used by each worker thread to copy the
61 // content from the master thread.
62 //
63 G4ThreadLocal G4IonTable::G4IonList* G4IonTable::fIonList = nullptr;
64 G4ThreadLocal std::vector<G4VIsotopeTable*>* G4IonTable::fIsotopeTableList = nullptr;
65 G4IonTable::G4IonList* G4IonTable::fIonListShadow = nullptr;
66 std::vector<G4VIsotopeTable*>* G4IonTable::fIsotopeTableListShadow = nullptr;
67
68 namespace lightions
69 {
70 static const G4ParticleDefinition* p_proton = nullptr;
71 static const G4ParticleDefinition* p_deuteron = nullptr;
72 static const G4ParticleDefinition* p_triton = nullptr;
73 static const G4ParticleDefinition* p_alpha = nullptr;
74 static const G4ParticleDefinition* p_He3 = nullptr;
75 void Init()
76 {
77 if (p_proton != nullptr) return;
78 p_proton = G4ParticleTable::GetParticleTable()->FindParticle("proton");
79 p_deuteron = G4ParticleTable::GetParticleTable()->FindParticle("deuteron");
80 p_triton = G4ParticleTable::GetParticleTable()->FindParticle("triton");
81 p_alpha = G4ParticleTable::GetParticleTable()->FindParticle("alpha");
82 p_He3 = G4ParticleTable::GetParticleTable()->FindParticle("He3");
83 }
84 } // namespace lightions
85
86 namespace antilightions
87 {
88 static const G4ParticleDefinition* p_proton = nullptr;
89 static const G4ParticleDefinition* p_deuteron = nullptr;
90 static const G4ParticleDefinition* p_triton = nullptr;
91 static const G4ParticleDefinition* p_alpha = nullptr;
92 static const G4ParticleDefinition* p_He3 = nullptr;
93 void Init()
94 {
95 if (p_proton != nullptr) return;
96 p_proton = G4ParticleTable::GetParticleTable()->FindParticle("anti_proton");
97 p_deuteron = G4ParticleTable::GetParticleTable()->FindParticle("anti_deuteron");
98 p_triton = G4ParticleTable::GetParticleTable()->FindParticle("anti_triton");
99 p_alpha = G4ParticleTable::GetParticleTable()->FindParticle("anti_alpha");
100 p_He3 = G4ParticleTable::GetParticleTable()->FindParticle("anti_He3");
101 }
102 } // namespace antilightions
103
104 #ifdef G4MULTITHREADED
105 G4Mutex G4IonTable::ionTableMutex = G4MUTEX_INITIALIZER;
106 #endif
107
108 // --------------------------------------------------------------------
109
110 G4IonTable::G4IonTable()
111 {
112 fIonList = new G4IonList();
113
114 // Set up the shadow pointer used by worker threads.
115 //
116 if (fIonListShadow == nullptr) {
117 fIonListShadow = fIonList;
118 }
119
120 fIsotopeTableList = new std::vector<G4VIsotopeTable*>;
121
122 // Set up the shadow pointer used by worker threads.
123 //
124 if (fIsotopeTableListShadow == nullptr) {
125 fIsotopeTableListShadow = fIsotopeTableList;
126 }
127
128 PrepareNuclideTable();
129 RegisterIsotopeTable(pNuclideTable);
130 }
131
132 G4IonTable::~G4IonTable()
133 {
134 // delete IsotopeTable if existing
135 if (fIsotopeTableList != nullptr) {
136 for (const auto fIsotopeTable : *fIsotopeTableList) {
137 if (fIsotopeTable != G4NuclideTable::GetNuclideTable()) {
138 delete fIsotopeTable;
139 }
140 }
141 fIsotopeTableList->clear();
142 delete fIsotopeTableList;
143 }
144 fIsotopeTableList = nullptr;
145
146 if (fIonList == nullptr) return;
147
148 // remove all contents in the Ion List
149 // No need to delete here because all particles are dynamic objects
150 fIonList->clear();
151 delete fIonList;
152 fIonList = nullptr;
153 }
154
155 G4IonTable* G4IonTable::GetIonTable()
156 {
157 return G4ParticleTable::GetParticleTable()->GetIonTable();
158 }
159
160 // Used by each worker thread to copy the content from the master thread
161 void G4IonTable::WorkerG4IonTable()
162 {
163 if (fIonList == nullptr) {
164 fIonList = new G4IonList();
165 }
166 else {
167 fIonList->clear();
168 }
169
170 for (const auto& it : *fIonListShadow) {
171 fIonList->insert(it);
172 }
173
174 // Do not copy Isotope Table to Worker thread
175 //
176 if (fIsotopeTableList == nullptr) {
177 fIsotopeTableList = new std::vector<G4VIsotopeTable*>;
178 for (const auto i : *fIsotopeTableListShadow) {
179 fIsotopeTableList->push_back(i);
180 }
181 }
182 }
183
184 void G4IonTable::InitializeLightIons()
185 {
186 lightions::Init();
187 antilightions::Init();
188 }
189
190 void G4IonTable::DestroyWorkerG4IonTable()
191 {
192 // delete IsotopeTable if existing
193 if (fIsotopeTableList != nullptr) {
194 for (auto fIsotopeTable : *fIsotopeTableList) {
195 if (fIsotopeTable != G4NuclideTable::GetNuclideTable()) {
196 delete fIsotopeTable;
197 }
198 }
199 fIsotopeTableList->clear();
200 delete fIsotopeTableList;
201 }
202 fIsotopeTableList = nullptr;
203
204 if (fIonList == nullptr) return;
205
206 // remove all contents in the Ion List
207 // No need to delete here because all particles are dynamic objects
208 fIonList->clear();
209 delete fIonList;
210 fIonList = nullptr;
211 }
212
213 G4ParticleDefinition* G4IonTable::CreateIon(G4int Z, G4int A, G4double E,
214 G4Ions::G4FloatLevelBase flb)
215 {
216 G4ParticleDefinition* ion = nullptr;
217
218 // check whether GenericIon has processes
219 G4ParticleDefinition* genericIon = G4ParticleTable::GetParticleTable()->GetGenericIon();
220 G4ProcessManager* pman = nullptr;
221 if (genericIon != nullptr) {
222 pman = genericIon->GetProcessManager();
223 }
224 if ((genericIon == nullptr) || (genericIon->GetParticleDefinitionID() < 0) || (pman == nullptr)) {
225 #ifdef G4VERBOSE
226 if (GetVerboseLevel() > 1) {
227 G4cout << "G4IonTable::CreateIon() : can not create ion of "
228 << " Z =" << Z << " A = " << A << " because GenericIon is not ready !!" << G4endl;
229 }
230 #endif
231 G4Exception("G4IonTable::CreateIon()", "PART105", JustWarning,
232 "Can not create ions because GenericIon is not ready");
233 return nullptr;
234 }
235
236 G4double life = 0.0;
237 G4DecayTable* decayTable = nullptr;
238 G4bool stable = true;
239 G4double mu = 0.0;
240 G4double Eex = 0.0;
241 G4int lvl = 0;
242 G4int J = 0;
243
244 const G4IsotopeProperty* fProperty = FindIsotope(Z, A, E, flb);
245 if (fProperty != nullptr) {
246 Eex = fProperty->GetEnergy();
247 lvl = fProperty->GetIsomerLevel();
248 J = fProperty->GetiSpin();
249 life = fProperty->GetLifeTime();
250 mu = fProperty->GetMagneticMoment();
251 decayTable = fProperty->GetDecayTable();
252 stable = (life <= 0.) || (decayTable == nullptr);
253 lvl = fProperty->GetIsomerLevel();
254 if (lvl < 0) lvl = 9;
255 }
256 else {
257 #ifdef G4VERBOSE
258 if (GetVerboseLevel() > 1) {
259 G4ExceptionDescription ed;
260 ed << "G4IonTable::CreateIon(): G4IsotopeProperty object is not found for"
261 << " Z = " << Z << " A = " << A << " E = " << E / keV << " (keV)";
262 if (flb != G4Ions::G4FloatLevelBase::no_Float) {
263 ed << " FloatingLevel +" << G4Ions::FloatLevelBaseChar(flb);
264 }
265 ed << ".\n"
266 << " Physics quantities such as life are not set for this ion.";
267 G4Exception("G4IonTable::CreateIon()", "PART70105", JustWarning, ed);
268 }
269 #endif
270 // excitation energy
271 Eex = E;
272 // lvl is assigned to 9 temporarily
273 if (Eex > 0.0) lvl = 9;
274 }
275
276 // Eex = G4NuclideTable::Round(Eex);
277 if (Eex == 0.0) lvl = 0;
278 // ion name
279 G4String name = "";
280 /////////////if (lvl<9) name = GetIonName(Z, A, lvl);
281 if (lvl == 0 && flb == G4Ions::G4FloatLevelBase::no_Float)
282 name = GetIonName(Z, A, lvl);
283 else
284 name = GetIonName(Z, A, Eex, flb);
285
286 // PDG encoding
287 G4int encoding = GetNucleusEncoding(Z, A, E, lvl);
288
289 // PDG mass
290 G4double mass = GetNucleusMass(Z, A) + Eex;
291
292 // PDG charge is set to one of nucleus
293 G4double charge = G4double(Z) * eplus;
294
295 // create an ion
296 // spin, parity, isospin values are fixed
297
298 // Request lock for particle table accesses. Some changes are inside
299 // this critical region.
300 //
301 // clang-format off
302 ion = new G4Ions( name, mass, 0.0*MeV, charge,
303 J, +1, 0,
304 0, 0, 0,
305 "nucleus", 0, A, encoding,
306 stable, life, decayTable, false,
307 "generic", 0,
308 Eex, lvl );
309 // clang-format on
310
311 // Release lock for particle table accesses.
312 //
313 ion->SetPDGMagneticMoment(mu);
314 static_cast<G4Ions*>(ion)->SetFloatLevelBase(flb);
315
316 // No Anti particle registered
317 ion->SetAntiPDGEncoding(0);
318
319 #ifdef G4VERBOSE
320 if (GetVerboseLevel() > 1) {
321 G4cout << "G4IonTable::CreateIon() : create ion of " << name << " " << Z << ", " << A
322 << " encoding=" << encoding;
323 if (E > 0.0) {
324 G4cout << " IsomerLVL=" << lvl << " excited energy=" << Eex / keV << "[keV]";
325 }
326 G4cout << G4endl;
327 }
328 #endif
329
330 // Add process manager to the ion
331 AddProcessManager(ion);
332
333 #ifdef G4MULTITHREADED
334 // Fill decay channels if this method is invoked from worker
335 if (G4Threading::IsWorkerThread()) {
336 if (!stable && (decayTable != nullptr)) {
337 G4int nCh = decayTable->entries();
338 for (G4int iCh = 0; iCh < nCh; ++iCh) {
339 decayTable->GetDecayChannel(iCh)->GetDaughter(0);
340 }
341 }
342 }
343 #endif
344
345 return ion;
346 }
347
348 G4ParticleDefinition* G4IonTable::CreateIon(G4int Z, G4int A, G4int LL, G4double E,
349 G4Ions::G4FloatLevelBase flb)
350 {
351 if (LL == 0) return CreateIon(Z, A, E, flb);
352
353 // create hyper nucleus
354 G4ParticleDefinition* ion = nullptr;
355
356 // check whether GenericIon has processes
357 G4ParticleDefinition* genericIon = G4ParticleTable::GetParticleTable()->GetGenericIon();
358 G4ProcessManager* pman = nullptr;
359 if (genericIon != nullptr) pman = genericIon->GetProcessManager();
360 if ((genericIon == nullptr) || (genericIon->GetParticleDefinitionID() < 0) || (pman == nullptr)) {
361 #ifdef G4VERBOSE
362 if (GetVerboseLevel() > 1) {
363 G4cout << "G4IonTable::CreateIon() : can not create ion of "
364 << " Z =" << Z << " A = " << A << " because GenericIon is not ready !!" << G4endl;
365 }
366 #endif
367 G4Exception("G4IonTable::CreateIon()", "PART105", JustWarning,
368 "Can not create ions because GenericIon is not ready");
369 return nullptr;
370 }
371
372 G4int J = 0;
373 G4double life = 0.0;
374 G4DecayTable* decayTable = nullptr;
375 G4bool stable = true;
376
377 // excitation energy
378 // G4double Eex = G4NuclideTable::Round(E);
379 G4double Eex = E;
380 G4double mass = GetNucleusMass(Z, A, LL) + Eex;
381 G4int lvl = 0;
382 // lvl is assigned to 9 temporarily
383 if (Eex > 0.0) lvl = 9;
384
385 // PDG encoding
386 G4int encoding = GetNucleusEncoding(Z, A, LL, E, lvl);
387
388 // PDG charge is set to one of nucleus
389 G4double charge = G4double(Z) * eplus;
390
391 // create an ion
392 // spin, parity, isospin values are fixed
393 //
394 // get ion name
395 G4String name = GetIonName(Z, A, LL, Eex, flb);
396
397 // clang-format off
398 ion = new G4Ions( name, mass, 0.0*MeV, charge,
399 J, +1, 0,
400 0, 0, 0,
401 "nucleus", 0, A, encoding,
402 stable, life, decayTable, false,
403 "generic", 0,
404 Eex, lvl );
405 // clang-format on
406
407 // Release lock for particle table accesses
408
409 G4double mu = 0.0; // magnetic moment
410 ion->SetPDGMagneticMoment(mu);
411 static_cast<G4Ions*>(ion)->SetFloatLevelBase(flb);
412
413 // No Anti particle registered
414 ion->SetAntiPDGEncoding(0);
415
416 #ifdef G4VERBOSE
417 if (GetVerboseLevel() > 1) {
418 G4cout << "G4IonTable::CreateIon() : create hyper ion of " << name << " " << Z << ", " << A
419 << ", " << LL << " encoding=" << encoding;
420 if (E > 0.0) {
421 G4cout << " IsomerLVL=" << lvl << " excited energy=" << Eex / keV << "[keV]";
422 }
423 G4cout << G4endl;
424 }
425 #endif
426
427 // Add process manager to the ion
428 AddProcessManager(ion);
429
430 return ion;
431 }
432
433 G4ParticleDefinition* G4IonTable::CreateIon(G4int Z, G4int A, G4int lvl)
434 {
435 // always create an ion for any lvl
436 return CreateIon(Z, A, 0.0, G4Ions::FloatLevelBase(lvl));
437 }
438
439 G4ParticleDefinition* G4IonTable::CreateIon(G4int Z, G4int A, G4int LL, G4int lvl)
440 {
441 return (LL == 0) ? CreateIon(Z, A, 0.0, G4Ions::G4FloatLevelBase(lvl))
442 : CreateIon(Z, A, LL, 0.0, G4Ions::G4FloatLevelBase::no_Float);
443 }
444
445 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4int lvl)
446 {
447 return GetIon(Z, A, 0.0, G4Ions::G4FloatLevelBase(lvl), 0);
448 }
449
450 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4int LL, G4int lvl)
451 {
452 return (LL == 0) ? GetIon(Z, A, 0.0, G4Ions::G4FloatLevelBase(lvl), 0)
453 : GetIon(Z, A, LL, 0.0, G4Ions::G4FloatLevelBase::no_Float, 0);
454 }
455
456 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4double E, G4int J)
457 {
458 return GetIon(Z, A, E, G4Ions::G4FloatLevelBase::no_Float, J);
459 }
460
461 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4double E, char flbChar, G4int J)
462 {
463 return GetIon(Z, A, E, G4Ions::FloatLevelBase(flbChar), J);
464 }
465
466 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4double E, G4Ions::G4FloatLevelBase flb,
467 G4int J)
468 {
469 if ((A < 1) || (Z <= 0) || (E < 0.0) || (A > 999) || (J < 0)) {
470 #ifdef G4VERBOSE
471 if (GetVerboseLevel() > 0) {
472 G4cout << "G4IonTable::GetIon() : illegal atomic number/mass"
473 << " Z =" << Z << " A = " << A << " E = " << E / keV << G4endl;
474 }
475 #endif
476 return nullptr;
477 }
478 auto flb1 = flb;
479
480 // Search ions with A, Z
481 G4ParticleDefinition* ion = FindIon(Z, A, E, flb, J);
482
483 // find out ground state floating level
484 if (ion == nullptr && E == 0.0) {
485 const G4IsotopeProperty* fProperty = FindIsotope(Z, A, E, flb);
486 if (nullptr != fProperty) {
487 flb1 = fProperty->GetFloatLevelBase();
488 if (flb != flb1) {
489 ion = FindIon(Z, A, E, flb1, J);
490 }
491 }
492 }
493
494 // create ion
495 #ifdef G4MULTITHREADED
496 if (ion == nullptr) {
497 if (G4Threading::IsWorkerThread()) {
498 G4MUTEXLOCK(&G4IonTable::ionTableMutex);
499 ion = FindIonInMaster(Z, A, E, flb1, J);
500 if (ion == nullptr) ion = CreateIon(Z, A, E, flb1);
501 InsertWorker(ion);
502 G4MUTEXUNLOCK(&G4IonTable::ionTableMutex);
503 }
504 else {
505 ion = CreateIon(Z, A, E, flb1);
506 }
507 }
508 #else
509 if (ion == nullptr) ion = CreateIon(Z, A, E, flb1);
510 #endif
511
512 return ion;
513 }
514
515 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4int LL, G4double E, G4int J)
516 {
517 return GetIon(Z, A, LL, E, G4Ions::G4FloatLevelBase::no_Float, J);
518 }
519
520 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4int LL, G4double E, char flbChar,
521 G4int J)
522 {
523 return GetIon(Z, A, LL, E, G4Ions::FloatLevelBase(flbChar), J);
524 }
525
526 G4ParticleDefinition* G4IonTable::GetIon(G4int Z, G4int A, G4int LL, G4double E,
527 G4Ions::G4FloatLevelBase flb, G4int J)
528 {
529 if (LL == 0) return GetIon(Z, A, E, flb, J);
530
531 if (A < 2 || Z < 0 || Z > A - LL || LL > A || A > 999) {
532 #ifdef G4VERBOSE
533 if (GetVerboseLevel() > 0) {
534 G4cout << "G4IonTable::GetIon() : illegal atomic number/mass"
535 << " Z =" << Z << " A = " << A << " L = " << LL << " E = " << E / keV << G4endl;
536 }
537 #endif
538 return nullptr;
539 }
540 if (A == 2) {
541 #ifdef G4VERBOSE
542 if (GetVerboseLevel() > 0) {
543 G4cout << "G4IonTable::GetIon() : No boud state for "
544 << " Z =" << Z << " A = " << A << " L = " << LL << " E = " << E / keV << G4endl;
545 }
546 #endif
547 return nullptr;
548 }
549
550 // Search ions with A, Z
551 G4ParticleDefinition* ion = FindIon(Z, A, LL, E, flb, J);
552
553 // create ion
554 #ifdef G4MULTITHREADED
555 if (ion == nullptr) {
556 if (G4Threading::IsWorkerThread()) {
557 G4MUTEXLOCK(&G4IonTable::ionTableMutex);
558 ion = FindIonInMaster(Z, A, LL, E, flb, J);
559 if (ion == nullptr) ion = CreateIon(Z, A, LL, E, flb);
560 InsertWorker(ion);
561 G4MUTEXUNLOCK(&G4IonTable::ionTableMutex);
562 }
563 else {
564 ion = CreateIon(Z, A, LL, E, flb);
565 }
566 }
567 #else
568 if (ion == nullptr) ion = CreateIon(Z, A, LL, E, flb);
569 #endif
570
571 return ion;
572 }
573
574 G4ParticleDefinition* G4IonTable::GetIon(G4int encoding)
575 {
576 G4int Z, A, LL, IsoLvl;
577 G4double E;
578 if (!GetNucleusByEncoding(encoding, Z, A, LL, E, IsoLvl)) {
579 #ifdef G4VERBOSE
580 if (GetVerboseLevel() > 0) {
581 G4cout << "G4IonTable::GetIon() : illegal encoding"
582 << " CODE:" << encoding << G4endl;
583 }
584 #endif
585 G4Exception("G4IonTable::GetIon()", "PART106", JustWarning, "illegal encoding for an ion");
586 return nullptr;
587 }
588 return GetIon(Z, A, LL, IsoLvl);
589 }
590
591 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4double E, G4int J)
592 {
593 return FindIon(Z, A, E, G4Ions::G4FloatLevelBase::no_Float, J);
594 }
595
596 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4double E, char flbChar, G4int J)
597 {
598 return FindIon(Z, A, E, G4Ions::FloatLevelBase(flbChar), J);
599 }
600
601 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4double E,
602 G4Ions::G4FloatLevelBase flb, G4int J)
603 {
604 if ((A < 1) || (Z <= 0) || (J < 0) || (E < 0.0) || (A > 999)) {
605 #ifdef G4VERBOSE
606 if (GetVerboseLevel() > 0) {
607 G4cout << "G4IonTable::FindIon(): illegal atomic number/mass"
608 << " or excitation level:" << G4endl << " Z =" << Z << " A = " << A
609 << " E = " << E / keV << G4endl;
610 }
611 #endif
612 G4Exception("G4IonTable::FindIon()", "PART107", JustWarning, "illegal atomic number/mass");
613 return nullptr;
614 }
615 // Search ions with A, Z ,E
616 // !! J is omitted now !!
617 const G4ParticleDefinition* ion = nullptr;
618 G4bool isFound = false;
619
620 // check if light ion
621 ion = GetLightIon(Z, A);
622 if (ion != nullptr && E == 0.0) {
623 // light ion
624 isFound = true;
625 }
626 else {
627 // -- loop over all particles in Ion table
628 G4int encoding = GetNucleusEncoding(Z, A);
629 const G4ParticleDefinition* ion1 = nullptr;
630 for (auto i = fIonList->find(encoding); i != fIonList->cend(); ++i) {
631 ion = i->second;
632 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
633 // excitation level
634 G4double anExcitaionEnergy = ((const G4Ions*)(ion))->GetExcitationEnergy();
635
636 if (std::fabs(E - anExcitaionEnergy) < pNuclideTable->GetLevelTolerance()) {
637 if (nullptr == ion1) ion1 = ion;
638 if (((const G4Ions*)(ion))->GetFloatLevelBase() == flb) {
639 isFound = true;
640 break;
641 }
642 }
643 }
644 // rerun search on ground level without check on floating
645 if (!isFound && E == 0.0 && nullptr != ion1) {
646 isFound = true;
647 ion = ion1;
648 }
649 }
650
651 if (isFound) {
652 return const_cast<G4ParticleDefinition*>(ion);
653 }
654
655 return nullptr;
656 }
657
658 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4int LL, G4double E, G4int J)
659 {
660 return (LL == 0) ? FindIon(Z, A, E, G4Ions::G4FloatLevelBase::no_Float, J)
661 : FindIon(Z, A, LL, E, G4Ions::G4FloatLevelBase::no_Float, J);
662 }
663
664 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4int LL, G4double E, char flbChar,
665 G4int J)
666 {
667 return FindIon(Z, A, LL, E, G4Ions::FloatLevelBase(flbChar), J);
668 }
669
670 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4int LL, G4double E,
671 G4Ions::G4FloatLevelBase flb, G4int J)
672 {
673 if (LL == 0) return FindIon(Z, A, E, flb, J);
674
675 if (A < 2 || Z < 0 || Z > A - LL || LL > A || A > 999) {
676 #ifdef G4VERBOSE
677 if (GetVerboseLevel() > 0) {
678 G4cout << "G4IonTable::FindIon(): illegal atomic number/mass"
679 << " or excitation level:" << G4endl << " Z =" << Z << " A = " << A << " L = " << LL
680 << " E = " << E / keV << G4endl;
681 }
682 #endif
683 G4Exception("G4IonTable::FindIon()", "PART107", JustWarning, "illegal atomic number/mass");
684 return nullptr;
685 }
686 // Search ions with A, Z ,E
687 // !! J is omitted now !!
688 const G4ParticleDefinition* ion = nullptr;
689 G4bool isFound = false;
690
691 // -- loop over all particles in Ion table
692 G4int encoding = GetNucleusEncoding(Z, A, LL, 0.0, 0);
693 for (auto i = fIonList->find(encoding); i != fIonList->cend(); ++i) {
694 ion = i->second;
695 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
696 if (ion->GetQuarkContent(3) != LL) break;
697 // excitation level
698 G4double anExcitaionEnergy = ((const G4Ions*)(ion))->GetExcitationEnergy();
699 if (std::fabs(E - anExcitaionEnergy) < pNuclideTable->GetLevelTolerance()) {
700 if (((const G4Ions*)(ion))->GetFloatLevelBase() == flb) {
701 isFound = true;
702 break;
703 }
704 }
705 }
706
707 if (isFound) {
708 return const_cast<G4ParticleDefinition*>(ion);
709 }
710
711 return nullptr;
712 }
713
714 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4int lvl)
715 {
716 return FindIon(Z, A, 0.0, G4Ions::FloatLevelBase(lvl), 0);
717 }
718
719 G4ParticleDefinition* G4IonTable::FindIon(G4int Z, G4int A, G4int LL, G4int lvl)
720 {
721 return (LL == 0) ? FindIon(Z, A, 0.0, G4Ions::FloatLevelBase(lvl), 0)
722 : FindIon(Z, A, LL, 0.0, G4Ions::G4FloatLevelBase::no_Float, 0);
723 }
724
725 G4int G4IonTable::GetNucleusEncoding(G4int Z, G4int A, G4double E, G4int lvl)
726 {
727 // PDG code for Ions
728 // Nuclear codes are given as 10-digit numbers +-100ZZZAAAI.
729 // For a nucleus consisting of np protons and nn neutrons
730 // A = np + nn and Z = np.
731 // I gives the isomer level, with I = 0 corresponding
732 // to the ground state and I >0 to excitations
733
734 if (Z == 1 && A == 1 && E == 0.0) return 2212; // proton
735
736 G4int encoding = 1000000000;
737 encoding += Z * 10000;
738 encoding += A * 10;
739 if (lvl > 0 && lvl < 10)
740 encoding += lvl; // isomer level
741 else if (E > 0.0)
742 encoding += 9; // isomer level
743
744 return encoding;
745 }
746
747 G4int G4IonTable::GetNucleusEncoding(G4int Z, G4int A, G4int LL, G4double E, G4int lvl)
748 {
749 // Get PDG code for Hyper-Nucleus Ions
750 // Nuclear codes are given as 10-digit numbers +-10LZZZAAAI.
751 // For a nucleus consisting of np protons and nn neutrons
752 // A = np + nn +nlambda and Z = np.
753 // LL = nlambda
754 // I gives the isomer level, with I = 0 corresponding
755 // to the ground state and I >0 to excitations
756
757 G4int encoding = GetNucleusEncoding(Z, A, E, lvl);
758 if (LL == 0) return encoding;
759 encoding += LL * 10000000;
760 if (Z == 1 && A == 1 && E == 0.0) encoding = 3122; // Lambda
761
762 return encoding;
763 }
764
765 G4bool G4IonTable::GetNucleusByEncoding(G4int encoding, G4int& Z, G4int& A, G4double& E, G4int& lvl)
766 {
767 if (encoding <= 0) return false; // anti particle
768
769 if (encoding == 2212) // proton
770 {
771 Z = 1;
772 A = 1;
773 E = 0.0;
774 lvl = 0;
775 return true;
776 }
777
778 encoding -= 1000000000;
779 Z = encoding / 10000;
780 encoding -= 10000 * Z;
781 A = encoding / 10;
782 lvl = encoding % 10;
783 return true;
784 }
785
786 G4bool G4IonTable::GetNucleusByEncoding(G4int encoding, G4int& Z, G4int& A, G4int& LL, G4double& E,
787 G4int& lvl)
788 {
789 if (encoding <= 0) return false; // anti particle
790
791 if (encoding == 3122) // Lambda
792 {
793 Z = 1;
794 A = 1;
795 LL = 1;
796 E = 0.0;
797 lvl = 0;
798 return true;
799 }
800
801 if (encoding % 10 != 0) {
802 // !!!not supported for excitation states !!!
803 return false;
804 }
805 if (encoding < 1000000000) {
806 // anti particle
807 return false;
808 }
809
810 encoding -= 1000000000;
811 LL = encoding / 10000000;
812 encoding -= 10000000 * LL;
813 Z = encoding / 10000;
814 encoding -= 10000 * Z;
815 A = encoding / 10;
816 lvl = encoding % 10;
817 return true;
818 }
819
820 G4String G4IonTable::GetIonName(G4int Z, G4int A, G4double E,
821 G4Ions::G4FloatLevelBase flb) const
822 {
823 G4String name = GetIonName(Z, A, 0);
824
825 // Excited energy or floating level
826 if (E > 0 || flb != G4Ions::G4FloatLevelBase::no_Float) {
827 std::ostringstream os;
828 os.setf(std::ios::fixed);
829 os.precision(3);
830 // Excited nucleus
831 os << '[' << E / keV;
832 if (flb != G4Ions::G4FloatLevelBase::no_Float) {
833 os << G4Ions::FloatLevelBaseChar(flb);
834 }
835 os << ']';
836 name += os.str();
837 }
838
839 return name;
840 }
841
842 G4String G4IonTable::GetIonName(G4int Z, G4int A, G4int LL, G4double E,
843 G4Ions::G4FloatLevelBase flb) const
844 {
845 if (LL == 0) return GetIonName(Z, A, E, flb);
846 G4String name = "";
847 for (G4int i = 0; i < LL; ++i) {
848 name += "L";
849 }
850 name += GetIonName(Z, A, E, flb);
851 return name;
852 }
853
854 G4String G4IonTable::GetIonName(G4int Z, G4int A, G4int lvl) const
855 {
856 std::ostringstream os;
857
858 // Atomic number
859 if ((0 < Z) && (Z <= numberOfElements)) {
860 os << elementName[Z - 1];
861 }
862 else {
863 os << "E" << Z << "-";
864 }
865 // Atomic Mass
866 os << A;
867
868 if (lvl > 0) {
869 // Isomer level for Excited nucelus
870 os << '[' << lvl << ']';
871 }
872 G4String name = os.str();
873 return name;
874 }
875
876 G4String G4IonTable::GetIonName(G4int Z, G4int A, G4int LL, G4int lvl) const
877 {
878 if (LL == 0) return GetIonName(Z, A, lvl);
879 G4String name = "";
880 for (G4int i = 0; i < LL; ++i) {
881 name += "L";
882 }
883 name += GetIonName(Z, A, lvl);
884 return name;
885 }
886
887 G4bool G4IonTable::IsIon(const G4ParticleDefinition* particle)
888 {
889 // Return true if the particle is ion
890 static const G4String nucleus("nucleus");
891 static const G4String proton("proton");
892
893 // Neutron is not ion
894 if ((particle->GetAtomicMass() > 0) && (particle->GetAtomicNumber() > 0)) {
895 return particle->GetBaryonNumber() > 0;
896 }
897
898 // Particles derived from G4Ions
899 if (particle->GetParticleType() == nucleus) return true;
900
901 // Proton (Hydrogen nucleus)
902 if (particle->GetParticleName() == proton) return true;
903
904 return false;
905 }
906
907 G4bool G4IonTable::IsAntiIon(const G4ParticleDefinition* particle)
908 {
909 // Return true if the particle is ion
910 static const G4String anti_nucleus("anti_nucleus");
911 static const G4String anti_proton("anti_proton");
912
913 // Anti_neutron is not ion
914 if ((particle->GetAtomicMass() > 0) && (particle->GetAtomicNumber() > 0)) {
915 return particle->GetBaryonNumber() < 0;
916 }
917
918 // Particles derived from G4Ions
919 if (particle->GetParticleType() == anti_nucleus) return true;
920
921 // Anti_proton (Anti_Hydrogen nucleus)
922 if (particle->GetParticleName() == anti_proton) return true;
923
924 return false;
925 }
926
927 G4bool G4IonTable::IsLightIon(const G4ParticleDefinition* particle) const
928 {
929 static const std::string names[] = {"proton", "alpha", "deuteron", "triton", "He3"};
930
931 // Return true if the particle is pre-defined ion
932 return std::find(names, names + 5, (particle->GetParticleName()).c_str()) != names + 5;
933 }
934
935 G4bool G4IonTable::IsLightAntiIon(const G4ParticleDefinition* particle) const
936 {
937 static const std::string names[] = {"anti_proton", "anti_alpha", "anti_deuteron", "anti_triton",
938 "anti_He3"};
939
940 // Return true if the particle is pre-defined ion
941 return std::find(names, names + 5, (particle->GetParticleName()).c_str()) != names + 5;
942 }
943
944 G4ParticleDefinition* G4IonTable::GetLightIon(G4int Z, G4int A) const
945 {
946 // Returns pointer to pre-defined ions
947 const G4ParticleDefinition* ion = nullptr;
948 if ((Z <= 2)) {
949 #ifndef G4MULTITHREADED
950 // In sequential use lazy-initialization
951 lightions::Init();
952 #endif
953 if ((Z == 1) && (A == 1)) {
954 ion = lightions::p_proton;
955 }
956 else if ((Z == 1) && (A == 2)) {
957 ion = lightions::p_deuteron;
958 }
959 else if ((Z == 1) && (A == 3)) {
960 ion = lightions::p_triton;
961 }
962 else if ((Z == 2) && (A == 4)) {
963 ion = lightions::p_alpha;
964 }
965 else if ((Z == 2) && (A == 3)) {
966 ion = lightions::p_He3;
967 }
968 }
969 return const_cast<G4ParticleDefinition*>(ion);
970 }
971
972 G4ParticleDefinition* G4IonTable::GetLightAntiIon(G4int Z, G4int A) const
973 {
974 // Returns pointer to pre-defined ions
975 const G4ParticleDefinition* ion = nullptr;
976 if ((Z <= 2)) {
977 #ifndef G4MULTITHREADED
978 // In sequential use lazy-initialization
979 antilightions::Init();
980 #endif
981 if ((Z == 1) && (A == 1)) {
982 ion = antilightions::p_proton;
983 }
984 else if ((Z == 1) && (A == 2)) {
985 ion = antilightions::p_deuteron;
986 }
987 else if ((Z == 1) && (A == 3)) {
988 ion = antilightions::p_triton;
989 }
990 else if ((Z == 2) && (A == 4)) {
991 ion = antilightions::p_alpha;
992 }
993 else if ((Z == 2) && (A == 3)) {
994 ion = antilightions::p_He3;
995 }
996 }
997 return const_cast<G4ParticleDefinition*>(ion);
998 }
999
1000 G4double G4IonTable::GetNucleusMass(G4int Z, G4int A, G4int LL, G4int lvl) const
1001 {
1002 if ((A < 1) || (Z < 0) || (LL < 0) || (lvl < 0) || (lvl > 9)) {
1003 #ifdef G4VERBOSE
1004 if (GetVerboseLevel() > 0) {
1005 G4cout << "G4IonTable::GetNucleusMass() : illegal atomic number/mass:" << G4endl
1006 << " Z =" << Z << " A = " << A << " L = " << LL << " lvl = " << lvl << G4endl;
1007 }
1008 #endif
1009 G4Exception("G4IonTable::GetNucleusMass()", "PART107", EventMustBeAborted,
1010 "illegal atomic number/mass");
1011 return -1.0;
1012 }
1013
1014 G4double mass;
1015 if (LL == 0) {
1016 // calculate nucleus mass
1017 const G4ParticleDefinition* ion = GetLightIon(Z, A);
1018
1019 if (ion != nullptr) {
1020 mass = ion->GetPDGMass();
1021 }
1022 else {
1023 // Use G4NucleiProperties::GetNuclearMass
1024 mass = G4NucleiProperties::GetNuclearMass(A, Z);
1025 }
1026
1027 // Isomer
1028 if (lvl > 0) {
1029 // -- loop over all particles in Ion table
1030 G4int encoding = GetNucleusEncoding(Z, A);
1031 G4bool isFound = false;
1032 for (auto i = fIonList->find(encoding); i != fIonList->cend(); ++i) {
1033 ion = i->second;
1034 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
1035 // Excitation level
1036 if (((const G4Ions*)(ion))->GetIsomerLevel() == lvl) {
1037 isFound = true;
1038 break;
1039 }
1040 }
1041 if (isFound) {
1042 // Return existing isomer mass
1043 mass = ion->GetPDGMass();
1044 }
1045 else {
1046 // Find isomer from IsotopeTable
1047 const G4IsotopeProperty* fProperty = FindIsotope(Z, A, lvl);
1048 if (fProperty != nullptr) mass += fProperty->GetEnergy();
1049 }
1050 }
1051 }
1052 else {
1053 mass = G4HyperNucleiProperties::GetNuclearMass(A, Z, LL);
1054 }
1055 return mass;
1056 }
1057
1058 G4double G4IonTable::GetIsomerMass(G4int Z, G4int A, G4int lvl) const
1059 {
1060 return GetNucleusMass(Z, A, 0, lvl);
1061 }
1062
1063 G4double G4IonTable::GetIonMass(G4int Z, G4int A, G4int LL, G4int lvl) const
1064 {
1065 return GetNucleusMass(Z, A, LL, lvl);
1066 }
1067
1068 void G4IonTable::clear()
1069 {
1070 if (G4ParticleTable::GetParticleTable()->GetReadiness()) {
1071 G4Exception("G4IonTable::clear()", "PART116", JustWarning,
1072 "No effects because readyToUse is true.");
1073 return;
1074 }
1075
1076 #ifdef G4VERBOSE
1077 if (GetVerboseLevel() > 2) {
1078 G4cout << "G4IonTable::Clear() : number of Ion registered = ";
1079 G4cout << fIonList->size() << G4endl;
1080 }
1081 #endif
1082 fIonList->clear();
1083 }
1084
1085 void G4IonTable::Insert(const G4ParticleDefinition* particle)
1086 {
1087 if (!IsIon(particle)) return;
1088 if (Contains(particle)) return;
1089
1090 G4int Z = particle->GetAtomicNumber();
1091 G4int A = particle->GetAtomicMass();
1092 G4int LL = particle->GetQuarkContent(3); // strangeness
1093 G4int encoding = GetNucleusEncoding(Z, A, LL); // encoding of the groud state
1094
1095 // Register the ion with its encoding of the ground state
1096 fIonListShadow->insert(std::pair<const G4int, const G4ParticleDefinition*>(encoding, particle));
1097 }
1098
1099 void G4IonTable::InsertWorker(const G4ParticleDefinition* particle)
1100 {
1101 if (particle == nullptr) return;
1102
1103 G4int Z = particle->GetAtomicNumber();
1104 G4int A = particle->GetAtomicMass();
1105 G4int LL = particle->GetQuarkContent(3); // strangeness
1106 G4int encoding = GetNucleusEncoding(Z, A, LL);
1107 G4bool found = false;
1108 if (encoding != 0) {
1109 for (auto i = fIonList->find(encoding); i != fIonList->cend(); ++i) {
1110 if (particle == i->second) {
1111 found = true;
1112 break;
1113 }
1114 }
1115 }
1116 if (found) return;
1117
1118 // Register the ion with its encoding of the gronud state
1119 fIonList->insert(std::pair<const G4int, const G4ParticleDefinition*>(encoding, particle));
1120 }
1121
1122 void G4IonTable::Remove(const G4ParticleDefinition* particle)
1123 {
1124 if (particle == nullptr) return;
1125 #ifdef G4MULTITHREADED
1126 if (G4Threading::IsWorkerThread()) {
1127 G4ExceptionDescription ed;
1128 ed << "Request of removing " << particle->GetParticleName()
1129 << " is ignored as it is invoked from a worker thread.";
1130 G4Exception("G4IonTable::Remove()", "PART10117", JustWarning, ed);
1131 return;
1132 }
1133 #endif
1134 if (G4ParticleTable::GetParticleTable()->GetReadiness()) {
1135 G4StateManager* pStateManager = G4StateManager::GetStateManager();
1136 G4ApplicationState currentState = pStateManager->GetCurrentState();
1137 if (currentState != G4State_PreInit) {
1138 G4String msg = "Request of removing ";
1139 msg += particle->GetParticleName();
1140 msg += " has No effects other than Pre_Init";
1141 G4Exception("G4IonTable::Remove()", "PART117", JustWarning, msg);
1142 return;
1143 }
1144
1145 #ifdef G4VERBOSE
1146 if (GetVerboseLevel() > 0) {
1147 G4cout << particle->GetParticleName() << " will be removed from the IonTable " << G4endl;
1148 }
1149 #endif
1150 }
1151
1152 if (IsIon(particle)) {
1153 G4int Z = particle->GetAtomicNumber();
1154 G4int A = particle->GetAtomicMass();
1155 G4int LL = particle->GetQuarkContent(3); // strangeness
1156 G4int encoding = GetNucleusEncoding(Z, A, LL);
1157 if (encoding != 0) {
1158 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1159 if (particle == i->second) {
1160 fIonListShadow->erase(i);
1161 break;
1162 }
1163 }
1164 }
1165 }
1166 else {
1167 #ifdef G4VERBOSE
1168 if (GetVerboseLevel() > 1) {
1169 G4cout << "G4IonTable::Remove :" << particle->GetParticleName() << " is not ions" << G4endl;
1170 }
1171 #endif
1172 }
1173 }
1174
1175 void G4IonTable::DumpTable(const G4String& particle_name) const
1176 {
1177 const G4ParticleDefinition* ion;
1178 for (const auto& idx : *fIonList) {
1179 ion = idx.second;
1180 if ((particle_name == "ALL") || (particle_name == "all")) {
1181 ion->DumpTable();
1182 }
1183 else if (particle_name == ion->GetParticleName()) {
1184 ion->DumpTable();
1185 }
1186 }
1187 }
1188
1189 // --------------------------------------------------------------------
1190 //
1191 // clang-format off
1192 const G4String G4IonTable::elementName[] =
1193 {
1194 "H", "He",
1195 "Li", "Be", "B", "C", "N", "O", "F", "Ne",
1196 "Na", "Mg", "Al", "Si", "P", "S", "Cl", "Ar",
1197 "K", "Ca", "Sc", "Ti", "V", "Cr", "Mn", "Fe", "Co", "Ni", "Cu", "Zn", "Ga", "Ge", "As", "Se", "Br", "Kr",
1198 "Rb", "Sr", "Y", "Zr", "Nb", "Mo","Tc", "Ru", "Rh", "Pd", "Ag", "Cd", "In", "Sn", "Sb", "Te", "I", "Xe",
1199 "Cs", "Ba",
1200 "La", "Ce", "Pr", "Nd", "Pm", "Sm", "Eu", "Gd", "Tb", "Dy", "Ho", "Er", "Tm", "Yb", "Lu",
1201 "Hf", "Ta", "W", "Re", "Os", "Ir", "Pt", "Au", "Hg", "Tl", "Pb", "Bi", "Po", "At", "Rn",
1202 "Fr", "Ra",
1203 "Ac", "Th", "Pa", "U", "Np", "Pu", "Am", "Cm", "Bk", "Cf", "Es", "Fm", "Md", "No", "Lr",
1204 "Rf", "Db", "Sg", "Bh", "Hs", "Mt", "Ds", "Rg", "Cn", "Nh", "Fl", "Mc", "Lv", "Ts", "Og"
1205 };
1206 // clang-format on
1207
1208 G4int G4IonTable::GetVerboseLevel() const
1209 {
1210 return G4ParticleTable::GetParticleTable()->GetVerboseLevel();
1211 }
1212
1213 void G4IonTable::AddProcessManager(G4ParticleDefinition* ion)
1214 {
1215 if (ion->IsGeneralIon()) {
1216 // Check whether GenericIon has processes
1217 G4ParticleDefinition* genericIon = G4ParticleTable::GetParticleTable()->GetGenericIon();
1218
1219 G4ProcessManager* pman = nullptr;
1220 if (genericIon != nullptr) pman = genericIon->GetProcessManager();
1221 if ((genericIon == nullptr) || (genericIon->GetParticleDefinitionID() < 0) || (pman == nullptr))
1222 {
1223 G4String msg = "G4IonTable::AddProcessManager(): cannot create ion of ";
1224 msg += ion->GetParticleName();
1225 msg += "\n because GenericIon is not available!!";
1226 G4Exception("G4IonTable::AddProcessManager()", "PART105", FatalException, msg);
1227 return;
1228 }
1229
1230 ion->SetParticleDefinitionID(genericIon->GetParticleDefinitionID());
1231 }
1232 else {
1233 // Is this a MuonicAtom ?
1234 auto muatom = dynamic_cast<G4MuonicAtom*>(ion);
1235
1236 if (muatom != nullptr) {
1237 #ifdef G4VERBOSE
1238 if (GetVerboseLevel() > 1) {
1239 G4cout << "G4IonTable::AddProcessManager(): "
1240 << "MuonicAtom dynamic_cast succeeded for " << ion->GetParticleName() << G4endl;
1241 }
1242 #endif
1243 // Check whether GenericMuonicAtom has processes
1244 G4ParticleDefinition* genericMA = G4ParticleTable::GetParticleTable()->GetGenericMuonicAtom();
1245
1246 G4ProcessManager* pman = nullptr;
1247 if (genericMA != nullptr) pman = genericMA->GetProcessManager();
1248 if ((genericMA == nullptr) || (genericMA->GetParticleDefinitionID() < 0) || (pman == nullptr))
1249 {
1250 G4String msg = "G4IonTable::AddProcessManager(): cannot create MuonicAtom ";
1251 msg += ion->GetParticleName();
1252 msg += "\n because GenericMuonicAtom is not available!!";
1253 G4Exception("G4IonTable::AddProcessManager()", "PART106", FatalException, msg);
1254 return;
1255 }
1256
1257 ion->SetParticleDefinitionID(genericMA->GetParticleDefinitionID());
1258 }
1259 else {
1260 G4String msg = "G4IonTable::AddProcessManager(): cannot create ";
1261 msg += ion->GetParticleName();
1262 msg += "\n because of unsupported particle type !!";
1263 G4Exception("G4IonTable::AddProcessManager()", "PART107", FatalException, msg);
1264 return;
1265 }
1266 }
1267 return;
1268 }
1269
1270 void G4IonTable::RegisterIsotopeTable(G4VIsotopeTable* table)
1271 {
1272 // check duplication
1273 G4String name = table->GetName();
1274 for (const auto fIsotopeTable : *fIsotopeTableList) {
1275 if (name == fIsotopeTable->GetName()) return;
1276 }
1277 // register
1278 fIsotopeTableList->push_back(table);
1279 }
1280
1281 G4VIsotopeTable* G4IonTable::GetIsotopeTable(std::size_t index) const
1282 {
1283 G4VIsotopeTable* fIsotopeTable = nullptr;
1284 if (index < fIsotopeTableList->size()) {
1285 fIsotopeTable = (*fIsotopeTableList)[index];
1286 }
1287 return fIsotopeTable;
1288 }
1289
1290 G4IsotopeProperty* G4IonTable::FindIsotope(G4int Z, G4int A, G4double E,
1291 G4Ions::G4FloatLevelBase flb) const
1292 {
1293 if (fIsotopeTableList == nullptr) return nullptr;
1294 if (fIsotopeTableList->empty()) return nullptr;
1295
1296 G4IsotopeProperty* property = nullptr;
1297
1298 for (std::size_t i = 0; i < fIsotopeTableList->size(); ++i) {
1299 G4VIsotopeTable* fIsotopeTable = (*fIsotopeTableList)[fIsotopeTableList->size() - i - 1];
1300 property = fIsotopeTable->GetIsotope(Z, A, E, flb);
1301 if (property != nullptr) break;
1302 }
1303
1304 return property;
1305 }
1306
1307 G4IsotopeProperty* G4IonTable::FindIsotope(G4int Z, G4int A, G4int lvl) const
1308 {
1309 if (fIsotopeTableList == nullptr) return nullptr;
1310 if (fIsotopeTableList->empty()) return nullptr;
1311
1312 G4IsotopeProperty* property = nullptr;
1313
1314 // iterate
1315 for (std::size_t i = 0; i < fIsotopeTableList->size(); ++i) {
1316 G4VIsotopeTable* fIsotopeTable = (*fIsotopeTableList)[fIsotopeTableList->size() - i - 1];
1317 property = fIsotopeTable->GetIsotope(Z, A, lvl);
1318 if (property != nullptr) break;
1319 }
1320
1321 return property;
1322 }
1323
1324 void G4IonTable::CreateAllIon()
1325 {
1326 PreloadNuclide();
1327 }
1328
1329 void G4IonTable::CreateAllIsomer()
1330 {
1331 PreloadNuclide();
1332 }
1333
1334 void G4IonTable::PrepareNuclideTable()
1335 {
1336 if (pNuclideTable == nullptr) pNuclideTable = G4NuclideTable::GetNuclideTable();
1337 }
1338
1339 void G4IonTable::PreloadNuclide()
1340 {
1341 if (isIsomerCreated || !G4Threading::IsMultithreadedApplication()) return;
1342
1343 pNuclideTable->GenerateNuclide();
1344
1345 for (std::size_t i = 0; i != pNuclideTable->entries(); ++i) {
1346 const G4IsotopeProperty* fProperty = pNuclideTable->GetIsotopeByIndex(i);
1347 G4int Z = fProperty->GetAtomicNumber();
1348 G4int A = fProperty->GetAtomicMass();
1349 G4double Eex = fProperty->GetEnergy();
1350 GetIon(Z, A, Eex);
1351 }
1352
1353 isIsomerCreated = true;
1354 }
1355
1356 G4ParticleDefinition* G4IonTable::GetParticle(G4int index) const
1357 {
1358 if ((index >= 0) && (index < Entries())) {
1359 auto idx = fIonList->cbegin();
1360 G4int counter = 0;
1361 while (idx != fIonList->cend()) // Loop checking, 09.08.2015, K.Kurashige
1362 {
1363 if (counter == index) {
1364 return const_cast<G4ParticleDefinition*>(idx->second);
1365 }
1366 ++counter;
1367 ++idx;
1368 }
1369 }
1370 #ifdef G4VERBOSE
1371 if (GetVerboseLevel() > 1) {
1372 G4cout << " G4IonTable::GetParticle"
1373 << " invalid index (=" << index << ")"
1374 << " entries = " << Entries() << G4endl;
1375 }
1376 #endif
1377 return nullptr;
1378 }
1379
1380 G4bool G4IonTable::Contains(const G4ParticleDefinition* particle) const
1381 {
1382 if (!IsIon(particle)) return false;
1383
1384 G4int Z = particle->GetAtomicNumber();
1385 G4int A = particle->GetAtomicMass();
1386 G4int LL = particle->GetQuarkContent(3); // strangeness
1387 G4int encoding = GetNucleusEncoding(Z, A, LL);
1388 G4bool found = false;
1389 if (encoding != 0) {
1390 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1391 if (particle == i->second) {
1392 found = true;
1393 break;
1394 }
1395 }
1396 }
1397 return found;
1398 }
1399
1400 G4int G4IonTable::Entries() const
1401 {
1402 return (G4int)fIonList->size();
1403 }
1404
1405 G4int G4IonTable::size() const
1406 {
1407 return (G4int)fIonList->size();
1408 }
1409
1410 G4ParticleDefinition* G4IonTable::FindIonInMaster(G4int Z, G4int A, G4double E,
1411 G4Ions::G4FloatLevelBase flb, G4int /*J*/)
1412 {
1413 // Search ions with A, Z ,E
1414 // !! J is omitted now !!
1415 const G4ParticleDefinition* ion = nullptr;
1416 G4bool isFound = false;
1417
1418 // -- loop over all particles in Ion table
1419 G4int encoding = GetNucleusEncoding(Z, A);
1420 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1421 ion = i->second;
1422 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
1423 // excitation level
1424 G4double anExcitaionEnergy = ((const G4Ions*)(ion))->GetExcitationEnergy();
1425 if (std::fabs(E - anExcitaionEnergy) < pNuclideTable->GetLevelTolerance()) {
1426 if (((const G4Ions*)(ion))->GetFloatLevelBase() == flb) {
1427 isFound = true;
1428 break;
1429 }
1430 }
1431 }
1432
1433 if (isFound) {
1434 return const_cast<G4ParticleDefinition*>(ion);
1435 }
1436
1437 return nullptr;
1438 }
1439
1440 G4ParticleDefinition* G4IonTable::FindIonInMaster(G4int Z, G4int A, G4int LL, G4double E,
1441 G4Ions::G4FloatLevelBase flb, G4int J)
1442 {
1443 if (LL == 0) return FindIon(Z, A, E, flb, J);
1444
1445 // Search ions with A, Z ,E
1446 // !! J is omitted now !!
1447 const G4ParticleDefinition* ion = nullptr;
1448 G4bool isFound = false;
1449
1450 // -- loop over all particles in Ion table
1451 G4int encoding = GetNucleusEncoding(Z, A, LL, 0.0, 0);
1452 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1453 ion = i->second;
1454 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
1455 if (ion->GetQuarkContent(3) != LL) break;
1456 // Excitation level
1457 G4double anExcitaionEnergy = ((const G4Ions*)(ion))->GetExcitationEnergy();
1458 if (std::fabs(E - anExcitaionEnergy) < pNuclideTable->GetLevelTolerance()) {
1459 if (((const G4Ions*)(ion))->GetFloatLevelBase() == flb) {
1460 isFound = true;
1461 break;
1462 }
1463 }
1464 }
1465
1466 if (isFound) {
1467 return const_cast<G4ParticleDefinition*>(ion);
1468 }
1469
1470 return nullptr;
1471 }
1472
1473 G4ParticleDefinition* G4IonTable::FindIonInMaster(G4int Z, G4int A, G4int lvl)
1474 {
1475 // Search ions with A, Z ,E
1476 // !! J is omitted now !!
1477 const G4ParticleDefinition* ion = nullptr;
1478 G4bool isFound = false;
1479
1480 // -- loop over all particles in Ion table
1481 G4int encoding = GetNucleusEncoding(Z, A);
1482 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1483 ion = i->second;
1484 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
1485 // Excitation level
1486 if (((const G4Ions*)(ion))->GetIsomerLevel() == lvl) {
1487 isFound = true;
1488 break;
1489 }
1490 }
1491
1492 if (isFound) {
1493 return const_cast<G4ParticleDefinition*>(ion);
1494 }
1495
1496 return nullptr;
1497 }
1498
1499 G4ParticleDefinition* G4IonTable::FindIonInMaster(G4int Z, G4int A, G4int LL, G4int lvl)
1500 {
1501 if (LL == 0) return FindIon(Z, A, lvl);
1502
1503 // Search ions with A, Z ,E, lvl
1504 const G4ParticleDefinition* ion = nullptr;
1505 G4bool isFound = false;
1506
1507 // -- loop over all particles in Ion table
1508 G4int encoding = GetNucleusEncoding(Z, A, LL);
1509 for (auto i = fIonListShadow->find(encoding); i != fIonListShadow->cend(); ++i) {
1510 ion = i->second;
1511 if ((ion->GetAtomicNumber() != Z) || (ion->GetAtomicMass() != A)) break;
1512 if (ion->GetQuarkContent(3) != LL) break;
1513 // excitation level
1514 if (((const G4Ions*)(ion))->GetIsomerLevel() == lvl) {
1515 isFound = true;
1516 break;
1517 }
1518 }
1519
1520 if (isFound) {
1521 return const_cast<G4ParticleDefinition*>(ion);
1522 }
1523
1524 return nullptr;
1525 }
1526
1527 G4double G4IonTable::GetLifeTime(const G4ParticleDefinition* particle) const
1528 {
1529 if ((particle->IsGeneralIon()) && (pNuclideTable == nullptr)) {
1530 G4Exception("G4IonTable::GetLifeTime()", "ParticleIon1001", FatalException,
1531 "Method is invoked before G4IonTable is initialized.");
1532 return 0.;
1533 }
1534 return particle->GetPDGLifeTime();
1535 }
1536
1537 G4double G4IonTable::GetLifeTime(G4int Z, G4int A, G4double E, char flbChar) const
1538 {
1539 return GetLifeTime(Z, A, E, G4Ions::FloatLevelBase(flbChar));
1540 }
1541
1542 G4double G4IonTable::GetLifeTime(G4int Z, G4int A, G4double E, G4Ions::G4FloatLevelBase flb) const
1543 {
1544 G4double life = -1001.0;
1545 const G4IsotopeProperty* fProperty = FindIsotope(Z, A, E, flb);
1546 if (fProperty != nullptr) life = fProperty->GetLifeTime();
1547 return life;
1548 }
1549
1550 G4ParticleDefinition* G4IonTable::GetMuonicAtom(G4Ions const* base)
1551 {
1552 if (base == nullptr || !IsIon(base)) {
1553 G4Exception("G4IonTable::GetMuonicAtom()", "PART987654321", FatalException,
1554 "Constructor argument is not a G4Ions");
1555 return nullptr;
1556 }
1557
1558 // We're assuming here that we get a base that is actually
1559 // constructed and unexcited ... strip excitations, Lambdas, and
1560 // isomers from the encoding
1561
1562 auto const Z = base->GetAtomicNumber();
1563 auto const A = base->GetAtomicMass();
1564 auto const baseenc = GetNucleusEncoding(Z, A);
1565 auto const encoding = baseenc + 1000000000;
1566
1567 // We have to do all the MT manipulations manually, because the
1568 // convenience functions assume a G4Ions with canonical PDG codes;
1569 // they recalculate the encoding from particle properties rather
1570 // than using the carried member function values. Thus, they will
1571 // do operations on the base ion, rather than the passed in
1572 // G4MuonicAtom
1573
1574 auto i = fIonList->find(encoding);
1575 if (i != fIonList->cend()) {
1576 return const_cast<G4ParticleDefinition*>(i->second);
1577 }
1578 // not in threadlocal list; check global list ...
1579 #ifdef G4MULTITHREADED
1580 if (G4Threading::IsWorkerThread()) {
1581 G4MUTEXLOCK(&G4IonTable::ionTableMutex);
1582 i = fIonListShadow->find(encoding);
1583 auto end = fIonListShadow->cend();
1584 G4MUTEXUNLOCK(&G4IonTable::ionTableMutex);
1585 if (i != end) {
1586 // we found it, stuff it into the threadlocal list
1587 fIonList->insert(*i);
1588 // and then return it ...
1589 return const_cast<G4ParticleDefinition*>(i->second);
1590 }
1591 }
1592 #endif
1593
1594 // not found in either list; create and potentially insert
1595 auto const name = "Mu" + GetIonName(Z, A);
1596
1597 G4MuonicAtom* muatom = G4MuonicAtomHelper::ConstructMuonicAtom(name, encoding, base);
1598
1599 // Not sure this is doing the right thing...
1600 AddProcessManager(muatom);
1601
1602 // Now, we have to push the muatom into the appropriate IonTables
1603 // first, recheck global list, in case another thread came along
1604 // before us and created this same muatom
1605
1606 #ifdef G4MULTITHREADED
1607 if (G4Threading::IsWorkerThread()) {
1608 G4MUTEXLOCK(&G4IonTable::ionTableMutex);
1609 // first, we need to make sure it hasn't been inserted by some
1610 // other thread
1611 auto j = fIonListShadow->find(encoding);
1612 if (j != fIonListShadow->cend()) {
1613 // oops ... someone else built a copy when we weren't looking;
1614 // cleanup our instantiation, and take a handle to the one in
1615 // the global list
1616 delete muatom;
1617 muatom = const_cast<G4MuonicAtom*>(static_cast<G4MuonicAtom const*>(j->second));
1618 }
1619 else {
1620 // otherwise, push onto the global list first
1621 fIonListShadow->insert(std::make_pair(encoding, muatom));
1622 }
1623 G4MUTEXUNLOCK(&G4IonTable::ionTableMutex);
1624 }
1625 #endif
1626 // in either case, push onto the the threadlocal list
1627 fIonList->insert(std::make_pair(encoding, muatom));
1628
1629 return muatom;
1630 }
1631
1632 G4ParticleDefinition* G4IonTable::GetMuonicAtom(G4int Z, G4int A)
1633 {
1634 // Need the cast because we need a G4Ions* to pass into the
1635 // function, but GetIon returns a G4ParticleDefinition*
1636 auto base = static_cast<G4Ions const*>(GetIon(Z, A, 0.0));
1637 return GetMuonicAtom(base);
1638 }
1639