#include <MyCfSource.hh>
Public Methods | |
| MyCfSource (int newIsotope=ISOTOPE, double newRmaxgen=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, and Cf252NeutronSpectrum.
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*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 const int 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;
00198 double phi;
00199 double etot;
00200 //
00201 avge = -1.33 + 119.6*pow(Z,0.3333)/A;
00202 phi = 2.51 - 1.13e-5*Z*Z*sqrt(A);
00203 etot = phi*Nneutron[n]+4.0;
00204 //
00205 // use the Brunson multiplicity (also via Valentine)
00206 //
00207 len = 1;
00208 float ng = 0.;
00209 funlux_(mspace,ng,len);
00210 Ngamma[n] = (int)ng;
00211 //cout << " " << Ngamma[n] << " photons" << endl;
00212 //
00213 // pick an energy for each gamma
00214 //
00215 double tote = 0.;
00216 for(int nn=0;nn<Ngamma[n]; nn++){
00217 int len = 1;
00218 funlux_(gspace,gammaE[n][nn],len);
00219 //cout << " 252Cf photon energy = " << gammaE[n][nn] << endl;
00220 tote += gammaE[n][nn];
00221
00222 const int len2 = 2;
00223 float rv[len2];
00224 ranlux_(rv,len2);
00225 //
00226 // 80% of gammas have T_1/2 = 0.01 ns and 20% have T_1/2 = 1 ns
00227 //
00228 double halflife;
00229 if(rv[0] < 0.8)
00230 halflife = 0.01;
00231 else
00232 halflife = 1.0;
00233 Tgamma[n][nn] = halflife*log(1/rv[1]);
00234 }
00235 //cout << " total energy = " << tote << endl;
00236 } // done looping over sources
00237
00238 // zero out remaining energy and position arrays
00239 for(int n=0;n<Nsources;n++){
00240 for(int nn=Nneutron[n];nn<maxNeutron;nn++){
00241 neutronE[n][nn] = 0.;
00242 Tneutron[n][nn] = 0.;
00243 }
00244 for(int nn=Ngamma[n];nn<maxGamma;nn++){
00245 gammaE[n][nn] = 0.;
00246 Tgamma[n][nn] = 0.;
00247 }
00248 }
00249 for(int n=Nsources;n<Ndim;n++){
00250 Nneutron[n] = 0;
00251 Ngamma[n] = 0;
00252 R[n] = 0.;
00253 zR[n] = 0.;
00254 phiR[n] = 0.;
00255 for(int nn=0;nn<maxNeutron;nn++){
00256 neutronE[n][nn] = 0.;
00257 Tneutron[n][nn] = 0.;
00258 }
00259 for(int nn=0;nn<maxGamma;nn++){
00260 gammaE[n][nn] = 0.;
00261 Tgamma[n][nn] = 0.;
00262 }
00263 }
00264
00265 } // done with Cf 252
00266
00267 if(Isotope == 255){
00268
00269 Nsources = 0;
00270 for(int n=Nsources;n<Ndim;n++){
00271 Nneutron[n] = 0;
00272 Ngamma[n] = 0;
00273 R[n] = 0.;
00274 zR[n] = 0.;
00275 phiR[n] = 0.;
00276 for(int nn=0;nn<maxNeutron;nn++){
00277 neutronE[n][nn] = 0.;
00278 Tneutron[n][nn] = 0.;
00279 }
00280 for(int nn=0;nn<maxGamma;nn++){
00281 gammaE[n][nn] = 0.;
00282 Tgamma[n][nn] = 0.;
00283 }
00284 }
00285 } // done with Cf 255
00286 }
|
|
||||||||||||
|
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