/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Public Attributes | Private Attributes
DsG4NeutronHPCapture Class Reference

#include <DsG4NeutronHPCapture.h>

Collaboration diagram for DsG4NeutronHPCapture:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 DsG4NeutronHPCapture ()
virtual ~DsG4NeutronHPCapture ()
G4HadFinalState * ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus)
void passNeutronCaptureInfoTool (INeutronCaptureInfo *)

Public Attributes

G4bool DebugMe

Private Attributes

G4double * xSec
G4NeutronHPChannel * theCapture
G4String dirName
G4int numEle
G4int it
G4HadFinalState theResult
G4HadFinalState * result
INeutronCaptureInfom_capinfo
G4DhNeutronCapture nCapture

Detailed Description

Definition at line 24 of file DsG4NeutronHPCapture.h.


Constructor & Destructor Documentation

DsG4NeutronHPCapture::DsG4NeutronHPCapture ( )

Definition at line 32 of file DsG4NeutronHPCapture.cc.

{
  DebugMe = false;
  if (DebugMe) G4cout  << "Begin DsG4NeutronCapture constructor " << G4endl;

  SetMinEnergy( 0.0 );
  SetMaxEnergy( 20.*MeV );
  
  char* envstr = getenv("G4NEUTRONHPDATA");

  if(!envstr) {
    throw G4HadronicException(__FILE__, __LINE__, 
                              "Please setenv G4NEUTRONHPDATA "
                              "to point to the neutron cross-section files.");
  }
  G4String tString = "/Capture/";
  dirName = envstr;
  dirName = dirName + tString;
  numEle = G4Element::GetNumberOfElements();

  G4cout << "Initializing DsG4NeutronCapture for " << numEle << " neutron HP channels" << G4endl;
  
  theCapture = new G4NeutronHPChannel[numEle];

  G4NeutronHPCaptureFS * theFS = new G4NeutronHPCaptureFS; 
  for (G4int i=0; i<numEle; i++) {
    
    if (DebugMe) G4cout << "initializing theCapture "<<i<<" "<< numEle
            << " name " << (*(G4Element::GetElementTable()))[i]->GetName()
            << G4endl; 
    // DEBUG:
    
    G4cout << "initializing theCapture "<<i<<" "<< numEle
         << ", name " << (*(G4Element::GetElementTable()))[i]->GetName() 
         << ", symbol " << (*(G4Element::GetElementTable()))[i]->GetSymbol() 
         << ", z   " << (*(G4Element::GetElementTable()))[i]->GetZ() 
         << ", n   " << (*(G4Element::GetElementTable()))[i]->GetN() 
         << ", a   " << (*(G4Element::GetElementTable()))[i]->GetA() 
         << ", nat " << (*(G4Element::GetElementTable()))[i]->GetNaturalAbandancesFlag() 
         << ", Niso   " << (*(G4Element::GetElementTable()))[i]->GetNumberOfIsotopes() 
         << G4endl;

    // Cast Z to nearest integer
    int elemZ = (int)((*(G4Element::GetElementTable()))[i]->GetZ() + 0.5);

    // To use new NNDC Gd Capture Gammas, change this to true
    // To use old Gd Capture Gammas from nuclear data tables (gdcap_spec.root),
    // but generated with an unbiased sampling method (Doc-5750), also need to change this to true
    // (see 'useOldNdtGdSpectra' in DsG4NNDCCaptureGammas.cc for more info)
    bool useNNDC_GdCapture = true; 

    if(elemZ == 64 && !useNNDC_GdCapture) { // Gd.
      // for Gd, DsG4GdNeutronHPCaptureFS is invoked.
      DsG4GdNeutronHPCaptureFS * theGdFS = new DsG4GdNeutronHPCaptureFS;
      theCapture[i].Init((*(G4Element::GetElementTable()))[i], dirName);
      if (DebugMe) G4cout << (*(G4Element::GetElementTable()))[i]->GetName() << G4endl;
      theCapture[i].Register(theGdFS);
      delete theGdFS;
    }
    else if(elemZ == 6 /*C*/
            || elemZ == 7 /*N*/
            || elemZ == 8 /*O*/
            || elemZ == 14 /*Si*/
            || elemZ == 15 /*P*/
            || elemZ == 16 /*S*/
            || elemZ == 24 /*Cr*/
            || elemZ == 25 /*Mn*/
            || elemZ == 26 /*Fe*/
            || elemZ == 28 /*Ni*/
            || (elemZ == 64 && useNNDC_GdCapture) /*Gd*/){
      // For these elements, we use hand-generated tables of correlated gammas.
      // These tables ensure conservation of energy
      theCapture[i].Init((*(G4Element::GetElementTable()))[i], dirName);
      DsG4NNDCNeutronHPCaptureFS * theNNDCFS = new DsG4NNDCNeutronHPCaptureFS();
      theCapture[i].Register(theNNDCFS);
      delete theNNDCFS;
    }
    else { 
      theCapture[i].Init((*(G4Element::GetElementTable()))[i], dirName);
      theCapture[i].Register(theFS);
    }
  }
  delete theFS;
}
DsG4NeutronHPCapture::~DsG4NeutronHPCapture ( ) [virtual]

