Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
Public Member Functions | List of all members
G4NeutronHPPhotonDist Class Reference

#include <G4NeutronHPPhotonDist.hh>

Public Member Functions

 G4NeutronHPPhotonDist ()
 
 ~G4NeutronHPPhotonDist ()
 
G4bool InitMean (std::ifstream &aDataFile)
 
void InitAngular (std::ifstream &aDataFile)
 
void InitEnergies (std::ifstream &aDataFile)
 
void InitPartials (std::ifstream &aDataFile)
 
G4ReactionProductVector * GetPhotons (G4double anEnergy)
 
G4double GetTargetMass ()
 
G4bool NeedsCascade ()
 
G4double GetLevelEnergy ()
 

Detailed Description

Definition at line 55 of file G4NeutronHPPhotonDist.hh.

Constructor & Destructor Documentation

◆ G4NeutronHPPhotonDist()

G4NeutronHPPhotonDist::G4NeutronHPPhotonDist ( )
inline

Definition at line 59 of file G4NeutronHPPhotonDist.hh.

60 : repFlag( 0 )
61 , targetMass( 0.0 )
62 , nDiscrete( 0 )
63 , isoFlag( 0 )
64 , tabulationType( 0 )
65 , nDiscrete2( 0 )
66 , nIso( 0 )
67 , nPartials( 0 )
68 , theInternalConversionFlag( 0 )
69 , nGammaEnergies( 0 )
70 , theBaseEnergy( 0.0 )
71 {
72
73 disType = 0;
74 energy = 0;
75 theYield = 0;
76 thePartialXsec = 0;
77 isPrimary = 0;
78 theShells = 0;
79 theGammas = 0;
80 nNeu = 0;
81 theLegendre = 0;
82 theAngular = 0;
83 distribution = 0;
84 probs = 0;
85 partials = 0;
86 actualMult = 0;
87
88 theLevelEnergies = 0;
89 theTransitionProbabilities = 0;
90 thePhotonTransitionFraction = 0;
91
92 }

◆ ~G4NeutronHPPhotonDist()

G4NeutronHPPhotonDist::~G4NeutronHPPhotonDist ( )
inline

Definition at line 94 of file G4NeutronHPPhotonDist.hh.

95 {
96 delete [] disType;
97 delete [] energy;
98 delete [] theYield;
99 delete [] thePartialXsec;
100 delete [] isPrimary;
101 delete [] theShells;
102 delete [] theGammas;
103 delete [] nNeu;
104 delete [] theAngular;
105 delete [] distribution;
106 delete [] probs;
107
108 if ( theLegendre != NULL )
109 {
110 for ( G4int i = 0 ; i < (nDiscrete2-nIso) ; i++ )
111 if ( theLegendre[i] != NULL ) delete[] theLegendre[i];
112
113 delete [] theLegendre;
114 }
115
116 if ( partials != 0 )
117 {
118 for ( G4int i = 0 ; i < nPartials ; i++ )
119 { delete partials[i]; }
120
121 delete [] partials;
122 }
123
124 delete [] actualMult;
125
126 // delete theLevelEnergies;
127 // delete theTransitionProbabilities;
128 // delete thePhotonTransitionFraction;
129// TKDB
130 delete [] theLevelEnergies;
131 delete [] theTransitionProbabilities;
132 delete [] thePhotonTransitionFraction;
133 }
G4int
int G4int
Definition: G4Types.hh:66

Member Function Documentation

◆ GetLevelEnergy()

G4double G4NeutronHPPhotonDist::GetLevelEnergy ( )
inline

Definition at line 149 of file G4NeutronHPPhotonDist.hh.

149{return theBaseEnergy;}

Referenced by G4NeutronHPInelasticCompFS::CompositeApply().

◆ GetPhotons()

G4ReactionProductVector * G4NeutronHPPhotonDist::GetPhotons ( G4double  anEnergy)

Definition at line 281 of file G4NeutronHPPhotonDist.cc.

