#include <MyCfSource.hh>
Public Methods | |
| MyCfSource (int newIsotope=k.Isotope, int newRmaxgen=k.Rmaxgen) | |
| MyCfSource constructor. More... | |
| ~MyCfSource () | |
| MyCfSource destructor. | |
| MyCfSource (const MyCfSource &CfSource) | |
| MyCfSource copy constructor. | |
| MyCfSource & | operator= (const MyCfSource &rhs) |
| MyCfSource overloaded = operator. | |
| int | GetNumCfSources () const |
| Returns the number of Cf sources. | |
| int | GetNumNeutron (int &n) const |
| Returns the neutron multiplicity of each source. | |
| int | GetNumGamma (int &n) const |
| Returns the prompt photon multiplicity of each source. | |
| double | GetCfNeutronEnergy (int &n, int &nn) const |
| Returns the energy of the neutrons produced by each Cf decay. More... | |
| double | GetCfNeutronTime (int &n, int &nn) const |
| double | GetCfGammaEnergy (int &n, int &nn) const |
| Returns the energy of the gammas produced by each Cf decay. More... | |
| double | GetCfGammaTime (int &n, int &nn) const |
| double | GetCfVertexR (int &n) const |
| Returns the radius of the neutron produced by each Cf decay. More... | |
| double | GetCfVertexCosTheta (int &n) const |
| Returns cos(theta) of the neutron produced by each Cf decay. More... | |
| double | GetCfVertexPhi (int &n) const |
| Returns the phi of the neutron produced by each Cf decay. More... | |
Static Public Methods | |
| float | Cf252NeutronSpectrum (const float &x) |
| The probability density of the prompt neutrons from the Cf decay as a function of neutron energy. | |
| float | Cf252GammaMultiplicity (const int &x) |
| The probability density of the prompt gammas from the Cf decay as a function of integer gamma multiplicity (exact). | |
| float | Cf252GammaMultiplicityFit (const float &x) |
| The probability density of the prompt gammas from the Cf decay as a function of float gamma multiplicity (fit). | |
| float | Cf252GammaSpectrum (const float &x) |
| The probability density of the prompt gammas from the Cf decay as a function of gamma energy. | |
The supported isotope is 252, which is also the default isotope given ReactorConstants. The neutrons are generated with the correct multiplicity and each neutron is given an energy by the Cf252NeutronSpectrum function, which the probability density of the produced neutrons as a function of neutron energy. The prompt photons are generated by the Cf252GammaMultiplicityFit and Cf252GammaSpectrum functions.
The number and locations of the Cf sources are read in from cf252_position.txt and parsed by MyCfSource to yield the source vertex R, cos(theta), and phi.
Definition at line 26 of file MyCfSource.hh.
|
||||||||||||
|
MyCfSource constructor. Called with the desired Cf isotope: 252, which is the default isotope given in ReactorConstants and the maximum allowed radius for generating particles in the simulation, Rmaxgen (default = 200 cm). If the source in cf252_position.txt is located with R > Rmaxgen the source is ignored and an error message is printed. Definition at line 19 of file MyCfSource.cpp. References Cf252GammaMultiplicityFit, Cf252GammaSpectrum, Cf252NeutronSpectrum, ReactorConstants::Ndim, and ReactorConstants::pi.
00020 {
00021
00022 Isotope = newIsotope;
00023 Rmaxgen = newRmaxgen;
00024
00025 // cout << "MyCfSource::MyCfSource(" << Isotope << "," << Rmaxgen
00026 // << ")" << endl;
00027
00028 // Cf252
00029 if(Isotope == 252){
00030
00031 // setup the probability density as a function of energy
00032 bool static first=1;
00033 float static fspace[200];
00034 float static mspace[200];
00035 float static gspace[200];
00036
00037 if(first){
00038 float xlo = 0.;
00039 float xhi = 50.;
00040 funlxp_(Cf252NeutronSpectrum,fspace,xlo,xhi);
00041 xhi = 25.;
00042 funlxp_(Cf252GammaMultiplicityFit,mspace,xlo,xhi);
00043 funlxp_(Cf252GammaSpectrum,gspace,xlo,xhi);
00044 first = 0;
00045 }
00046
00047 // input the source number and position
00048 fstream cf252_position("data/cf252_position.txt", ios::in);
00049 //cout << "Cf input file: data/cf252_position.txt" << endl;
00050
00051 Nsources = 0;
00052 char paramLine[256];
00053 string dummy;
00054 bool done = false;
00055 bool numSet = false;
00056 int positionNum = 0;
00057
00058 while(!done){
00059
00060 // '#' specifies a comment line -- ignore it
00061 if(cf252_position.peek() == '#'){
00062
00063 cf252_position.getline(paramLine, 256);
00064 continue;
00065 }
00066
00067 // 'E'ND or 'e'nd specifies end of file -- set done flag to true
00068 if(cf252_position.peek() == 'e' || cf252_position.peek() == 'E'){
00069
00070 if(Nsources == 0) cout << " Cf252 ERROR: "
00071 << "number of cf252 sources = 0" << endl;
00072 done = true;
00073 continue;
00074 }
00075
00076 // 'N'SOURCES = or 'n'sources = specifies number of sources
00077 if(cf252_position.peek() == 'n' || cf252_position.peek() == 'N'){
00078
00079 if(numSet){
00080
00081 cf252_position.getline( paramLine, 256 );
00082 continue;
00083 }
00084 cf252_position >> dummy;
00085 while(cf252_position.peek() == ' ' || cf252_position.peek() == '=')
00086 cf252_position.ignore(1);
00087
00088 cf252_position >> Nsources;
00089 if(Nsources == 0){
00090 cout << " Cf252 ERROR: number of cf252 sources = 0" << endl;
00091 done = true;
00092 continue;
00093 }
00094
00095 numSet = true;
00096 continue;
00097 }
00098
00099 // 'P'OSITION = or 'p'osition = specifies position of the sourse(s)
00100 if(cf252_position.peek() == 'p' || cf252_position.peek() == 'P'){
00101
00102 if(Nsources == 0 || positionNum >= Nsources){
00103
00104 cout << " Cf252 ERROR: too many positions " << positionNum+1
00105 << " given for the number of sources " << Nsources << endl;
00106 cf252_position.getline( paramLine, 256 );
00107 positionNum++;
00108 continue;
00109 }
00110 cf252_position >> dummy;
00111 while(cf252_position.peek() == ' ' || cf252_position.peek() == '=')
00112 cf252_position.ignore(1);
00113
00114 cf252_position >> R[positionNum] >> zR[positionNum]
00115 >> phiR[positionNum];
00116
00117 // check for unphysical source positions
00118 if(R[positionNum] > Rmaxgen){
00119 cout << " Cf252 ERROR: R[" << positionNum << "] of "
00120 << R[positionNum] << " > Rmaxgen" << endl;
00121 done = true;
00122 continue;
00123 }
00124 if(zR[positionNum] < -1. || zR[positionNum] > 1.){
00125 cout << " Cf252 ERROR: CosTheta["
00126 << positionNum << "] of " << zR[positionNum] << endl;
00127 done = true;
00128 continue;
00129 }
00130 if(phiR[positionNum] > 2*k.pi){
00131 cout << " Cf252 ERROR: Phi[" << positionNum << "] of "
00132 << phiR[positionNum] << " > 2*pi" << endl;
00133 done = true;
00134 continue;
00135 }
00136
00137 positionNum++;
00138 continue;
00139 }
00140
00141 // ignore extra whitespace
00142 if(cf252_position.peek() == ' ' || cf252_position.peek() == '\n'){
00143
00144 cf252_position.ignore(1);
00145 continue;
00146 }
00147 }
00148
00149 if(positionNum < Nsources){
00150 cout << " Cf252 ERROR: not enough positions " << positionNum
00151 << " given for the number of sources " << Nsources << endl;
00152 }
00153 //
00154 // neutron multiplicity distribution (via Los Alamo manual)
00155 //
00156 double p[maxNeutron+1] = {0.002,0.028,0.155,0.428,
00157 0.732,0.917,0.983,0.998,1.000};
00158 //
00159 for(int n=0;n<Nsources;n++){
00160 //cout << "Source " << n+1 << endl << " Cf Source Position = "
00161 // << R[n] << ", " << zR[n] << ", " << phiR[n] << endl;
00162
00163 // pick a neutron multiplicity
00164 bool passed = false;
00165 int len;
00166 while(!passed){
00167 len = 1;
00168 float r[len];
00169 ranlux_(r,len);
00170
00171 for(int i=0;i<maxNeutron+1; i++){
00172
00173 if(p[i] > r[0] && !passed){
00174 Nneutron[n] = i;
00175 if(Nneutron[n])passed = true;
00176 }
00177 }
00178 }
00179 //cout << " " << Nneutron[n] << " neutrons" << endl;
00180 //
00181 // pick an energy for each neutron
00182 //
00183 for(int nn=0;nn<Nneutron[n]; nn++){
00184 len = 1;
00185 funlux_(fspace,neutronE[n][nn],len);
00186 //cout << " Cf neutron energy = " << neutronE[n][nn] << endl;
00187 Tneutron[n][nn] = 0.;
00188 }
00189 /*
00190 The total energy in the prompt gammas is a function of the
00191 material and the number of prompt neutrons produced; these
00192 formulae come from T. Valentine's summary
00193 */
00194 float Z = 98.;
00195 float A = 252.;
00196 //
00197 double avge = -1.33 + 119.6*pow(Z,0.3333)/A;
00198 double phi = 2.51 - 1.13e-5*Z*Z*sqrt(A);
00199 double etot = phi*Nneutron[n]+4.0;
00200 //
00201 // use the Brunson multiplicity (also via Valentine)
00202 //
00203 len = 1;
00204 float ng = 0.;
00205 funlux_(mspace,ng,len);
00206 Ngamma[n] = (int)ng;
00207 //cout << " " << Ngamma[n] << " photons" << endl;
00208 //
00209 // pick an energy for each gamma
00210 //
00211 double tote = 0.;
00212 for(int nn=0;nn<Ngamma[n]; nn++){
00213 len = 1;
00214 funlux_(gspace,gammaE[n][nn],len);
00215 //cout << " 252Cf photon energy = " << gammaE[n][nn] << endl;
00216 tote += gammaE[n][nn];
00217
00218 len = 2;
00219 float rv[len];
00220 ranlux_(rv,len);
00221 //
00222 // 80% of gammas have T_1/2 = 0.01 ns and 20% have T_1/2 = 1 ns
00223 //
00224 double halflife;
00225 if(rv[0] < 0.8)
00226 halflife = 0.01;
00227 else
00228 halflife = 1.0;
00229 Tgamma[n][nn] = halflife*log(1/rv[1]);
00230 }
00231 //cout << " total energy = " << tote << endl;
00232 } // done looping over sources
00233
00234 // zero out remaining energy and position arrays
00235 for(int n=0;n<Nsources;n++){
00236 for(int nn=Nneutron[n];nn<maxNeutron;nn++){
00237 neutronE[n][nn] = 0.;
00238 Tneutron[n][nn] = 0.;
00239 }
00240 for(int nn=Ngamma[n];nn<maxGamma;nn++){
00241 gammaE[n][nn] = 0.;
00242 Tgamma[n][nn] = 0.;
00243 }
00244 }
00245 for(int n=Nsources;n<k.Ndim;n++){
00246 Nneutron[n] = 0;
00247 Ngamma[n] = 0;
00248 R[n] = 0.;
00249 zR[n] = 0.;
00250 phiR[n] = 0.;
00251 for(int nn=0;nn<maxNeutron;nn++){
00252 neutronE[n][nn] = 0.;
00253 Tneutron[n][nn] = 0.;
00254 }
00255 for(int nn=0;nn<maxGamma;nn++){
00256 gammaE[n][nn] = 0.;
00257 Tgamma[n][nn] = 0.;
00258 }
00259 }
00260
00261 } // done with Cf 252
00262
00263 if(Isotope == 255){
00264
00265 Nsources = 0;
00266 for(int n=Nsources;n<k.Ndim;n++){
00267 Nneutron[n] = 0;
00268 Ngamma[n] = 0;
00269 R[n] = 0.;
00270 zR[n] = 0.;
00271 phiR[n] = 0.;
00272 for(int nn=0;nn<maxNeutron;nn++){
00273 neutronE[n][nn] = 0.;
00274 Tneutron[n][nn] = 0.;
00275 }
00276 for(int nn=0;nn<maxGamma;nn++){
00277 gammaE[n][nn] = 0.;
00278 Tgamma[n][nn] = 0.;
00279 }
00280 }
00281 } // done with Cf 255
00282 }
|
|
||||||||||||
|
Returns the energy of the gammas produced by each Cf decay. Called with the integer index for each gamma in the gammaE array, from 0 to the total number of sources and 0 to the total number of gammas/source. Definition at line 85 of file MyCfSource.hh. Referenced by ReactorEvent::ReactorEvent.
00085 {return gammaE[n][nn];}
|
|
||||||||||||
|
Returns the energy of the neutrons produced by each Cf decay. Called with the integer index for each neutron in the neutronE array, from 0 to the total number of sources and 0 to the total number of neutrons/source. Definition at line 78 of file MyCfSource.hh. Referenced by ReactorEvent::ReactorEvent.
00078 {return neutronE[n][nn];}
|
|
|
Returns cos(theta) of the neutron produced by each Cf decay. Called with the integer index for each neutron in the zR array, from 0 to the total number of sources. Definition at line 96 of file MyCfSource.hh. Referenced by ReactorEvent::ReactorEvent.
00096 {return zR[n];}
|
|
|
Returns the phi of the neutron produced by each Cf decay. Called with the integer index for each neutron in the phiR array, from 0 to the total number of sources. Definition at line 101 of file MyCfSource.hh. Referenced by ReactorEvent::ReactorEvent.
00101 {return phiR[n];}
|
|
|
Returns the radius of the neutron produced by each Cf decay. Called with the integer index for each neutron in the R array, from 0 to the total number of sources. Definition at line 91 of file MyCfSource.hh. Referenced by ReactorEvent::ReactorEvent.
00091 {return R[n];}
|
1.2.14 written by Dimitri van Heesch,
© 1997-2002