Definition at line 117 of file DsG4NeutronHPCapture.cc.

{
  delete [] theCapture;
}

Member Function Documentation

G4HadFinalState * DsG4NeutronHPCapture::ApplyYourself ( const G4HadProjectile &  aTrack,
G4Nucleus &  aTargetNucleus 
)

Definition at line 124 of file DsG4NeutronHPCapture.cc.

{
  result = new G4HadFinalState(); 

  // Initialize
  G4int gammaNum = 0; 
  G4double capGammaE[20] = {0}; 
  G4double capTime = 0; 
  G4double capGammaEsum = 0;
  
  //if(getenv("NeutronHPCapture")) 
  if ( DebugMe ) {G4cout <<" ### DsG4NeutronHPCapture called"<< G4endl;} 

  // get cross-sections for elements in current material
  const G4Material * theMaterial = aTrack.GetMaterial();
  G4int n = theMaterial->GetNumberOfElements();
  G4int index = theMaterial->GetElement(0)->GetIndex();
  if(n!=1) {
    xSec = new G4double[n];
    G4double sum=0;
    G4int i;
    const G4double * NumAtomsPerVolume = theMaterial->GetVecNbOfAtomsPerVolume();
    G4double rWeight;    
    G4NeutronHPThermalBoost aThermalE;
    for (i=0; i<n; i++) {
      index = theMaterial->GetElement(i)->GetIndex();
      rWeight = NumAtomsPerVolume[i];
      xSec[i] = theCapture[index].GetXsec(
                                          aThermalE.GetThermalEnergy(aTrack,
                                                                     theMaterial->GetElement(i),
                                                                     theMaterial->GetTemperature()));

      xSec[i] *= rWeight;
      sum+=xSec[i];
    }
    
    // determine the target nucleus
    G4double random = G4UniformRand();
    G4double running = 0;
    for (i=0; i<n; i++) {
      running += xSec[i];
      index = theMaterial->GetElement(i)->GetIndex();
      if(random<=running/sum) break;
    }
    if(i==n) i=std::max(0, n-1);
    delete [] xSec;
  }

  capTime = aTrack.GetGlobalTime()/(1000*ns);
  
  G4String targetname = (*(G4Element::GetElementTable()))[index]->GetName();
  G4double targetZ = (*(G4Element::GetElementTable()))[index]->GetZ();
  G4double targetA = (*(G4Element::GetElementTable()))[index]->GetN(); 
  
  result = new G4HadFinalState();

  // Allow up to 100 tries to get physically meaningful capture gamma energy and number of gammas
  G4int tries = 100;  
  G4bool done = false;  
  while (!done) { 

    gammaNum = 0; 
    capGammaEsum = 0;

    result = theCapture[index].ApplyYourself(aTrack);
      
    G4int num = result->GetNumberOfSecondaries();

    if (DebugMe) G4cout << "DDEBUG: number of secondaries: " << num << G4endl;
      
    G4String secname;
    G4double seckine;
    for(int ii=0;ii<num;ii++) {
      secname=(result->GetSecondary(ii))->GetParticle()->GetDefinition()->GetParticleName();
      seckine=(result->GetSecondary(ii))->GetParticle()->GetKineticEnergy()/MeV;  

      if (DebugMe) G4cout << "   DDEBUG   name: " << secname << G4endl; 
      // Attention: the recoiling target is one of the secondaries.
      if(secname == "gamma") {
        gammaNum++;
        capGammaE[gammaNum-1] = seckine;
        nCapture.PushCapGammaE(seckine); // Push info into n-cap tool
      }
      if(secname == "gamma") {
        capGammaEsum += seckine;
      }    
    }//end of loop over secondaries
      
    // djaffe: decide if gammas are physically reasonable
      //         Captures on H and Gd give valid results.
      //         Enforce # of gammas and total energy for captures on carbon.
      //         For captures on other nuclei, just avoid zero or very small gamma energy.
      //         This needs to be fixed in the future.
      int iz = (int)(targetZ + 0.5); // convert to integer

      switch (iz) {
      case 1: // hydrogen
        done = true ; // nH always OK
        break;
      case 64: // Gd
        done = true ; // nGd always OK
        break;
        //case 6: // carbon
        //      done = ((gammaNum == 1) || (gammaNum == 2)) && std::abs(capGammaEsum-4.946*MeV)<0.1*MeV ;
        //      break;
      default: // everything else
        if ( capGammaEsum > 0.01*MeV ) {
          G4bool zero = false;
          for (int ii=0; ii<gammaNum; ii++) { if (capGammaE[ii] < 0.001*MeV) zero = true; }
          done = !zero;
        }
        break;
      }
      --tries; // don't try forever
      if (tries == 0)   {
        done = true; 
        G4cout << " DsG4NeutronHPCapture: GIVING UP. Could not achieve acceptable final state gammas  " << G4endl;
        G4cout << " DsG4NeutronHPCapture: Z,A " <<  targetZ <<","<<targetA
                 <<" "<<targetname << " N(gamma)="
                 << gammaNum << " E(gamma)= " ;                             
        for (int ii=0;ii<gammaNum;ii++){ G4cout << capGammaE[ii] <<", " ;} 
        G4cout << G4endl; 
      }
    } // !done
    
    
    /* recording the capture information --- Wei Wang, Aug 14, 2008 */
    nCapture.SetCapTargetZ(targetZ);
    nCapture.SetCapTargetA(targetA);
    nCapture.SetCapTime(capTime);
    nCapture.SetCapGammaN(gammaNum);
    
    m_capinfo->addCapture(nCapture);
    /* capture information recorded */

    return result; 
}
void DsG4NeutronHPCapture::passNeutronCaptureInfoTool ( INeutronCaptureInfo p_capinfo)

Definition at line 263 of file DsG4NeutronHPCapture.cc.

{
  m_capinfo = p_capinfo;
}

Member Data Documentation

Definition at line 34 of file DsG4NeutronHPCapture.h.

G4double* DsG4NeutronHPCapture::xSec [private]

Definition at line 40 of file DsG4NeutronHPCapture.h.

G4NeutronHPChannel* DsG4NeutronHPCapture::theCapture [private]

Definition at line 41 of file DsG4NeutronHPCapture.h.

G4String DsG4NeutronHPCapture::dirName [private]

Definition at line 42 of file DsG4NeutronHPCapture.h.

Definition at line 43 of file DsG4NeutronHPCapture.h.

G4int DsG4NeutronHPCapture::it [private]

Definition at line 44 of file DsG4NeutronHPCapture.h.

G4HadFinalState DsG4NeutronHPCapture::theResult [private]

Definition at line 46 of file DsG4NeutronHPCapture.h.

G4HadFinalState* DsG4NeutronHPCapture::result [private]

Definition at line 47 of file DsG4NeutronHPCapture.h.

Definition at line 52 of file DsG4NeutronHPCapture.h.

Definition at line 53 of file DsG4NeutronHPCapture.h.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:17:58 for DetSim by doxygen 1.7.4