PROGRAM f_rst_Summary ! A utility program in the Restore4 family. ! Reads a f.rst file OUTPUT from Restore4 (e.g., fUSM_NI_090.0Ma.rst), ! adds up the model increments of fault-offset (during time-steps so far completed) ! towards meeting goals for fault offset established by the original f.rst input file, ! and writes these as 3 new columns, to the right of the f.rst INPUT table. ! The resulting f_rst_Summary file (e.g., fUSM_NI_090.0Ma_rst_Summary.txt) ! can be opened with Excel or other spreadsheet for plotting model vs. target offsets. ! Note 2 structural differences from the similar companion program p_rst_Summary: !(1) Paleomagnetic data gives 2 targets per site (latitude anomaly, and ! vertical-axis rotation), but each fault-offset-component datum gives only one. !(2) Paleomagnetic data always span to the present, so modeling of such data are either ! COMPLETE (if the model report time is older than the datum in question), or ! INCOMPLETE (if the model report time is younger). However, for fault-offsets, there ! is a third group: WAITING (those for which the end of fault slip is older than the ! model report time). ! By Peter Bird, UCLA, April 2021 !- - - - - - - - - - - - - - - - - - - IMPLICIT NONE CHARACTER(1) :: tab CHARACTER(5) :: c5, c5_code CHARACTER(6) :: c6, dummy_f_code CHARACTER(7) :: c7 CHARACTER(8) :: c8 CHARACTER(10) :: c10 CHARACTER(47) :: c47_citation CHARACTER(50) :: c50 CHARACTER(80) :: output_filename, f_rst_filename CHARACTER(134) :: f_rst_format ! to read fault-offset data CHARACTER(134) :: f_rst_titles ! COULD be used to write fault-offset data--BUT, may contain unwanted Neotectonic-rate column-headers; ! ALSO, potentially-usable entries in the left-hand side would need to be separated by tabs. CHARACTER(134) :: line INTEGER :: f_rst_count ! number of fault-offset-components in dataset INTEGER :: delta_t_denominator, i, ios, j, high_time, length REAL*8 :: delta_t_Ma ! inferred timestep REAL*8 :: highest_t_Ma ! greatest geologic age reached so far, in Ma REAL*8 :: delta_t_numerator, duration_in_my, fraction, overlap_in_my, t1, t2, t3, t4 !-------------------------------------------------------------------------------------------------------- !ARRAYS USED TO HOLD DATA FOUND IN THE INPUT FILE: CHARACTER*6, DIMENSION(:), ALLOCATABLE :: f_code ! e.g., "F1180N"; includes offset-component byte ! (1:f_rst_count = fault-offset index) CHARACTER*50, DIMENSION(:), ALLOCATABLE :: f_fault_name ! (1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:), ALLOCATABLE :: f_target_offset_km ! (1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:), ALLOCATABLE :: f_offset_sigma_km ! (1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:), ALLOCATABLE :: f_MaxAge_Ma ! (1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:), ALLOCATABLE :: f_MinAge_Ma ! (1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:,:),ALLOCATABLE :: f_goal_rate_mmpa ! in units of km/m.y. = mm/a. !(1:num_timesteps = timestep index; 1:f_rst_count = fault-offset index) REAL*8, DIMENSION(:,:),ALLOCATABLE :: f_model_rate_mmpa ! in units of km/m.y. = mm/a. !(1:num_timesteps = timestep index; 1:f_rst_count = fault-offset index) REAL*8, DIMENSION(0:2) :: f_errors ! 3 norms of the list of errors in offset (each offet normalized by its sigma before combining): !(0:2 = L0,L1,L2 norm). !ARRAYS BEING PREPARED FOR OUTPUT AS EXTRAL COLUMNS OF THE OUTPUT FILE: INTEGER, DIMENSION(:), ALLOCATABLE :: summary_class012 ! 0 = WAITING; 1 = INCOMPLETE; 2 = COMPLETE REAL*8, DIMENSION(:), ALLOCATABLE :: summary_goal_offset_km REAL*8, DIMENSION(:), ALLOCATABLE :: summary_goal_sigma_km REAL*8, DIMENSION(:), ALLOCATABLE :: summary_model_offset_km !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (*, "(' PROGRAM f_rst_Summary:')") WRITE (*, "(' A utility program in the Restore4 family.')") WRITE (*, "(' Reads a f.rst file OUTPUT from Restore4 (e.g., fUSM_NI_090.0Ma.rst),')") WRITE (*, "(' adds up the model increments of fault-offset (during completed time-steps)')") WRITE (*, "(' towards meeting goals for fault offset in the original f.rst input file,')") WRITE (*, "(' and writes these as 3 new columns, to the right of the f.rst INPUT table.')") WRITE (*, "(' The resulting f_rst_Summary file (e.g., fUSM_NI_090.0Ma_rst_Summary.txt)')") WRITE (*, "(' can be opened with Excel or other spreadsheet to plot model vs. target offsets.')") WRITE (*, "(' Note 2 structural differences from the similar companion program p_rst_Summary:')") WRITE (*, "(' (1) Paleomagnetic data gives 2 targets per site (latitude anomaly, and')") WRITE (*, "(' vertical-axis rotation), but each fault-offset datum gives only one.')") WRITE (*, "(' (2) Paleomagnetic data span to the present, so modeling of such data are either')") WRITE (*, "(' COMPLETE (if the model report time is older than the datum in question), or')") WRITE (*, "(' INCOMPLETE (if report time is younger). However, for fault-offsets, there')") WRITE (*, "(' is a third group: WAITING (those for which the end of fault slip is older')") WRITE (*, "(' than the model report time).')") WRITE (*, "(' By Peter Bird, UCLA, April 2021')") CALL Pause() 10 WRITE (*, *) WRITE (*, "(' Enter filename of f.rst file OUTPUT from Restore: ')") READ (*, "(A)") f_rst_filename OPEN (UNIT = 1, FILE = TRIM(f_rst_filename), STATUS = "OLD", IOSTAT = ios) IF (ios /= 0) THEN WRITE (*, "(' ERROR: File ', A, ' not found (in current folder). Try again...')") TRIM(f_rst_filename) CALL Pause() GO TO 10 END IF !First time through, just count number of offset-components (f_rst_count) and greatest timestep achieved (top_time): highest_t_Ma = 0.0D0 ! just initializing before search delta_t_Ma = 0.0D0 f_rst_count = 0 delta_t_numerator = 0.0D0 delta_t_denominator = 0 READ (1, "(A)") f_rst_format READ (1, "(A)") f_rst_titles scanning: DO READ (1, "(A)", IOSTAT = ios) line IF (ios /= 0) EXIT scanning IF (line(1:1) == '*') THEN ! progress toward offset in one timestep line(1:1) = ' ' READ (line, *) t1, t2 delta_t_numerator = delta_t_numerator + (t2 - t1) delta_t_denominator = delta_t_denominator + 1 highest_t_Ma = MAX(highest_t_Ma, t2) ELSE ! header line with original INPUT line of the INPUT f_rst file READ (line, f_rst_format, IOSTAT = ios) dummy_f_code IF (ios /= 0) EXIT scanning ! perhaps input file has empty lines at end? f_rst_count = f_rst_count + 1 END IF END DO scanning CLOSE (1) ! f_rst input file (which was OUTPUT from Restore) !Infer timestep: IF (delta_t_denominator > 0) THEN delta_t_Ma = delta_t_numerator / delta_t_denominator ELSE WRITE (*, *) WRITE (*, "(' ERROR: This seems to be(?) a f.rst INPUT file, not a f.rst OUTPUT file.')") CALL Pause() STOP END IF !Infer high_time (INTEGER): high_time = NINT(highest_t_Ma / delta_t_Ma) WRITE (*, *) WRITE (*, "(' This file describes ', I5, ' fault-offset-components within the current F-E grid.')") f_rst_count WRITE (*, "(' This file contains results of ', I5, ' timesteps, extending back to ', F8.3, ' Ma.')") high_time, highest_t_Ma CALL Pause() !DIMENSION all arrays needed to record input data: ALLOCATE ( f_code(f_rst_count) ) ALLOCATE ( f_fault_name(f_rst_count) ) ALLOCATE ( f_target_offset_km(f_rst_count) ) ALLOCATE ( f_offset_sigma_km(f_rst_count) ) ALLOCATE ( f_MaxAge_Ma(f_rst_count) ) ALLOCATE ( f_MinAge_Ma(f_rst_count) ) ALLOCATE ( f_goal_rate_mmpa(high_time, f_rst_count) ) ALLOCATE ( f_model_rate_mmpa(high_time, f_rst_count) ) !Initialize all arrays (to simplify debugging): f_code = ' ' f_fault_name = ' ' f_target_offset_km = 0.0D0 f_offset_sigma_km = 0.0D0 f_MaxAge_Ma = 0.0D0 f_MinAge_Ma = 0.0D0 f_goal_rate_mmpa = 0.0D0 f_model_rate_mmpa = 0.0D0 !Re-read the input f.rst file, and memorize contents: OPEN (UNIT = 1, FILE = TRIM(f_rst_filename), STATUS = "OLD", IOSTAT = ios) READ (1, "(A)") f_rst_format READ (1, "(A)") f_rst_titles i = 0 ! index of fault-offset-component, 1:f_rst_count. Incremented whenever a header-line is read. memorizing: DO READ (1, "(A)", IOSTAT = ios) line IF (ios /= 0) EXIT memorizing IF (line(1:1) == '*') THEN ! rate of progress toward offset in one timestep, in mm/a = km/m.y. line(1:1) = ' ' READ (line, *) t1, t2, t3, t4 j = NINT(t2 / delta_t_Ma) f_model_rate_mmpa(j, i) = t3 f_goal_rate_mmpa(j, i) = t4 ELSE ! header line with original INPUT line of the INPUT f_rst file READ (line, f_rst_format, IOSTAT = ios) f_code(i+1), f_fault_name(i+1), & & f_target_offset_km(i+1), f_offset_sigma_km(i+1), & & f_MaxAge_Ma(i+1), f_MinAge_Ma(i+1) IF (ios /= 0) EXIT memorizing ! perhaps input file has empty lines at end? i = i + 1 END IF END DO memorizing CLOSE (1) ! f_rst input file (which was OUTPUT from Restore); now DONE with this. !Create extra arrays to describe Summary results. Note that MANY of these are incomplete/prorated. !Therefore we have to have logically distinct, differently-named arrays, even though the values !may be the same for those sites where the simulation is complete. ALLOCATE ( summary_class012(f_rst_count) ) ALLOCATE ( summary_goal_offset_km(f_rst_count) ) ALLOCATE ( summary_goal_sigma_km(f_rst_count) ) ALLOCATE ( summary_model_offset_km(f_rst_count) ) summary_class012 = -1 summary_goal_offset_km = 0.0D0 summary_goal_sigma_km = 0.0D0 summary_model_offset_km = 0.0D0 !Fill in values of all output arrays: DO i = 1, f_rst_count IF (f_MaxAge_Ma(i) <= highest_t_Ma) THEN summary_class012(i) = 2 ! COMPLETE ELSE IF (f_MinAge_Ma(i) < highest_t_Ma) THEN summary_class012(i) = 1 ! INCOMPLETE ELSE summary_class012(i) = 0 ! WAITING END IF summary_goal_offset_km(i) = 0.0D0 ! just initializing, before sum over timesteps... summary_model_offset_km(i) = 0.0D0 ! just initializing, before sum over timesteps... DO j = 1, high_time summary_goal_offset_km(i) = summary_goal_offset_km(i) + delta_t_Ma * f_goal_rate_mmpa(j, i) summary_model_offset_km(i) = summary_model_offset_km(i) + delta_t_Ma * f_model_rate_mmpa(j, i) END DO IF (summary_class012(i) == 2) THEN ! COMPLETE; so, use original sigmas: summary_goal_sigma_km(i) = f_offset_sigma_km(i) ELSE IF (summary_class012(i) == 1) THEN ! INCOMPLETE; use pro-rated sigmas: overlap_in_my = MIN((highest_t_Ma - 0.0D0), & ! This is the time-duration we have modeled (within the fault-offset time bracket): & (f_MaxAge_Ma(i) - f_MinAge_Ma(i)), & & (highest_t_Ma - f_MinAge_Ma(i))) duration_in_my = f_MaxAge_Ma(i) - f_MinAge_Ma(i) IF (duration_in_my > 0.0D0) THEN fraction = overlap_in_my / duration_in_my ELSE ! bad data??? fraction = 1.0D0 ENDIF summary_goal_sigma_km(i) = fraction * f_offset_sigma_km(i) ELSE ! WAITING summary_goal_sigma_km(i) = 0.0D0 ! (prorated, as above) END IF END DO !OUTPUT all results, using tabs for ease in opening in Excel, and grouping COMPLETE vs. INCOMPLETE vs. WAITING fault-offset-components: output_filename = TRIM(f_rst_filename) length = LEN_TRIM(f_rst_filename) output_filename((length-3):(length-3)) = '_' ! changing the '.' of ".rst" to the '_' of "_rst" output_filename = TRIM(output_filename) // "_Summary.txt" OPEN (UNIT = 2, FILE = TRIM(output_filename)) ! unconditional OPEN; overwrites any existing file !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -- - - WRITE (2, "('COMPLETE:')") tab = CHAR(9) ! from ASCII table line = "Code" // tab // "Fault name (c50)" // tab // "Offset(km)" // tab // "Sigma(km)" // tab // "Began(Ma)" // tab // "Ended(Ma)" // tab // & & tab // & ! (extra tab to separate input from output sections of table(s) & "Goal_km" // tab // "Sigma_km" // tab // "Model_km" // tab // & & tab // & ! (extra tab to separate the 2 output sections of table(s): linear, and log-ready. & "POS.Goal" // tab // "POS.Sigma" // tab // "POS.Model" WRITE (2, "(A)") TRIM(line) DO i = 1, f_rst_count IF (summary_class012(i) == 2) THEN ! COMPLETE line = f_code(i) WRITE (c50, "(A)") TRIM(f_fault_name(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c50)) WRITE (c10, "(F10.3)") f_target_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") f_offset_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MaxAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MinAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns (input data, vs. output summary) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") summary_goal_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_goal_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_model_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns within summary: linear vs. log-ready: !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_sigma_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_model_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (2, "(A)") TRIM(line) END IF END DO !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -- - - WRITE (2, "('INCOMPLETE:')") tab = CHAR(9) ! from ASCII table line = "Code" // tab // "Fault name (c50)" // tab // "Offset(km)" // tab // "Sigma(km)" // tab // "Began(Ma)" // tab // "Ended(Ma)" // tab // & & tab // & ! (extra tab to separate input from output sections of table(s) & "Goal_km" // tab // "Sigma_km" // tab // "Model_km" // tab // & & tab // & ! (extra tab to separate the 2 output sections of table(s): linear, and log-ready. & "POS.Goal" // tab // "POS.Sigma" // tab // "POS.Model" WRITE (2, "(A)") TRIM(line) DO i = 1, f_rst_count IF (summary_class012(i) == 1) THEN ! INCOMPLETE line = f_code(i) WRITE (c50, "(A)") TRIM(f_fault_name(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c50)) WRITE (c10, "(F10.3)") f_target_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") f_offset_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MaxAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MinAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns (input data, vs. output summary) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") summary_goal_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_goal_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_model_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns within summary: linear vs. log-ready: !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_sigma_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_model_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (2, "(A)") TRIM(line) END IF END DO !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -- - - WRITE (2, "('WAITING:')") tab = CHAR(9) ! from ASCII table line = "Code" // tab // "Fault name (c50)" // tab // "Offset(km)" // tab // "Sigma(km)" // tab // "Began(Ma)" // tab // "Ended(Ma)" // tab // & & tab // & ! (extra tab to separate input from output sections of table(s) & "Goal_km" // tab // "Sigma_km" // tab // "Model_km" // tab // & & tab // & ! (extra tab to separate the 2 output sections of table(s): linear, and log-ready. & "POS.Goal" // tab // "POS.Sigma" // tab // "POS.Model" WRITE (2, "(A)") TRIM(line) DO i = 1, f_rst_count IF (summary_class012(i) == 0) THEN ! WAITING line = f_code(i) WRITE (c50, "(A)") TRIM(f_fault_name(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c50)) WRITE (c10, "(F10.3)") f_target_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") f_offset_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MaxAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.1)") f_MinAge_Ma(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns (input data, vs. output summary) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") summary_goal_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_goal_sigma_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") summary_model_offset_km(i) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) line = TRIM(line) // tab ! extra tab to separate 2 groups of columns within summary: linear vs. log-ready: !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_goal_sigma_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) WRITE (c10, "(F10.3)") MAX(0.01D0, summary_model_offset_km(i)) line = TRIM(line) // tab // TRIM(ADJUSTL(c10)) !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - WRITE (2, "(A)") TRIM(line) END IF END DO !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -- - - CLOSE (2) WRITE (*, *) WRITE (*, "(' Job completed.')") CALL Pause() CONTAINS SUBROUTINE Pause() IMPLICIT NONE WRITE (*,"(' Press [Enter]...'\)") READ (*,*) END SUBROUTINE Pause END PROGRAM f_rst_Summary