Quasi-MC mode is now compatible with statistical mode

MintonGroup · Nov 28, 2023 · 2ce51c9 · 2ce51c9
1 parent a33cf19
commit 2ce51c9
Show file tree

Hide file tree

Showing 2 changed files with 11 additions and 4 deletions.
diff --git a/ctem/ctem/driver.py b/ctem/ctem/driver.py
@@ -169,7 +169,10 @@ def __init__(self, param_file="ctem.in", isnew=True):
                     rclist[:,0] = np.exp(interp(np.log(rclist[:,0])))
 
                     #Convert age in Ga to "interval time"
-                    rclist[:,5] = (self.user['interval'] * self.user['numintervals']) - craterproduction.Tscale(rclist[:,5], 'NPF_Moon')
+                    if (self.user['runtype'].upper() == 'STATISTICAL'):
+                        rclist[:,5] = (self.user['interval']) - craterproduction.Tscale(rclist[:,5], 'NPF_Moon')
+                    else:
+                        rclist[:,5] = (self.user['interval'] * self.user['numintervals']) - craterproduction.Tscale(rclist[:,5], 'NPF_Moon')
                     rclist = rclist[rclist[:,5].argsort()]
 
                     #Export to dat file

diff --git a/src/crater/crater_populate.f90 b/src/crater/crater_populate.f90
@@ -146,7 +146,11 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt
    clock = 0.0_DP
    finterval = 1.0_DP / real(ntotcrat,kind=DP)
    if (user%doregotrack) then
-      maxage = user%interval * user%numintervals
+      if (user%runtype .eq. 'STATISTICAL') then
+         maxage = user%interval
+      else
+         maxage = user%interval * user%numintervals
+      end if
       if (maxage < 0._DP ) then
          write(*,*) "MAJOR ERROR: Negative age!"
          stop
@@ -183,13 +187,13 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt
          end if
          if (crater%timestamp < 2330._DP) then
             if (oldGa > 0._DP) then 
-               if (user%numintervals .eq. 1) then
+               if ((user%numintervals .eq. 1) .or. (user%runtype .eq. 'STATISTICAL')) then
                   crater%timestampGa = util_t_from_scale(maxage-crater%timestamp,agemin,oldGa)
                else
                   crater%timestampGa = util_t_from_scale(maxage-crater%timestamp,1e-10_DP,oldGa)
                end if
             else
-               if (user%numintervals .eq. 1) then
+               if ((user%numintervals .eq. 1) .or. (user%runtype .eq. 'STATISTICAL')) then
                   crater%timestampGa = util_t_from_scale(maxage-crater%timestamp,agemin,maxageGa)
                else
                   crater%timestampGa = util_t_from_scale(maxage-crater%timestamp,1e-10_DP,maxageGa)