| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

DsG4GdCaptureGammas Class Reference

#include <DsG4GdCaptureGammas.h>

List of all members.


Public Member Functions

 DsG4GdCaptureGammas ()
 ~DsG4GdCaptureGammas ()
G4ReactionProductVector * GetGammas ()
std::vector< double > GetEnergy ()

Private Attributes

TFile * specfile
TH1F * capspec

Detailed Description

Definition at line 18 of file DsG4GdCaptureGammas.h.


Constructor & Destructor Documentation

DsG4GdCaptureGammas::DsG4GdCaptureGammas (  ) 

Definition at line 25 of file DsG4GdCaptureGammas.cc.

00026 {
00027     //First back up the current TDirectory. This is hard-wired, but at 
00028     //this point we assume that an output root file has been opened at the
00029     //very beginning of the run. I will try to find a prettier fix  
00030     //later. jianglai 05/09/2006
00031     TDirectory *dirKeep = gDirectory;
00032 
00033     // The experimental spectrum is saved in gdcap_spec.root.
00034     // This file should be placed in the data/ directory.
00035     // It will get copied to the InstallArea/${tag}/data directory when
00036     // the package is 'make'd.
00037     // The InstallArea/${tag}/data which should be linked to with the $DAYA_DATA_DIR
00038     // environment variable.
00039 
00040     std::string specfilename = "gdcap_spec.root";
00041     const char* prefix = getenv("DAYA_DATA_DIR");
00042     if(prefix) {
00043         if(strlen(prefix)>0) { 
00044             std::string newname = prefix;
00045             newname += "/";
00046             newname += specfilename;
00047             specfilename = newname;
00048         }
00049     }
00050 
00051     specfile = TFile::Open(specfilename.c_str());
00052 
00053     if(!specfile)
00054         G4cout << "Can not open gdcap_spec.root. "
00055                << "Make sure DAYA_DATA_DIR is set to it's location." << G4endl;
00056     if(!(specfile->IsOpen()))
00057         G4cout << "Can not open gdcap_spec.root. "
00058                << "Make sure DAYA_DATA_DIR is set to it's location." << G4endl;
00059 
00060     capspec = (TH1F*)(specfile->Get("h1"));
00061 
00062     // return to the TDirectory used before opening TFile
00063     if (dirKeep) {
00064         dirKeep->cd();
00065     }  
00066 }

DsG4GdCaptureGammas::~DsG4GdCaptureGammas (  ) 

Definition at line 68 of file DsG4GdCaptureGammas.cc.

00069 {
00070   /*
00071     if(specfile && specfile->IsOpen()) {
00072         specfile->Close();
00073         G4cout <<"gdcap_spec.root closed"<<G4endl;
00074     }
00075     else 
00076         G4cout <<"gdcap_spec.root is not closed."<<G4endl;
00077   */
00078 }


Member Function Documentation

G4ReactionProductVector * DsG4GdCaptureGammas::GetGammas (  ) 

Definition at line 80 of file DsG4GdCaptureGammas.cc.

00081 {
00082     G4ReactionProductVector* theGammas = new G4ReactionProductVector;
00083     vector<double> energy = GetEnergy();
00084     for(unsigned int i=0; i<energy.size(); i++) {
00085 
00086         G4ReactionProduct* theOne = new G4ReactionProduct;
00087         theOne->SetDefinition( G4Gamma::Gamma() );
00088 
00089         // Get the gammas direction. 
00090         // Isotropic emission.
00091         G4double costheta = 2.*G4UniformRand()-1;
00092         G4double theta = acos(costheta);
00093         G4double phi = twopi*G4UniformRand();
00094         G4double sinth = sin(theta);
00095         theOne->SetTotalEnergy( energy[i] );
00096         G4ThreeVector temp(energy[i]*sinth*cos(phi), 
00097                            energy[i]*sinth*sin(phi),
00098                            energy[i]*cos(theta) );
00099         theOne->SetMomentum( temp ) ;
00100         theGammas->push_back(theOne);
00101     }
00102     return theGammas;
00103 }

vector< double > DsG4GdCaptureGammas::GetEnergy (  ) 

Definition at line 157 of file DsG4GdCaptureGammas.cc.

00158 {
00159     vector<double> Energy;
00160     G4double TotalEnergy;
00161     if (G4UniformRand()<0.815)
00162         TotalEnergy = 7.93; // captured by Gd157
00163     else
00164         TotalEnergy = 8.53; // captured by Gd155
00165     G4double energy_sum = 0.0;
00166     G4int    gamma_num = 0;
00167     G4double gamma_energy[100]={0};
00168     G4double energy_tmp;
00169   
00170   sample:
00171   
00172 
00173     //Use the CLHEP random engine to sample the distribution instead
00174     //Jianglai 10/01/2006
00175     //G4double energy=capspec->GetRandom();// sample from the experimental spectrum.
00176     G4double energy = getRandomCLHEP(capspec);
00177   
00178     energy_sum = energy_sum+energy;
00179     gamma_num = gamma_num+1;
00180     gamma_energy[gamma_num] = energy;
00181 
00182     if(energy_sum>TotalEnergy) {
00183         if((energy_sum-TotalEnergy)<0.8) {
00184             energy_tmp = gamma_energy[gamma_num-1]+energy-(energy_sum-TotalEnergy);
00185             if(energy_tmp>7.4) {
00186                 energy_sum = energy_sum-energy;
00187                 gamma_num = gamma_num-1;
00188                 goto sample;
00189             }
00190             gamma_num = gamma_num-1;
00191             gamma_energy[gamma_num]=energy_tmp;
00192             goto goon;
00193         }
00194         energy_sum = energy_sum-energy;
00195         gamma_num = gamma_num-1;
00196         goto sample;
00197     }
00198     goto sample;
00199 
00200   goon:
00201   
00202     // set the sequence
00203     for(int i=gamma_num;i>1;i--) {
00204         for(int j=1;j<i;j++) {
00205             if(gamma_energy[j]<gamma_energy[j+1]) {
00206                 energy_tmp = gamma_energy[j+1];
00207                 gamma_energy[j+1] = gamma_energy[j];
00208                 gamma_energy[j] = energy_tmp;
00209             }
00210         }
00211     }
00212     if(gamma_num<=3)
00213         goto  ok;
00214   
00215     if(gamma_energy[gamma_num]<0.8) {
00216         energy_tmp = gamma_energy[gamma_num]+gamma_energy[gamma_num-1];
00217         if(energy_tmp>7.4)
00218             goto ok;
00219         gamma_energy[gamma_num-1]=energy_tmp;
00220         gamma_num=gamma_num-1;
00221     }
00222    
00223   ok:
00224   
00225     for(int k=1;k<=gamma_num;k++)
00226         Energy.push_back(gamma_energy[k]); 
00227    
00228     return Energy;
00229 }


Member Data Documentation

TFile* DsG4GdCaptureGammas::specfile [private]

Definition at line 29 of file DsG4GdCaptureGammas.h.

TH1F* DsG4GdCaptureGammas::capspec [private]

Definition at line 30 of file DsG4GdCaptureGammas.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 Mon Apr 11 20:53:27 2011 for DetSim by doxygen 1.4.7