G4NeutronHPCapture.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // neutron_hp -- source file
00027 // J.P. Wellisch, Nov-1996
00028 // A prototype of the low energy neutron transport model.
00029 //
00030 // 070523 bug fix for G4FPE_DEBUG on by A. Howard ( and T. Koi)
00031 //
00032 #include "G4NeutronHPCapture.hh"
00033 #include "G4SystemOfUnits.hh"
00034 #include "G4NeutronHPCaptureFS.hh"
00035 #include "G4NeutronHPDeExGammas.hh"
00036 #include "G4ParticleTable.hh"
00037 #include "G4IonTable.hh"
00038 
00039 #include "G4NeutronHPManager.hh"
00040 
00041   G4NeutronHPCapture::G4NeutronHPCapture()
00042    :G4HadronicInteraction("NeutronHPCapture")
00043   {
00044     SetMinEnergy( 0.0 );
00045     SetMaxEnergy( 20.*MeV );
00046 //    G4cout << "Capture : start of construction!!!!!!!!"<<G4endl;
00047     if(!getenv("G4NEUTRONHPDATA")) 
00048        throw G4HadronicException(__FILE__, __LINE__, "Please setenv G4NEUTRONHPDATA to point to the neutron cross-section files.");
00049     dirName = getenv("G4NEUTRONHPDATA");
00050     G4String tString = "/Capture";
00051     dirName = dirName + tString;
00052     numEle = G4Element::GetNumberOfElements();
00053 //    G4cout << "+++++++++++++++++++++++++++++++++++++++++++++++++"<<G4endl;
00054 //    G4cout <<"Disname="<<dirName<<" numEle="<<numEle<<G4endl;
00055     //theCapture = new G4NeutronHPChannel[numEle];
00056 //    G4cout <<"G4NeutronHPChannel constructed"<<G4endl;
00057     G4NeutronHPCaptureFS * theFS = new G4NeutronHPCaptureFS;
00058     //for (G4int i=0; i<numEle; i++)
00059     //{
00060 //  //    G4cout << "initializing theCapture "<<i<<" "<< numEle<<G4endl;
00061     //  theCapture[i].Init((*(G4Element::GetElementTable()))[i], dirName);
00062     //  theCapture[i].Register(theFS);
00063     //}
00064     for ( G4int i = 0 ; i < numEle ; i++ ) 
00065     {
00066        theCapture.push_back( new G4NeutronHPChannel );
00067        (*theCapture[i]).Init((*(G4Element::GetElementTable()))[i], dirName);
00068        (*theCapture[i]).Register(theFS);
00069     }
00070     delete theFS;
00071 //    G4cout << "-------------------------------------------------"<<G4endl;
00072 //    G4cout << "Leaving G4NeutronHPCapture::G4NeutronHPCapture"<<G4endl;
00073   }
00074   
00075   G4NeutronHPCapture::~G4NeutronHPCapture()
00076   {
00077     //delete [] theCapture;
00078 //    G4cout << "Leaving G4NeutronHPCapture::~G4NeutronHPCapture"<<G4endl;
00079      for ( std::vector<G4NeutronHPChannel*>::iterator 
00080            ite = theCapture.begin() ; ite != theCapture.end() ; ite++ )
00081      {
00082         delete *ite;
00083      }
00084      theCapture.clear();
00085   }
00086   
00087   #include "G4NeutronHPThermalBoost.hh"
00088   G4HadFinalState * G4NeutronHPCapture::ApplyYourself(const G4HadProjectile& aTrack, G4Nucleus& aNucleus )
00089   {
00090 
00091     if ( numEle < (G4int)G4Element::GetNumberOfElements() ) addChannelForNewElement();
00092 
00093     G4NeutronHPManager::GetInstance()->OpenReactionWhiteBoard();
00094     if(getenv("NeutronHPCapture")) G4cout <<" ####### G4NeutronHPCapture called"<<G4endl;
00095     const G4Material * theMaterial = aTrack.GetMaterial();
00096     G4int n = theMaterial->GetNumberOfElements();
00097     G4int index = theMaterial->GetElement(0)->GetIndex();
00098     if(n!=1)
00099     {
00100       xSec = new G4double[n];
00101       G4double sum=0;
00102       G4int i;
00103       const G4double * NumAtomsPerVolume = theMaterial->GetVecNbOfAtomsPerVolume();
00104       G4double rWeight;    
00105       G4NeutronHPThermalBoost aThermalE;
00106       for (i=0; i<n; i++)
00107       {
00108         index = theMaterial->GetElement(i)->GetIndex();
00109         rWeight = NumAtomsPerVolume[i];
00110         //xSec[i] = theCapture[index].GetXsec(aThermalE.GetThermalEnergy(aTrack,
00111         xSec[i] = (*theCapture[index]).GetXsec(aThermalE.GetThermalEnergy(aTrack,
00112                                                                      theMaterial->GetElement(i),
00113                                                                      theMaterial->GetTemperature()));
00114         xSec[i] *= rWeight;
00115         sum+=xSec[i];
00116       }
00117       G4double random = G4UniformRand();
00118       G4double running = 0;
00119       for (i=0; i<n; i++)
00120       {
00121         running += xSec[i];
00122         index = theMaterial->GetElement(i)->GetIndex();
00123         //if(random<=running/sum) break;
00124         if( sum == 0 || random <= running/sum ) break;
00125       }
00126       if(i==n) i=std::max(0, n-1);
00127       delete [] xSec;
00128     }
00129 
00130     //return theCapture[index].ApplyYourself(aTrack);
00131     //G4HadFinalState* result = theCapture[index].ApplyYourself(aTrack);
00132     G4HadFinalState* result = (*theCapture[index]).ApplyYourself(aTrack);
00133     aNucleus.SetParameters(G4NeutronHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargA(),G4NeutronHPManager::GetInstance()->GetReactionWhiteBoard()->GetTargZ());
00134     G4NeutronHPManager::GetInstance()->CloseReactionWhiteBoard();
00135     return result; 
00136   }
00137 
00138 const std::pair<G4double, G4double> G4NeutronHPCapture::GetFatalEnergyCheckLevels() const
00139 {
00140    //return std::pair<G4double, G4double>(10*perCent,10*GeV);
00141    return std::pair<G4double, G4double>(10*perCent,DBL_MAX);
00142 }
00143 
00144 void G4NeutronHPCapture::addChannelForNewElement()
00145 {
00146    G4NeutronHPCaptureFS* theFS = new G4NeutronHPCaptureFS;
00147    for ( G4int i = numEle ; i < (G4int)G4Element::GetNumberOfElements() ; i++ ) 
00148    {
00149       G4cout << "G4NeutronHPCapture Prepairing Data for the new element of " << (*(G4Element::GetElementTable()))[i]->GetName() << G4endl;
00150       theCapture.push_back( new G4NeutronHPChannel );
00151       (*theCapture[i]).Init((*(G4Element::GetElementTable()))[i], dirName);
00152       (*theCapture[i]).Register(theFS);
00153    }
00154    delete theFS;
00155    numEle = (G4int)G4Element::GetNumberOfElements();
00156 }

Generated on Mon May 27 17:48:59 2013 for Geant4 by  doxygen 1.4.7