282{
283
284 //G4cout << "G4NeutronHPPhotonDist::GetPhotons repFlag " << repFlag << G4endl;
285 // the partial cross-section case is not in this yet. @@@@ << 070601 TK add partial
286 G4int i, ii, iii;
287 G4int nSecondaries = 0;
288 G4ReactionProductVector * thePhotons = new G4ReactionProductVector;
289 if(repFlag==1)
290 {
291 G4double current=0;
292 for(i=0; i<nDiscrete; i++)
293 {
294 current = theYield[i].GetY(anEnergy);
295 actualMult[i] = G4Poisson(current); // max cut-off still missing @@@
296 if(nDiscrete==1&&current<1.0001)
297 {
298 actualMult[i] = static_cast<G4int>(current);
299 if(current<1)
300 {
301 actualMult[i] = 0;
302 if(G4UniformRand()<current) actualMult[i] = 1;
303 }
304 }
305 nSecondaries += actualMult[i];
306 }
307 //G4cout << "nSecondaries " << nSecondaries << " anEnergy " << anEnergy/eV << G4endl;
308 for(i=0;i<nSecondaries;i++)
309 {
310 G4ReactionProduct * theOne = new G4ReactionProduct;
311 theOne->SetDefinition(G4Gamma::Gamma());
312 thePhotons->push_back(theOne);
313 }
314 G4int count=0;
315
316/*
317G4double totalCascadeEnergy = 0.;
318G4double lastCascadeEnergy = 0.;
319G4double eGamm = 0;
320G4int maxEnergyIndex = 0;
321*/
322 //Gcout << "nDiscrete " << nDiscrete << " nPartials " << nPartials << G4endl;
323//3456
324 if ( nDiscrete == 1 && nPartials == 1 )
325 {
326 if ( actualMult[ 0 ] > 0 )
327 {
328 if ( disType[0] == 1 ) // continuum
329 {
330
331/*
332 for(ii=0; ii< actualMult[0]; ii++)
333 {
334
335 G4double sum=0, run=0;
336 for(iii=0; iii<nPartials; iii++) sum+=probs[iii].GetY(anEnergy);
337 G4double random = G4UniformRand();
338 G4int theP = 0;
339 for(iii=0; iii<nPartials; iii++)
340 {
341 run+=probs[iii].GetY(anEnergy);
342 theP = iii;
343 if(random<run/sum) break;
344 }
345 if(theP==nPartials) theP=nPartials-1; // das sortiert J aus.
346 sum=0;
347 G4NeutronHPVector * temp;
348 temp = partials[theP]->GetY(anEnergy); //@@@ look at, seems fishy
349 // Looking for TotalCascdeEnergy or LastMaxEnergy
350 if (ii == 0)
351 {
352 maxEnergyIndex = temp->GetVectorLength()-1;
353 totalCascadeEnergy = temp->GetX(maxEnergyIndex);
354 lastCascadeEnergy = totalCascadeEnergy;
355 }
356 lastCascadeEnergy -= eGamm;
357 if (ii != actualMult[i]-1) eGamm = temp->SampleWithMax(lastCascadeEnergy);
358 else eGamm = lastCascadeEnergy;
359 thePhotons->operator[](count)->SetKineticEnergy(eGamm);
360 delete temp;
361
362 }
363*/
364 G4NeutronHPVector * temp;
365 temp = partials[ 0 ]->GetY(anEnergy); //@@@ look at, seems fishy
366 G4double maximumE = temp->GetX( temp->GetVectorLength()-1 ); // This is an assumption.
367
368 //G4cout << "start " << actualMult[ 0 ] << " maximumE " << maximumE/eV << G4endl;
369
370 std::vector< G4double > photons_e_best( actualMult[ 0 ] , 0.0 );
371 G4double best = DBL_MAX;
372 G4int maxTry = 1000;
373 for ( G4int j = 0 ; j < maxTry ; j++ )
374 {
375 std::vector< G4double > photons_e( actualMult[ 0 ] , 0.0 );
376 for ( std::vector< G4double >::iterator
377 it = photons_e.begin() ; it < photons_e.end() ; it++ )
378 {
379 *it = temp->Sample();
380 }
381 if ( std::accumulate( photons_e.begin() , photons_e.end() , 0.0 ) > maximumE )
382 {
383 if ( std::accumulate( photons_e.begin() , photons_e.end() , 0.0 ) < best )
384 photons_e_best = photons_e;
385 continue;
386 }
387 else
388 {
389 for ( std::vector< G4double >::iterator
390 it = photons_e.begin() ; it < photons_e.end() ; it++ )
391 {
392 thePhotons->operator[](count)->SetKineticEnergy( *it );
393 }
394 //G4cout << "OK " << actualMult[0] << " j " << j << " total photons E "
395 // << std::accumulate( photons_e.begin() , photons_e.end() , 0.0 )/eV << " ratio " << std::accumulate( photons_e.begin() , photons_e.end() , 0.0 ) / maximumE
396 // << G4endl;
397
398 break;
399 }
400 G4cout << "NeutronHPPhotonDist could not find fitted energy set for multiplicity of " << actualMult[0] << "." << G4endl;
401 G4cout << "NeutronHPPhotonDist will use the best set." << G4endl;
402 for ( std::vector< G4double >::iterator
403 it = photons_e_best.begin() ; it < photons_e_best.end() ; it++ )
404 {
405 thePhotons->operator[](count)->SetKineticEnergy( *it );
406 }
407 //G4cout << "Not Good " << actualMult[0] << " j " << j << " total photons E "
408 // << best/eV << " ratio " << best / maximumE
409 // << G4endl;
410 }
411 // TKDB
412 delete temp;
413 }
414 else // discrete
415 {
416 thePhotons->operator[](count)->SetKineticEnergy(energy[i]);
417 }
418 count++;
419 if(count > nSecondaries) throw G4HadronicException(__FILE__, __LINE__, "G4NeutronHPPhotonDist::GetPhotons inconsistancy");
420 }
421
422 }
423 else
424 {
425 for(i=0; i<nDiscrete; i++)
426 {
427 for(ii=0; ii< actualMult[i]; ii++)
428 {
429 if(disType[i]==1) // continuum
430 {
431 G4double sum=0, run=0;
432 for(iii=0; iii<nPartials; iii++) sum+=probs[iii].GetY(anEnergy);
433 G4double random = G4UniformRand();
434 G4int theP = 0;
435 for(iii=0; iii<nPartials; iii++)
436 {
437 run+=probs[iii].GetY(anEnergy);
438 theP = iii;
439 if(random<run/sum) break;
440 }
441 if(theP==nPartials) theP=nPartials-1; // das sortiert J aus.
442 sum=0;
443 G4NeutronHPVector * temp;
444 temp = partials[theP]->GetY(anEnergy); //@@@ look at, seems fishy
445 G4double eGamm = temp->Sample();
446 thePhotons->operator[](count)->SetKineticEnergy(eGamm);
447 delete temp;
448 }
449 else // discrete
450 {
451 thePhotons->operator[](count)->SetKineticEnergy(energy[i]);
452 }
453 count++;
454 if(count > nSecondaries) throw G4HadronicException(__FILE__, __LINE__, "G4NeutronHPPhotonDist::GetPhotons inconsistancy");
455 }
456 }
457 }
458 // now do the angular distributions...
459 if( isoFlag == 1)
460 {
461 for (i=0; i< nSecondaries; i++)
462 {
463 G4double costheta = 2.*G4UniformRand()-1;
464 G4double theta = std::acos(costheta);
465 G4double phi = twopi*G4UniformRand();
466 G4double sinth = std::sin(theta);
467 G4double en = thePhotons->operator[](i)->GetTotalEnergy();
468 G4ThreeVector temp(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
469 thePhotons->operator[](i)->SetMomentum( temp ) ;
470 // G4cout << "Isotropic distribution in PhotonDist"<<temp<<G4endl;
471 }
472 }
473 else
474 {
475 for(i=0; i<nSecondaries; i++)
476 {
477 G4double currentEnergy = thePhotons->operator[](i)->GetTotalEnergy();
478 for(ii=0; ii<nDiscrete2; ii++)
479 {
480 if (std::abs(currentEnergy-theGammas[ii])<0.1*keV) break;
481 }
482 if(ii==nDiscrete2) ii--; // fix for what seems an (file12 vs file 14) inconsistancy found in the ENDF 7N14 data. @@
483 if(ii<nIso)
484 {
485 // isotropic distribution
486 G4double theta = pi*G4UniformRand();
487 G4double phi = twopi*G4UniformRand();
488 G4double sinth = std::sin(theta);
489 G4double en = thePhotons->operator[](i)->GetTotalEnergy();
490 G4ThreeVector tempVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
491 thePhotons->operator[](i)->SetMomentum( tempVector ) ;
492 }
493 else if(tabulationType==1)
494 {
495 // legendre polynomials
496 G4int it(0);
497 for (iii=0; iii<nNeu[ii-nIso]; iii++) // find the neutron energy
498 {
499 it = iii;
500 if(theLegendre[ii-nIso][iii].GetEnergy()>anEnergy)
501 break;
502 }
503 G4NeutronHPLegendreStore aStore(2);
504 aStore.SetCoeff(1, &(theLegendre[ii-nIso][it]));
505 //aStore.SetCoeff(0, &(theLegendre[ii-nIso][it-1]));
506 //TKDB 110512
507 if ( it > 0 )
508 {
509 aStore.SetCoeff(0, &(theLegendre[ii-nIso][it-1]));
510 }
511 else
512 {
513 aStore.SetCoeff(0, &(theLegendre[ii-nIso][it]));
514 }
515 G4double cosTh = aStore.SampleMax(anEnergy);
516 G4double theta = std::acos(cosTh);
517 G4double phi = twopi*G4UniformRand();
518 G4double sinth = std::sin(theta);
519 G4double en = thePhotons->operator[](i)->GetTotalEnergy();
520 G4ThreeVector tempVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
521 thePhotons->operator[](i)->SetMomentum( tempVector ) ;
522 }
523 else
524 {
525 // tabulation of probabilities.
526 G4int it(0);
527 for (iii=0; iii<nNeu[ii-nIso]; iii++) // find the neutron energy
528 {
529 it = iii;
530 if(theAngular[ii-nIso][iii].GetEnergy()>anEnergy)
531 break;
532 }
533 G4double costh = theAngular[ii-nIso][it].GetCosTh(); // no interpolation yet @@
534 G4double theta = std::acos(costh);
535 G4double phi = twopi*G4UniformRand();
536 G4double sinth = std::sin(theta);
537 G4double en = thePhotons->operator[](i)->GetTotalEnergy();
538 G4ThreeVector tmpVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*costh );
539 thePhotons->operator[](i)->SetMomentum( tmpVector ) ;
540 }
541 }
542 }
543 }
544 else if(repFlag == 2)
545 {
546 G4double * running = new G4double[nGammaEnergies];
547 running[0]=theTransitionProbabilities[0];
548 //G4int i; //declaration at 284th
549 for(i=1; i<nGammaEnergies; i++)
550 {
551 running[i]=running[i-1]+theTransitionProbabilities[i];
552 }
553 G4double random = G4UniformRand();
554 G4int it=0;
555 for(i=0; i<nGammaEnergies; i++)
556 {
557 it = i;
558 if(random < running[i]/running[nGammaEnergies-1]) break;
559 }
560 delete [] running;
561 G4double totalEnergy = theBaseEnergy - theLevelEnergies[it];
562 G4ReactionProduct * theOne = new G4ReactionProduct;
563 theOne->SetDefinition(G4Gamma::Gamma());
564 random = G4UniformRand();
565 if(theInternalConversionFlag==2 && random>thePhotonTransitionFraction[it])
566 {
567 theOne->SetDefinition(G4Electron::Electron());
568 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
569 //But never enter at least with G4NDL3.13
570 totalEnergy += G4Electron::Electron()->GetPDGMass(); //proposed correction: add this line for electron
571 }
572 theOne->SetTotalEnergy(totalEnergy);
573 if( isoFlag == 1 )
574 {
575 G4double costheta = 2.*G4UniformRand()-1;
576 G4double theta = std::acos(costheta);
577 G4double phi = twopi*G4UniformRand();
578 G4double sinth = std::sin(theta);
579 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
580 //G4double en = theOne->GetTotalEnergy();
581 G4double en = theOne->GetTotalMomentum();
582 //But never cause real effect at least with G4NDL3.13 TK
583 G4ThreeVector temp(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
584 theOne->SetMomentum( temp ) ;
585 }
586 else
587 {
588 G4double currentEnergy = theOne->GetTotalEnergy();
589 for(ii=0; ii<nDiscrete2; ii++)
590 {
591 if (std::abs(currentEnergy-theGammas[ii])<0.1*keV) break;
592 }
593 if(ii==nDiscrete2) ii--; // fix for what seems an (file12 vs file 14) inconsistancy found in the ENDF 7N14 data. @@
594 if(ii<nIso)
595 {
596 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
597 // isotropic distribution
598 //G4double theta = pi*G4UniformRand();
599 G4double theta = std::acos(2.*G4UniformRand()-1.);
600 //But this is alos never cause real effect at least with G4NDL3.13 TK not repFlag == 2 AND isoFlag != 1
601 G4double phi = twopi*G4UniformRand();
602 G4double sinth = std::sin(theta);
603 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
604 //G4double en = theOne->GetTotalEnergy();
605 G4double en = theOne->GetTotalMomentum();
606 //But never cause real effect at least with G4NDL3.13 TK
607 G4ThreeVector tempVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
608 theOne->SetMomentum( tempVector ) ;
609 }
610 else if(tabulationType==1)
611 {
612 // legendre polynomials
613 G4int itt(0);
614 for (iii=0; iii<nNeu[ii-nIso]; iii++) // find the neutron energy
615 {
616 itt = iii;
617 if(theLegendre[ii-nIso][iii].GetEnergy()>anEnergy)
618 break;
619 }
620 G4NeutronHPLegendreStore aStore(2);
621 aStore.SetCoeff(1, &(theLegendre[ii-nIso][itt]));
622 //aStore.SetCoeff(0, &(theLegendre[ii-nIso][it-1]));
623 //TKDB 110512
624 if ( itt > 0 )
625 {
626 aStore.SetCoeff(0, &(theLegendre[ii-nIso][itt-1]));
627 }
628 else
629 {
630 aStore.SetCoeff(0, &(theLegendre[ii-nIso][itt]));
631 }
632 G4double cosTh = aStore.SampleMax(anEnergy);
633 G4double theta = std::acos(cosTh);
634 G4double phi = twopi*G4UniformRand();
635 G4double sinth = std::sin(theta);
636 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
637 //G4double en = theOne->GetTotalEnergy();
638 G4double en = theOne->GetTotalMomentum();
639 //But never cause real effect at least with G4NDL3.13 TK
640 G4ThreeVector tempVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*std::cos(theta) );
641 theOne->SetMomentum( tempVector ) ;
642 }
643 else
644 {
645 // tabulation of probabilities.
646 G4int itt(0);
647 for (iii=0; iii<nNeu[ii-nIso]; iii++) // find the neutron energy
648 {
649 itt = iii;
650 if(theAngular[ii-nIso][iii].GetEnergy()>anEnergy)
651 break;
652 }
653 G4double costh = theAngular[ii-nIso][itt].GetCosTh(); // no interpolation yet @@
654 G4double theta = std::acos(costh);
655 G4double phi = twopi*G4UniformRand();
656 G4double sinth = std::sin(theta);
657 //Bug reported Chao Zhang ([email protected]), Dongming Mei([email protected]) Feb. 25, 2009
658 //G4double en = theOne->GetTotalEnergy();
659 G4double en = theOne->GetTotalMomentum();
660 //But never cause real effect at least with G4NDL3.13 TK
661 G4ThreeVector tmpVector(en*sinth*std::cos(phi), en*sinth*std::sin(phi), en*costh );
662 theOne->SetMomentum( tmpVector ) ;
663 }
664 }
665 thePhotons->push_back(theOne);
666 }
667 else if( repFlag==0 )
668 {
669
670// TK add
671 if ( thePartialXsec == 0 )
672 {
673 //G4cout << "repFlag is 0, but no PartialXsec data" << G4endl;
674 //G4cout << "This is not support yet." << G4endl;
675 return thePhotons;
676 }
677
678// Partial Case
679
680 G4ReactionProduct * theOne = new G4ReactionProduct;
681 theOne->SetDefinition( G4Gamma::Gamma() );
682 thePhotons->push_back( theOne );
683
684// Energy
685
686 //G4cout << "Partial Case nDiscrete " << nDiscrete << G4endl;
687 G4double sum = 0.0;
688 std::vector < G4double > dif( nDiscrete , 0.0 );
689 for ( G4int j = 0 ; j < nDiscrete ; j++ )
690 {
691 G4double x = thePartialXsec[ j ].GetXsec( anEnergy ); // x in barn
692 if ( x > 0 )
693 {
694 sum += x;
695 }
696 dif [ j ] = sum;
697 //G4cout << "j " << j << ", x " << x << ", dif " << dif [ j ] << G4endl;
698 }
699
700 G4double rand = G4UniformRand();
701
702 G4int iphoton = 0;
703 for ( G4int j = 0 ; j < nDiscrete ; j++ )
704 {
705 G4double y = rand*sum;
706 if ( dif [ j ] > y )
707 {
708 iphoton = j;
709 break;
710 }
711 }
712 //G4cout << "iphoton " << iphoton << G4endl;
713 //G4cout << "photon energy " << theGammas[ iphoton ] /eV << G4endl;
714
715// Angle
716 G4double cosTheta = 0.0; // mu
717
718 if ( isoFlag == 1 )
719 {
720
721// Isotropic Case
722
723 cosTheta = 2.*G4UniformRand()-1;
724
725 }
726 else
727 {
728
729 if ( iphoton < nIso )
730 {
731
732// still Isotropic
733
734 cosTheta = 2.*G4UniformRand()-1;
735
736 }
737 else
738 {
739
740 //G4cout << "Not Isotropic and isoFlag " << isoFlag << G4endl;
741 //G4cout << "tabulationType " << tabulationType << G4endl;
742 //G4cout << "nDiscrete2 " << nDiscrete2 << G4endl;
743 //G4cout << "nIso " << nIso << G4endl;
744 //G4cout << "size of nNeu " << nDiscrete2-nIso << G4endl;
745 //G4cout << "nNeu[iphoton-nIso] " << nNeu[iphoton-nIso] << G4endl;
746
747 if ( tabulationType == 1 )
748 {
749// legendre polynomials
750
751 G4int iangle = 0;
752 for ( G4int j = 0 ; j < nNeu [ iphoton - nIso ] ; j++ )
753 {
754 iangle = j;
755 if ( theLegendre[ iphoton - nIso ][ j ].GetEnergy() > anEnergy ) break;
756 }
757
758 G4NeutronHPLegendreStore aStore( 2 );
759 aStore.SetCoeff( 1 , &( theLegendre[ iphoton - nIso ][ iangle ] ) );
760 aStore.SetCoeff( 0 , &( theLegendre[ iphoton - nIso ][ iangle - 1 ] ) );
761
762 cosTheta = aStore.SampleMax( anEnergy );
763
764 }
765 else if ( tabulationType == 2 )
766 {
767
768// tabulation of probabilities.
769
770 G4int iangle = 0;
771 for ( G4int j = 0 ; j < nNeu [ iphoton - nIso ] ; j++ )
772 {
773 iangle = j;
774 if ( theAngular[ iphoton - nIso ][ j ].GetEnergy() > anEnergy ) break;
775 }
776
777 cosTheta = theAngular[iphoton-nIso][ iangle ].GetCosTh(); // no interpolation yet @@
778
779 }
780 }
781 }
782
783// Set
784 G4double phi = twopi*G4UniformRand();
785 G4double theta = std::acos( cosTheta );
786 G4double sinTheta = std::sin( theta );
787
788 G4double photonE = theGammas[ iphoton ];
789 G4ThreeVector direction ( sinTheta*std::cos( phi ) , sinTheta * std::sin( phi ) , cosTheta );
790 G4ThreeVector photonP = photonE * direction;
791 thePhotons->operator[]( 0 )->SetMomentum( photonP ) ;
792
793 }
794 else
795 {
796 delete thePhotons;
797 thePhotons = 0; // no gamma data available; some work needed @@@@@@@
798 }
799 return thePhotons;
800}
G4Poisson
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:50
G4ReactionProductVector
std::vector< G4ReactionProduct * > G4ReactionProductVector
Definition: G4ReactionProductVector.hh:43
G4double
double G4double
Definition: G4Types.hh:64
G4endl
#define G4endl
Definition: G4ios.hh:52
G4cout
G4DLLIMPORT std::ostream G4cout
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:53
G4Electron::Electron
static G4Electron * Electron()
Definition: G4Electron.cc:94
G4Gamma::Gamma
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4NeutronHPAngularP::GetCosTh
G4double GetCosTh(G4int l)
Definition: G4NeutronHPAngularP.hh:80
G4NeutronHPPartial::GetY
G4double GetY(G4int i, G4int j)
Definition: G4NeutronHPPartial.hh:116
G4NeutronHPVector::GetVectorLength
G4int GetVectorLength() const
Definition: G4NeutronHPVector.hh:193
G4NeutronHPVector::GetX
G4double GetX(G4int i) const
Definition: G4NeutronHPVector.hh:121
G4NeutronHPVector::GetXsec
G4double GetXsec(G4int i)
Definition: G4NeutronHPVector.hh:120
G4NeutronHPVector::GetY
G4double GetY(G4double x)
Definition: G4NeutronHPVector.hh:192
G4ReactionProduct::SetMomentum
void SetMomentum(const G4double x, const G4double y, const G4double z)
Definition: G4ReactionProduct.cc:166
G4ReactionProduct::SetTotalEnergy
void SetTotalEnergy(const G4double en)
Definition: G4ReactionProduct.hh:138
G4ReactionProduct::GetTotalMomentum
G4double GetTotalMomentum() const
Definition: G4ReactionProduct.hh:123
G4ReactionProduct::GetTotalEnergy
G4double GetTotalEnergy() const
Definition: G4ReactionProduct.hh:126
G4ReactionProduct::SetDefinition
void SetDefinition(G4ParticleDefinition *aParticleDefinition)
Definition: G4ReactionProduct.cc:155
G4INCL::Math::pi
const G4double pi
Definition: G4INCLGlobals.hh:70
DBL_MAX
#define DBL_MAX
Definition: templates.hh:83

