#include <iomanip>
#include "globals.hh"
#include "G4SystemOfUnits.hh"
#include "G3toG4.hh"
#include "G3EleTable.hh"
#include "G3MatTable.hh"
#include "G4Material.hh"
#include "G4Isotope.hh"
Go to the source code of this file.
Functions | |
void | PG4gsmixt (G4String *tokens) |
void | G4gsmixt (G4int imate, G4String name, G4double *a, G4double *z, G4double dens, G4int nlmat, G4double *wmat) |
void G4gsmixt | ( | G4int | imate, | |
G4String | name, | |||
G4double * | a, | |||
G4double * | z, | |||
G4double | dens, | |||
G4int | nlmat, | |||
G4double * | wmat | |||
) |
Definition at line 73 of file G4gsmixt.cc.
References G4Material::AddElement(), FatalException, G3Ele, G3Mat, G4Exception(), G3EleTable::GetEle(), and G3MatTable::put().
00075 { 00076 // in Geant3: 00077 // After a call with ratios by number (negative number of elements), 00078 // the ratio array is changed to the ratio by weight, so all successive 00079 // calls with the same array must specify the number of elements as 00080 // positive 00081 G4int i=0; 00082 if (nlmat<0) { 00083 // in case of proportions given in atom counts (nlmat<0), 00084 // the wmat[i] are converted to weight fractions 00085 G4double aMol = 0.; 00086 for (i=0; i<std::abs(nlmat); i++) { 00087 // total molecular weight 00088 aMol += wmat[i]*a[i]; 00089 } 00090 if (aMol == 0.) { 00091 G4String text = "G4mixt: Total molecular weight in " + name + " = 0."; 00092 G4Exception("G4gsmixt()", "G3toG40016", FatalException, text); 00093 return; 00094 } 00095 for (i=0; i<std::abs(nlmat); i++) { 00096 // weight fractions 00097 wmat[i] = wmat[i]*a[i]/aMol; 00098 } 00099 } 00100 00101 // create material with given number of components 00102 // (elements) 00103 00104 G4Material* material 00105 = new G4Material(name, dens*g/cm3, std::abs(nlmat)); 00106 for (i=0; i<std::abs(nlmat); i++) { 00107 // add units 00108 // G4Element* element = G4Element(z[i], a[i]*g/mole, name); 00109 G4Element* element = G3Ele.GetEle(z[i]); 00110 material->AddElement(element, wmat[i]); 00111 } 00112 00113 // add the material to the List 00114 G3Mat.put(imate, material); 00115 }
void PG4gsmixt | ( | G4String * | tokens | ) |
Definition at line 43 of file G4gsmixt.cc.
References G4String::data(), G3fillParams(), G4gsmixt(), Ipar, PTgsmixt, Rpar, and Spar.
Referenced by G3CLEval().
00044 { 00045 // fill the parameter containers 00046 G3fillParams(tokens,PTgsmixt); 00047 00048 // interpret the parameters 00049 G4String name = Spar[0].data(); 00050 G4int imate = Ipar[0]; 00051 G4int nlmat = Ipar[1]; 00052 //G4double dens = Rpar[0]*g/cm3; 00053 G4double dens = Rpar[0]; 00054 G4double *a = Rpar + 1; 00055 G4double *z = Rpar + 1+std::abs(nlmat); 00056 G4double *wmat = Rpar + 1 + 2*std::abs(nlmat); 00057 /* 00058 for (int i=0; i<std::abs(nlmat); i++){ 00059 //Rpar[i]=Rpar[i]*g/mole; 00060 Rpar[i]=Rpar[i]; 00061 }; 00062 */ 00063 G4gsmixt(imate,name,a,z,dens,nlmat,wmat); 00064 }