Hello All,
I’m running some GRAS Reverse simulations to get TID values, and got stuck with apparent very low values when running a simulation for GEO with a 5 years mission. My initial checkup/validation was to use the JUICE case.
In the Forward case of JUICE I got values very close to the ones I found in the JUICE example for FMC:
Target_000 = 41.719 krad
Target_001 = 41.086 krad
Target_011 = 44.847 krad
Target_111 = 45.653 krad
However, using the RMC method, I’m arriving into very low values, in the order of ~7krad. Not sure if the way I’m accessing the results is correct though. I access the root files via a root script and go to the total_dose_vs_primary_kine:
From there I get the Total Dose via the integral, with the following piece of code:
const int num_components = 4;
const char *Component[num_components];
Component[0] = "Target_000";
Component[1] = "Target_001";
Component[2] = "Target_002";
Component[3] = "Target_003";
double Components [num_components];
double Dose_TrpE [num_components];
for(int i=0;i<num_components;i++){
Components[i] = i+1;
int Ntrpe,Nsolarp,Ngcrz1;
int N1,N2;
char File_TrpE[80],File_SolarP[80],File_GCRZ1[80];
char File1[80],File2[80];
Ntrpe = sprintf(File_TrpE,"%s_TrpE_JUICE_RMC.root",Component[i]);
N1 = sprintf(File1,"Spectrum1/TID_%s_dose",Component[i]);
N2 = sprintf(File2,"Spectrum1/TID_%s_total_dose_vs_primary_kine",Component[i]);
TFile *ftrpe = TFile::Open(File_TrpE);
TH1D *h1trpe, *h2trpe;
ftrpe->GetObject(File2, h2trpe);
Dose_TrpE[i] = h2trpe->Integral();
cout<<Component[i]<<" = "<<1e-3*Dose_TrpE[i]<<endl;
}
I multiplied by 1e-3 to get in krad, since my call function was like this:
/gras/analysis/dose/addModule TID_{TargetName}
/gras/analysis/dose/TID_{TargetName}/bookTuples true
/gras/analysis/dose/TID_{TargetName}/addVolume {TargetName}
/gras/analysis/dose/TID_{TargetName}/setUnit rad
/gras/analysis/dose/TID_{TargetName}/initialise
Am I doing anything wrong in this process? I checked and rechecked my source and my physics. The source I simply used the ones from the JUICE example, and set the normalization as:
/gras/analysis/adjoint/addUserDefinedPrimSpectrum e- 1.680080E+17 1/cm2 Arb
The full source file I’m using:
/gps/ene/type Arb
/gps/hist/type arb
/gps/hist/point 9.200000E-04 2.248400E+11
/gps/hist/point 1.200000E-03 1.922500E+11
/gps/hist/point 1.600000E-03 1.633400E+11
/gps/hist/point 2.100000E-03 1.383300E+11
/gps/hist/point 2.700000E-03 1.178500E+11
/gps/hist/point 3.500000E-03 9.848500E+10
/gps/hist/point 4.500000E-03 8.004300E+10
/gps/hist/point 5.900000E-03 6.136800E+10
/gps/hist/point 7.700000E-03 4.490700E+10
/gps/hist/point 1.000000E-02 3.108000E+10
/gps/hist/point 1.300000E-02 2.004000E+10
/gps/hist/point 1.700000E-02 1.195900E+10
/gps/hist/point 3.000000E-02 3.525600E+09
/gps/hist/point 6.100000E-02 1.240200E+09
/gps/hist/point 8.900000E-02 4.845400E+08
/gps/hist/point 1.300000E-01 1.592000E+08
/gps/hist/point 1.800000E-01 7.051400E+07
/gps/hist/point 2.700000E-01 2.380900E+07
/gps/hist/point 4.000000E-01 7.978200E+06
/gps/hist/point 6.100000E-01 1.832300E+06
/gps/hist/point 9.100000E-01 5.152000E+05
/gps/hist/point 1.280000E+00 1.692100E+05
/gps/hist/point 1.990000E+00 3.264700E+04
/gps/hist/point 2.440000E+00 1.324500E+04
/gps/hist/point 3.070000E+00 4.555100E+03
/gps/hist/point 3.970000E+00 1.018900E+03
/gps/hist/point 5.200000E+00 0.000000E+00
-> #Normalization
/gras/analysis/adjoint/addUserDefinedPrimSpectrum e- 1.680080E+17 1/cm2 Arb
-> #Definition of the adjoint source
##################################
/adjoint/DefineAdjSourceOnExtSurfaceOfAVolume {TargetName}
/adjoint/SetAdjSourceEmin 1.0000e-05 MeV
/adjoint/SetAdjSourceEmax 5.2000e+00 MeV
-> # Definition of the external source
###################################
#/adjoint/DefineSphericalExtSource 0 0 0 40 cm #
/adjoint/SetExtSourceEmax 5.2000e+00 MeV
And my physics:
/gras/physics/addPhysics rmc_em_standard
/gras/physics/addPhysics stopping
/adjoint_physics/SetEminForAdjointModels 1.0000e-05 MeV
/adjoint_physics/SetEmaxForAdjointModels 5.2000e+00 MeV
/adjoint_physics/UseBremsstrahlung true
/adjoint_physics/UseCompton true
/adjoint_physics/UseAdjointCompton false #default
/adjoint_physics/UseMS true
/adjoint_physics/UsePEEffect true
/adjoint_physics/UseAdjointPEEffect false #default
/adjoint_physics/UseProtonIonisation true
/adjoint_physics/SetSplittingFactorForPrimFwdGamma 10 #default
/run/initialize
/run/setCut 1.0000e-02 mm
/gras/physics/describe
It’s a long message, but I figured it would be easier for anyone could can help me.
Thanks in advance.