Referenced by G4NeutronHPCaptureFS::ApplyYourself(), G4NeutronHPInelasticBaseFS::BaseApply(), G4NeutronHPInelasticCompFS::CompositeApply(), and G4NeutronHPFSFissionFS::GetPhotons().

◆ GetTargetMass()

G4double G4NeutronHPPhotonDist::GetTargetMass ( )
inline

Definition at line 145 of file G4NeutronHPPhotonDist.hh.

145{return targetMass;}

Referenced by G4NeutronHPCaptureFS::Init().

◆ InitAngular()

void G4NeutronHPPhotonDist::InitAngular ( std::ifstream &  aDataFile)

Definition at line 121 of file G4NeutronHPPhotonDist.cc.

122{
123
124 G4int i, ii;
125 //angular distributions
126 aDataFile >> isoFlag;
127 if (isoFlag != 1)
128 {
129if ( repFlag == 2 ) G4cout << "G4NeutronHPPhotonDist: repFlag == 2 && isoFlag != 1 is unexpected! If you use G4ND3.x, then please report to Geant4 Hyper News. Thanks." << G4endl;
130 aDataFile >> tabulationType >> nDiscrete2 >> nIso;
131//080731
132 if ( theGammas != NULL && nDiscrete2 != nDiscrete ) G4cout << "080731c G4NeutronHPPhotonDist nDiscrete2 != nDiscrete, It looks like something wrong in your NDL files. Please update the latest. If you still have this messages after the update, then please report to Geant4 Hyper News." << G4endl;
133
134 // The order of cross section (InitPartials) and distribution (InitAngular here) data are different, we have to re-coordinate consistent data order.
135 std::vector < G4double > vct_gammas_par;
136 std::vector < G4double > vct_shells_par;
137 std::vector < G4int > vct_primary_par;
138 std::vector < G4int > vct_distype_par;
139 std::vector < G4NeutronHPVector* > vct_pXS_par;
140 if ( theGammas != NULL )
141 {
142 //copy the cross section data
143 for ( i = 0 ; i < nDiscrete ; i++ )
144 {
145 vct_gammas_par.push_back( theGammas[ i ] );
146 vct_shells_par.push_back( theShells[ i ] );
147 vct_primary_par.push_back( isPrimary[ i ] );
148 vct_distype_par.push_back( disType[ i ] );
149 G4NeutronHPVector* hpv = new G4NeutronHPVector;
150 *hpv = thePartialXsec[ i ];
151 vct_pXS_par.push_back( hpv );
152 }
153 }
154 if ( theGammas == NULL ) theGammas = new G4double[nDiscrete2];
155 if ( theShells == NULL ) theShells = new G4double[nDiscrete2];
156//080731
157
158 for (i=0; i< nIso; i++) // isotropic photons
159 {
160 aDataFile >> theGammas[i] >> theShells[i];
161 theGammas[i]*=eV;
162 theShells[i]*=eV;
163 }
164 nNeu = new G4int [nDiscrete2-nIso];
165 if(tabulationType==1)theLegendre=new G4NeutronHPLegendreTable *[nDiscrete2-nIso];
166 if(tabulationType==2)theAngular =new G4NeutronHPAngularP *[nDiscrete2-nIso];
167 for(i=nIso; i< nDiscrete2; i++)
168 {
169 if(tabulationType==1)
170 {
171 aDataFile >> theGammas[i] >> theShells[i] >> nNeu[i-nIso];
172 theGammas[i]*=eV;
173 theShells[i]*=eV;
174 theLegendre[i-nIso]=new G4NeutronHPLegendreTable[nNeu[i-nIso]];
175 theLegendreManager.Init(aDataFile);
176 for (ii=0; ii<nNeu[i-nIso]; ii++)
177 {
178 theLegendre[i-nIso][ii].Init(aDataFile);
179 }
180 }
181 else if(tabulationType==2)
182 {
183 aDataFile >> theGammas[i] >> theShells[i] >> nNeu[i-nIso];
184 theGammas[i]*=eV;
185 theShells[i]*=eV;
186 theAngular[i-nIso]=new G4NeutronHPAngularP[nNeu[i-nIso]];
187 for (ii=0; ii<nNeu[i-nIso]; ii++)
188 {
189 theAngular[i-nIso][ii].Init(aDataFile);
190 }
191 }
192 else
193 {
194 G4cout << "tabulation type: tabulationType"<<G4endl;
195 throw G4HadronicException(__FILE__, __LINE__, "cannot deal with this tabulation type for angular distributions.");
196 }
197 }
198//080731
199 if ( vct_gammas_par.size() > 0 )
200 {
201 //Reordering cross section data to corrsponding distribution data
202 for ( i = 0 ; i < nDiscrete ; i++ )
203 {
204 for ( G4int j = 0 ; j < nDiscrete ; j++ )
205 {
206 // Checking gamma and shell to identification
207 if ( theGammas[ i ] == vct_gammas_par [ j ] && theShells [ i ] == vct_shells_par[ j ] )
208 {
209 isPrimary [ i ] = vct_primary_par [ j ];
210 disType [ i ] = vct_distype_par [ j ];
211 thePartialXsec[ i ] = ( *( vct_pXS_par[ j ] ) );
212 }
213 }
214 }
215 //Garbage collection
216 for ( std::vector < G4NeutronHPVector* >::iterator
217 it = vct_pXS_par.begin() ; it != vct_pXS_par.end() ; it++ )
218 {
219 delete *it;
220 }
221 }
222//080731
223 }
224}
G4InterpolationManager::Init
void Init(G4int aScheme, G4int aRange)
Definition: G4InterpolationManager.hh:85
G4NeutronHPAngularP::Init
void Init(std::ifstream &aDataFile)
Definition: G4NeutronHPAngularP.hh:52
G4NeutronHPLegendreTable::Init
void Init(std::ifstream &aDataFile)
Definition: G4NeutronHPLegendreTable.hh:67

