diff --git a/python/ctem/ctem/__init__.py b/python/ctem/ctem/__init__.py index 4774c6af..b9aa95ac 100644 --- a/python/ctem/ctem/__init__.py +++ b/python/ctem/ctem/__init__.py @@ -1,3 +1,4 @@ # from ctem.ctem_io_readers import * # from ctem.ctem_io_writers import * # from ctem.ctem_driver import * +from ctem.frontend import * \ No newline at end of file diff --git a/python/ctem/ctem/frontend.py b/python/ctem/ctem/frontend.py new file mode 100644 index 00000000..6340d5ec --- /dev/null +++ b/python/ctem/ctem/frontend.py @@ -0,0 +1,46 @@ +import numpy +import os +import subprocess +import shutil +from ctem import io + +class Simulation: + """ + This is a class that defines the basic CTEM simulation object + """ + def __init__(self, param_file="ctem.in"): + currentdir = os.getcwd() + self.parameters = { + 'restart': None, + 'runtype': None, + 'popupconsole': None, + 'saveshaded': None, + 'saverego': None, + 'savepres': None, + 'savetruelist': None, + 'seedn': 1, + 'totalimpacts': 0, + 'ncount': 0, + 'curyear': 0.0, + 'fracdone': 1.0, + 'masstot': 0.0, + 'interval': 0.0, + 'numintervals': 0, + 'pix': -1.0, + 'gridsize': -1, + 'seed': 0, + 'maxcrat': 1.0, + 'shadedminhdefault': 1, + 'shadedmaxhdefault': 1, + 'shadedminh': 0.0, + 'shadedmaxh': 0.0, + 'workingdir': currentdir, + 'ctemfile': param_file, + 'datfile': 'ctem.dat', + 'impfile': None, + 'sfdcompare': None, + 'sfdfile': None + } + self.parameters = io.read_param(self.parameters) + + return \ No newline at end of file diff --git a/python/ctem/ctem/io.py b/python/ctem/ctem/io.py new file mode 100644 index 00000000..664d8bb8 --- /dev/null +++ b/python/ctem/ctem/io.py @@ -0,0 +1,87 @@ +import numpy +import os + +def real2float(realstr): + """ + Converts a Fortran-generated ASCII string of a real value into a numpy float type. Handles cases where double precision + numbers in exponential notation use 'd' or 'D' instead of 'e' or 'E' + + Parameters + ---------- + realstr : string + Fortran-generated ASCII string of a real value. + + Returns + ------- + : float + The converted floating point value of the input string + """ + return float(realstr.replace('d', 'E').replace('D', 'E')) + +def read_param(parameters): + # Read and parse ctem.in file + inputfile = os.path.join(parameters['workingdir'], parameters['ctemfile']) + + # Read ctem.in file + print('Reading input file ' + parameters['ctemfile']) + with open(inputfile, 'r') as fp: + lines = fp.readlines() + + # Process file text + for line in lines: + fields = line.split() + if len(fields) > 0: + if ('pix' == fields[0].lower()): parameters['pix'] = real2float(fields[1]) + if ('gridsize' == fields[0].lower()): parameters['gridsize'] = int(fields[1]) + if ('seed' == fields[0].lower()): parameters['seed'] = int(fields[1]) + if ('sfdfile' == fields[0].lower()): parameters['sfdfile'] = fields[1] + if ('impfile' == fields[0].lower()): parameters['impfile'] = fields[1] + if ('maxcrat' == fields[0].lower()): parameters['maxcrat'] = real2float(fields[1]) + if ('sfdcompare' == fields[0].lower()): parameters['sfdcompare'] = fields[1] + if ('interval' == fields[0].lower()): parameters['interval'] = real2float(fields[1]) + if ('numintervals' == fields[0].lower()): parameters['numintervals'] = int(fields[1]) + if ('popupconsole' == fields[0].lower()): parameters['popupconsole'] = fields[1] + if ('saveshaded' == fields[0].lower()): parameters['saveshaded'] = fields[1] + if ('saverego' == fields[0].lower()): parameters['saverego'] = fields[1] + if ('savepres' == fields[0].lower()): parameters['savepres'] = fields[1] + if ('savetruelist' == fields[0].lower()): parameters['savetruelist'] = fields[1] + if ('runtype' == fields[0].lower()): parameters['runtype'] = fields[1] + if ('restart' == fields[0].lower()): parameters['restart'] = fields[1] + if ('shadedminh' == fields[0].lower()): + parameters['shadedminh'] = real2float(fields[1]) + parameters['shadedminhdefault'] = 0 + if ('shadedmaxh' == fields[0].lower()): + parameters['shadedmaxh'] = real2float(fields[1]) + parameters['shadedmaxhdefault'] = 0 + + # Test values for further processing + if (parameters['interval'] <= 0.0): + print('Invalid value for or missing variable INTERVAL in ' + inputfile) + if (parameters['numintervals'] <= 0): + print('Invalid value for or missing variable NUMINTERVALS in ' + inputfile) + if (parameters['pix'] <= 0.0): + print('Invalid value for or missing variable PIX in ' + inputfile) + if (parameters['gridsize'] <= 0): + print('Invalid value for or missing variable GRIDSIZE in ' + inputfile) + if (parameters['seed'] == 0): + print('Invalid value for or missing variable SEED in ' + inputfile) + if (parameters['sfdfile'] is None): + print('Invalid value for or missing variable SFDFILE in ' + inputfile) + if (parameters['impfile'] is None): + print('Invalid value for or missing variable IMPFILE in ' + inputfile) + if (parameters['popupconsole'] is None): + print('Invalid value for or missing variable POPUPCONSOLE in ' + inputfile) + if (parameters['saveshaded'] is None): + print('Invalid value for or missing variable SAVESHADED in ' + inputfile) + if (parameters['saverego'] is None): + print('Invalid value for or missing variable SAVEREGO in ' + inputfile) + if (parameters['savepres'] is None): + print('Invalid value for or missing variable SAVEPRES in ' + inputfile) + if (parameters['savetruelist'] is None): + print('Invalid value for or missing variable SAVETRUELIST in ' + inputfile) + if (parameters['runtype'] is None): + print('Invalid value for or missing variable RUNTYPE in ' + inputfile) + if (parameters['restart'] is None): + print('Invalid value for or missing variable RESTART in ' + inputfile) + + return parameters diff --git a/python/ctem/tests/testio/ctem.in b/python/ctem/tests/testio/ctem.in new file mode 100755 index 00000000..e57310dc --- /dev/null +++ b/python/ctem/tests/testio/ctem.in @@ -0,0 +1,70 @@ +! CTEM Input file + + +! Testing input. These are used to perform non-Monte Carlo tests. +testflag T ! Set to T to create a single crater with user-defined impactor properties +testimp 5.934e3 ! 93km crater ! Diameter of test impactor (m) +testvel 15.0e3 ! Velocity of test crater (m/s) +testang 90.0 ! Impact angle of test crater (deg) - Default 90.0 +testxoffset 0.0e0 ! x-axis offset of crater center from grid center (m) - Default 0.0 +testyoffset 0.0e0 ! y-axis offset of crater center from grid center (m) - Default 0.0 +tallyonly F ! Tally the craters without generating any craters +testtally F + +! IDL driver in uts +interval 1.0 +numintervals 1 ! Total number of intervals +restart F ! Restart a previous run +impfile NPFextrap.dat ! Impactor SFD rate file (col 1: Dimp (m), col 2: ! impactors > D (m**(-2) y**(-1)) +popupconsole F ! Pop up console window every output interval +saveshaded T ! Output shaded relief images +saverego F ! Output regolith map images +savepres F ! Output simplified console display images (presentation-compatible images) +savetruelist T ! Save the true cumulative crater distribution for each interval (large file size) +runtype statistical ! single: craters accumulate in successive intervals + ! statistical: surface is reset between intervals +sfdcompare FassettCounts.txt + +! CTEM required inputs +seed 23790 ! Random number generator seed +gridsize 1000 ! Size of grid in pixels +numlayers 10 ! Number of perched layers +pix 0.3468773099817e3 ! Pixel size (m) - Copernicus DEM compariso +mat rock ! Material (rock or ice) +! Bedrock scaling parameters +mu_b 0.55e0 ! Experimentally derived parameter for bedrock crater scaling law +kv_b 0.20e0 ! Experimentally derived parameter for bedrock crater scaling law +trho_b 2250.0e0 ! Target density (bedrock) (kg/m**3) +ybar_b 0.0e6 ! Bedrock strength (Pa) +! Regolith scaling parameters +mu_r 0.55e0 ! Experimentally derived parameter for regolith crater scaling law +kv_r 0.20e0 ! Experimentally derived parameter for regolith crater scaling law +trho_r 2250.0e0 ! Target density (regolith) (kg/m**3) +ybar_r 0.00e6 ! Regolith strength (Pa) +! Body parameters +gaccel 1.62e0 ! Gravitational acceleration at target (m/s**2) +trad 1737.35e3 ! Target radius (m) +prho 2500.0e0 ! Projectile density (kg/m**3) +sfdfile production.dat ! Impactor SFD file (col 1: Dimp (m), col 2: ! impactors > D +velfile lunar-MBA-impactor-velocities.dat ! Impactor velocity dist file + +! Seismic shaking input (required if seismic shaking is set to T, otherwise ignored) +doseismic F ! Perform seismic shaking calculations with each impact - Default F + +! Optional inputF These have internally set default values that work reasonable well. Comment them out with +deplimit 9e99 ! Depth limit for craters (m) - Default is to ignore. +maxcrat 1.00e0 ! Fraction of gridsize that maximum crater can be - Default 1.0 +killatmaxcrater F ! Stop the run if a crater larger than the maximum is produced - Default F +basinimp 35.0e3 ! Size of impactor to switch to lunar basin scaling law - Default is to ignore +docollapse T ! Do slope collapse - Default T +dosoftening F ! Do ejecta softening - Default T +doangle T ! Vary the impact angle. Set to F to have only vertical impacts - Default T +Kd1 0.0001 +psi 2.000 +fe 2.00 +ejecta_truncation 2.0 +dorays F +superdomain F + +dorealistic T + diff --git a/python/ctem/tests/testio/testio.py b/python/ctem/tests/testio/testio.py new file mode 100644 index 00000000..ae329313 --- /dev/null +++ b/python/ctem/tests/testio/testio.py @@ -0,0 +1,7 @@ +import ctem + +d = dir(ctem) + +sim = ctem.Simulation() +for k, v in sim.parameters.items(): + print(k, v) \ No newline at end of file