From 79d94be8d0c8e24ac73b72207f5cc0a1396fdd64 Mon Sep 17 00:00:00 2001 From: David A Minton Date: Fri, 15 Jul 2022 16:09:08 -0400 Subject: [PATCH] Cleaned up file io for rego tracking. File names now conform to the standard CTEM style, removed redundant set of top files, and added rego files (age, melt, comp, and stack) to the redirected save files when doregotrack is turned on. --- python/ctem/ctem/driver.py | 7 ++ src/crater/crater_populate.f90 | 2 +- src/globals/module_globals.f90 | 3 +- src/io/io_read_regotrack.f90 | 51 ++++++++------- src/io/io_write_regotrack.f90 | 115 +++++++-------------------------- 5 files changed, 58 insertions(+), 120 deletions(-) diff --git a/python/ctem/ctem/driver.py b/python/ctem/ctem/driver.py index 52174cc8..afdb2aa6 100644 --- a/python/ctem/ctem/driver.py +++ b/python/ctem/ctem/driver.py @@ -59,6 +59,11 @@ def __init__(self, param_file="ctem.in", isnew=True): 'ejc': 'surface_ejc.dat', 'pos': 'surface_pos.dat', 'time': 'surface_time.dat', + 'stack': 'surface_stacknum.dat', + 'rego': 'surface_rego.dat', + 'melt': 'surface_melt.dat', + 'comp': 'surface_comp.dat', + 'age' : 'surface_age.dat', 'ocum': 'ocumulative.dat', 'odist': 'odistribution.dat', 'pdist': 'pdistribution.dat', @@ -309,6 +314,8 @@ def process_output(self): if (self.user['savetruelist'].upper() == 'T'): self.redirect_outputs(['tcum'], 'dist') self.redirect_outputs(['impmass'], 'misc') + if (self.user['saverego'].upper() == 'T') : + self.redirect_outputs(['stack','rego','age','melt','comp'], 'rego') diff --git a/src/crater/crater_populate.f90 b/src/crater/crater_populate.f90 index c3c779c3..50ac5760 100644 --- a/src/crater/crater_populate.f90 +++ b/src/crater/crater_populate.f90 @@ -89,7 +89,7 @@ subroutine crater_populate(user,surf,crater,domain,prod,production_list,vdist,nt ! doregotrack & age simulation test real(DP) :: melt, clock, age, thick real(SP),dimension(user%gridsize, user%gridsize) :: agetop - real(SP),dimension(60) :: agetot + real(SP),dimension(MAXAGEBINS) :: agetot type(regolisttype),pointer :: current => null() real(DP) :: age_resolution diff --git a/src/globals/module_globals.f90 b/src/globals/module_globals.f90 index f29b7600..6da3fb25 100644 --- a/src/globals/module_globals.f90 +++ b/src/globals/module_globals.f90 @@ -63,6 +63,7 @@ module module_globals real(DP),parameter :: PF = 5.0e9 ! The shock pressure exceeding the hungoit elastic limit of common geological materials in Pa. (5 GPa) real(DP),parameter :: RAD_GP = 20.0_DP ! The maximum radial position of producing impact glass spherules within a transient crater (unit of crater radii, crater%rad) + type regodatatype real(SP),dimension(MAXAGEBINS) :: age real(DP) :: thickness @@ -275,7 +276,7 @@ module module_globals character(*),parameter :: TIMEFILE = 'surface_time.dat' character(*),parameter :: EJCOVFILE = 'surface_ejc.dat' character(*),parameter :: DEMFILE = 'surface_dem.dat' -character(*),parameter :: REGOFILE = 'surface_regotop.dat' +character(*),parameter :: REGOFILE = 'surface_rego.dat' character(*),parameter :: MELTFILE = 'surface_melt.dat' character(*),parameter :: COMPFILE = 'surface_comp.dat' character(*),parameter :: STACKNUMFILE = 'surface_stacknum.dat' diff --git a/src/io/io_read_regotrack.f90 b/src/io/io_read_regotrack.f90 index ed93300c..3acdbb14 100644 --- a/src/io/io_read_regotrack.f90 +++ b/src/io/io_read_regotrack.f90 @@ -28,14 +28,14 @@ subroutine io_read_regotrack(user,surf) ! Internals integer(I4B), parameter :: LUN=7 - integer(I4B), parameter :: LUM=8 - integer(I4B), parameter :: LUP=9 - integer(I4B), parameter :: LUC=10 - integer(I4B), parameter :: LUA=11 + integer(I4B), parameter :: FMELT = 10 + integer(I4B), parameter :: FREGO = 11 + integer(I4B), parameter :: FCOMP = 12 + integer(I4B), parameter :: FAGE = 13 real(DP),dimension(user%gridsize,user%gridsize) :: regotop,melt,comp - real(SP),dimension(user%gridsize,user%gridsize,60) :: age + real(SP),dimension(user%gridsize,user%gridsize,MAXAGEBINS) :: age integer(I4B),dimension(user%gridsize,user%gridsize) :: stacks_num - real(DP), dimension(:),allocatable :: regotopi,melti,compi,agei + real(DP), dimension(:), allocatable :: regotopi,melti,compi,agei type(regodatatype) :: newsurfi integer(I4B) :: ioerr,i,j,k,q,itmp integer(kind=8) :: recsize @@ -46,34 +46,34 @@ subroutine io_read_regotrack(user,surf) ! Open a file for obtaining the number of stacks that is stored in each linked list ioerr = 0 recsize = sizeof(itmp) * user%gridsize * user%gridsize - open(LUP,file=STACKNUMFILE,status='old',form='unformatted',recl=recsize,access='direct',iostat=ioerr) + open(LUN,file=STACKNUMFILE,status='old',form='unformatted',recl=recsize,access='direct',iostat=ioerr) if (ioerr/=0) then write(*,*) 'Error! Cannot read file ',trim(adjustl(STACKNUMFILE)) stop end if - read(LUP,rec=1) stacks_num - close(LUP) + read(LUN,rec=1) stacks_num + close(LUN) ! Open files for regolith thickness/melt fraction - open(LUN,file=REGOFILE,status='old',form='unformatted',iostat=ioerr) + open(FREGO,file=REGOFILE,status='old',form='unformatted',iostat=ioerr) if (ioerr/=0) then write(*,*) 'Error! Cannot read file ',trim(adjustl(REGOFILE)) stop end if - open(LUC,file=COMPFILE,status='old',form='unformatted',iostat=ioerr) + open(FCOMP,file=COMPFILE,status='old',form='unformatted',iostat=ioerr) if (ioerr/=0) then write(*,*) 'Error! Cannot read file ',trim(adjustl(COMPFILE)) stop end if - open(LUM,file=MELTFILE,status='old',form='unformatted',iostat=ioerr) + open(FMELT,file=MELTFILE,status='old',form='unformatted',iostat=ioerr) if (ioerr/=0) then write(*,*) 'Error! Cannot read file ',trim(adjustl(MELTFILE)) stop end if - open(LUA,file=AGEFILE,status='old',form='unformatted',iostat=ioerr) + open(FAGE,file=AGEFILE,status='old',form='unformatted',iostat=ioerr) if (ioerr/=0) then write(*,*) 'Error! Cannot read file ',trim(adjustl(MELTFILE)) stop @@ -89,18 +89,18 @@ subroutine io_read_regotrack(user,surf) allocate(regotopi(stacks_num(i,j))) allocate(compi(stacks_num(i,j))) allocate(melti(stacks_num(i,j))) - allocate(agei(60 * stacks_num(i,j))) + allocate(agei(MAXAGEBINS * stacks_num(i,j))) do k=1,stacks_num(i,j) - read(LUN) regotop(i,j) + read(FREGO) regotop(i,j) regotopi(k) = regotop(i,j) - read(LUC) comp(i,j) + read(FCOMP) comp(i,j) compi(k) = comp(i,j) - read(LUM) melt(i,j) + read(FMELT) melt(i,j) melti(k) = melt(i,j) - read(LUA) age(i,j,:) - do q=1,60 - agei(60*k - (60-q)) = age(i,j,q) + read(FAGE) age(i,j,:) + do q=1,MAXAGEBINS + agei(MAXAGEBINS*k - (MAXAGEBINS-q)) = age(i,j,q) end do end do @@ -108,8 +108,8 @@ subroutine io_read_regotrack(user,surf) newsurfi%thickness = regotopi(k) newsurfi%comp = compi(k) newsurfi%meltfrac = melti(k) - do q=1,60 - newsurfi%age(q) = agei(60*k-(60-q)) + do q=1,MAXAGEBINS + newsurfi%age(q) = agei(MAXAGEBINS*k-(MAXAGEBINS-q)) end do call util_push(surf(i,j)%regolayer,newsurfi) end do @@ -118,8 +118,9 @@ subroutine io_read_regotrack(user,surf) end do end do - close(LUN) - close(LUC) - close(LUM) + close(FMELT) + close(FREGO) + close(FCOMP) + close(FAGE) return end subroutine io_read_regotrack diff --git a/src/io/io_write_regotrack.f90 b/src/io/io_write_regotrack.f90 index a6f3c2e1..5f522321 100644 --- a/src/io/io_write_regotrack.f90 +++ b/src/io/io_write_regotrack.f90 @@ -27,14 +27,12 @@ subroutine io_write_regotrack(user,surf) ! Regotrack Internals integer(I4B) :: i,j,k - integer(I4B), parameter :: LUN=7 - integer(I4B), parameter :: LUM=8 - integer(I4B), parameter :: LUC=9 - integer(I4B), parameter :: LUA=10 + integer(I4B), parameter :: LUN = 7 + integer(I4B), parameter :: FMELT = 10 + integer(I4B), parameter :: FREGO = 11 + integer(I4B), parameter :: FCOMP = 12 + integer(I4B), parameter :: FAGE = 13 type(regolisttype),pointer :: current => null() - real(DP),dimension(user%gridsize,user%gridsize) :: regotop,comp,melt - real(SP),dimension(user%gridsize,user%gridsize) :: agetop - real(SP),dimension(user%gridsize,user%gridsize,60) :: age integer(I4B),dimension(user%gridsize,user%gridsize) :: stacks_num integer(kind=8) :: recsize real(DP) :: dtmp @@ -42,104 +40,35 @@ subroutine io_write_regotrack(user,surf) integer(I4B) :: itmp real(DP),dimension(user%gridsize,user%gridsize) :: comptop, rego real(DP),dimension(:),allocatable :: marehisto - real(DP) :: mare, z - ! Mixing - !real(DP),parameter :: zmix = 0.0_DP - !real(DP) :: z, zmare - - ! Output multiple "comphisto" files - character(len=255) :: fname - character(len=255), parameter :: clockfile = 'tic-toc.dat' - integer(I4B) :: tictoc - logical :: exist - ! Executable code - ! Output mulitple "comphisto" files - inquire(file=clockfile, exist=exist) - if (exist) then - open(LUN,file=clockfile,status='old') - read(LUN,*) tictoc - else - write(*,*) clockfile,' is missing!' - end if - tictoc = tictoc + 1 - close(LUN) - open(LUN,file=clockfile,status='replace') - write(LUN,*) tictoc - close(LUN) - - write(fname,'(a,i4.4)') 'surface_melt',tictoc - open(LUN,file=fname,status='replace',form='unformatted') - write(fname,'(a,i4.4)') 'surface_rego',tictoc - open(LUM,file=fname,status='replace',form='unformatted') - write(fname,'(a,i4.4)') 'surface_comp',tictoc - open(LUC,file=fname,status='replace',form='unformatted') - write(fname,'(a,i4.4)') 'surface_age',tictoc - open(LUA,file=fname,status='replace',form='unformatted') + open(FMELT,file=MELTFILE,status='replace',form='unformatted') + open(FREGO,file=REGOFILE,status='replace',form='unformatted') + open(FCOMP,file=COMPFILE,status='replace',form='unformatted') + open(FAGE,file=AGEFILE,status='replace',form='unformatted') + stacks_num(:,:) = 0 do j=1,user%gridsize do i=1,user%gridsize - stacks_num(i,j) = 0 current => surf(i,j)%regolayer - comptop(i,j) = current%regodata%comp - rego(i,j) = current%regodata%thickness - agetop(i,j) = current%regodata%age(1) do - if (.not. associated(current)) exit - stacks_num(i,j) = stacks_num(i,j) + 1 - regotop(i,j) = current%regodata%thickness - comp(i,j) = current%regodata%comp - melt(i,j) = current%regodata%meltfrac - age(i,j,:)= current%regodata%age(:) - write(LUM) regotop(i,j) - write(LUC) comp(i,j) - write(LUN) melt(i,j) - write(LUA) age(i,j,:) - current => current%next + if (.not. associated(current)) exit ! We've reached the bottom of the linked list + stacks_num(i,j) = stacks_num(i,j) + 1 + write(FMELT) current%regodata%meltfrac + write(FREGO) current%regodata%comp + write(FCOMP) current%regodata%meltfrac + write(FAGE) current%regodata%age(:) + current => current%next end do end do end do - close(LUN) - close(LUM) - close(LUC) - close(LUA) - - recsize = sizeof(dtmp) * user%gridsize * user%gridsize - open(LUN,file='agetop.dat',status='replace',form='unformatted',recl=recsize,access='direct') - write(LUN,rec=1) agetop - close(LUN) - - recsize = sizeof(dtmp) * user%gridsize * user%gridsize - open(LUN,file='comptop.dat',status='replace',form='unformatted',recl=recsize,access='direct') - write(LUN,rec=1) comptop - close(LUN) - - allocate(marehisto(user%gridsize)) - !open(LUN,file='comphisto',status='replace') ! Output comphisto one time - write(fname,'(a,i4.4)') 'comphisto',tictoc - open(LUN,file=fname,status='replace') - - do i=1,user%gridsize - marehisto(i) = 0.0_DP - mare = 0.0_DP - do j=1,user%gridsize - mare = mare + comptop(i,j) - end do - marehisto(i) = mare/real(user%gridsize) - write(LUN,*) real(i-user%gridsize/2)*user%pix/1000.0,marehisto(i)*100.0 - end do - close(LUN) - deallocate(marehisto) - - recsize = sizeof(dtmp) * user%gridsize * user%gridsize - open(LUN,file='regotop.dat',status='replace',form='unformatted',recl=recsize,access='direct') - write(LUN,rec=1) rego - close(LUN) + close(FMELT) + close(FREGO) + close(FCOMP) + close(FAGE) recsize = sizeof(itmp) * user%gridsize * user%gridsize - write(fname,'(a,i4.4)') 'surface_stacknum',tictoc - open(LUN,file=fname,status='replace',form='unformatted',recl=recsize,access='direct') + open(LUN,file=STACKNUMFILE,status='replace',form='unformatted',recl=recsize,access='direct') write(LUN,rec=1) stacks_num close(LUN)