Referenced by G4FissionLibrary::Init(), G4NeutronHPCaptureFS::Init(), G4NeutronHPFSFissionFS::Init(), G4NeutronHPInelasticCompFS::Init(), and G4NeutronHPInelasticBaseFS::Init().

◆ InitEnergies()

void G4NeutronHPPhotonDist::InitEnergies ( std::ifstream &  aDataFile)

Definition at line 227 of file G4NeutronHPPhotonDist.cc.

228{
229 G4int i, energyDistributionsNeeded = 0;
230 for (i=0; i<nDiscrete; i++)
231 {
232 if( disType[i]==1) energyDistributionsNeeded =1;
233 }
234 if(!energyDistributionsNeeded) return;
235 aDataFile >> nPartials;
236 distribution = new G4int[nPartials];
237 probs = new G4NeutronHPVector[nPartials];
238 partials = new G4NeutronHPPartial * [nPartials];
239 G4int nen;
240 G4int dummy;
241 for (i=0; i<nPartials; i++)
242 {
243 aDataFile >> dummy;
244 probs[i].Init(aDataFile, eV);
245 aDataFile >> nen;
246 partials[i] = new G4NeutronHPPartial(nen);
247 partials[i]->InitInterpolation(aDataFile);
248 partials[i]->Init(aDataFile);
249 }
250}
G4NeutronHPPartial::Init
void Init(std::ifstream &aDataFile)
Definition: G4NeutronHPPartial.hh:71
G4NeutronHPPartial::InitInterpolation
void InitInterpolation(G4int i, std::ifstream &aDataFile)
Definition: G4NeutronHPPartial.hh:61
G4NeutronHPVector::Init
void Init(std::ifstream &aDataFile, G4int total, G4double ux=1., G4double uy=1.)
Definition: G4NeutronHPVector.hh:215

Referenced by G4FissionLibrary::Init(), G4NeutronHPCaptureFS::Init(), G4NeutronHPFSFissionFS::Init(), G4NeutronHPInelasticCompFS::Init(), and G4NeutronHPInelasticBaseFS::Init().

◆ InitMean()

G4bool G4NeutronHPPhotonDist::InitMean ( std::ifstream &  aDataFile)

Definition at line 54 of file G4NeutronHPPhotonDist.cc.

55{
56
57 G4bool result = true;
58 if(aDataFile >> repFlag)
59 {
60
61 aDataFile >> targetMass;
62 if(repFlag==1)
63 {
64 // multiplicities
65 aDataFile >> nDiscrete;
66 disType = new G4int[nDiscrete];
67 energy = new G4double[nDiscrete];
68 actualMult = new G4int[nDiscrete];
69 theYield = new G4NeutronHPVector[nDiscrete];
70 for (G4int i=0; i<nDiscrete; i++)
71 {
72 aDataFile >> disType[i]>>energy[i];
73 energy[i]*=eV;
74 theYield[i].Init(aDataFile, eV);
75 }
76 }
77 else if(repFlag == 2)
78 {
79 aDataFile >> theInternalConversionFlag;
80 aDataFile >> theBaseEnergy;
81 theBaseEnergy*=eV;
82 aDataFile >> theInternalConversionFlag;
83 // theInternalConversionFlag == 1 No IC, theInternalConversionFlag == 2 with IC
84 aDataFile >> nGammaEnergies;
85 theLevelEnergies = new G4double[nGammaEnergies];
86 theTransitionProbabilities = new G4double[nGammaEnergies];
87 if(theInternalConversionFlag == 2) thePhotonTransitionFraction = new G4double[nGammaEnergies];
88 for(G4int ii=0; ii<nGammaEnergies; ii++)
89 {
90 if(theInternalConversionFlag == 1)
91 {
92 aDataFile >> theLevelEnergies[ii] >> theTransitionProbabilities[ii];
93 theLevelEnergies[ii]*=eV;
94 }
95 else if(theInternalConversionFlag == 2)
96 {
97 aDataFile >> theLevelEnergies[ii] >> theTransitionProbabilities[ii] >> thePhotonTransitionFraction[ii];
98 theLevelEnergies[ii]*=eV;
99 }
100 else
101 {
102 throw G4HadronicException(__FILE__, __LINE__, "G4NeutronHPPhotonDist: Unknown conversion flag");
103 }
104 }
105 // Note, that this is equivalent to using the 'Gamma' classes.
106 // throw G4HadronicException(__FILE__, __LINE__, "G4NeutronHPPhotonDist: Transition probability array not sampled for the moment.");
107 }
108 else
109 {
110 G4cout << "Data representation in G4NeutronHPPhotonDist: "<<repFlag<<G4endl;
111 throw G4HadronicException(__FILE__, __LINE__, "G4NeutronHPPhotonDist: This data representation is not implemented.");
112 }
113 }
114 else
115 {
116 result = false;
117 }
118 return result;
119}
G4bool
bool G4bool
Definition: G4Types.hh:67

Referenced by G4FissionLibrary::Init(), G4NeutronHPCaptureFS::Init(), G4NeutronHPFSFissionFS::Init(), G4NeutronHPInelasticCompFS::Init(), and G4NeutronHPInelasticBaseFS::Init().

◆ InitPartials()

void G4NeutronHPPhotonDist::InitPartials ( std::ifstream &  aDataFile)

Definition at line 252 of file G4NeutronHPPhotonDist.cc.

253{
254
255 //G4cout << "G4NeutronHPPhotonDist::InitPartials " << G4endl;
256 aDataFile >> nDiscrete >> targetMass;
257 if(nDiscrete != 1)
258 {
259 theTotalXsec.Init(aDataFile, eV);
260 }
261 G4int i;
262 theGammas = new G4double[nDiscrete];
263 theShells = new G4double[nDiscrete];
264 isPrimary = new G4int[nDiscrete];
265 disType = new G4int[nDiscrete];
266 thePartialXsec = new G4NeutronHPVector[nDiscrete];
267 for(i=0; i<nDiscrete; i++)
268 {
269 aDataFile>>theGammas[i]>>theShells[i]>>isPrimary[i]>>disType[i];
270 theGammas[i]*=eV;
271 theShells[i]*=eV;
272 thePartialXsec[i].Init(aDataFile, eV);
273 }
274
275 //G4cout << "G4NeutronHPPhotonDist::InitPartials Test " << G4endl;
276 //G4cout << "G4NeutronHPPhotonDist::InitPartials nDiscrete " << nDiscrete << G4endl;
277 //G4NeutronHPVector* aHP = new G4NeutronHPVector;
278 //aHP->Check(1);
279}

Referenced by G4NeutronHPInelasticCompFS::Init(), and G4NeutronHPInelasticBaseFS::Init().

◆ NeedsCascade()

G4bool G4NeutronHPPhotonDist::NeedsCascade ( )
inline

Definition at line 147 of file G4NeutronHPPhotonDist.hh.

147{return repFlag==2;}
The documentation for this class was generated from the following files: