"Drift corrrection" for gplates; copies latest positions of each plate to a near-present time and replaces all flowlines based on mid-ocean ridges

drift_correct.py

import sys
import os
import shutil
import fileinput
import pygplates

#Update rotation.rot file to copy latest rotation to a near-present day age
Adjust_rotation = True

#Redraw all flowlines based on mid-ocean ridge features
Fix_flowlines = True

#name of file with mid-ocean ridge features
mor_file = 'mid-ocean ridges.gpml'

#name of file with flowlines features (note: will be erased and replaced)
flow_file = 'flowlines.gpml'

#rotations at or before this age won't be copied,
#and flowlines will be drawn back only to this age
max_age = 2000.

#rotations will be copied to this age
min_age = 1.0


if Adjust_rotation:
    print('Adjusting rotation...')
    if os.path.exists('rotation.rot_backup3'):
        shutil.copy('rotation.rot_backup3', 'rotation.rot_backup4')    
    if os.path.exists('rotation.rot_backup2'):
        shutil.copy('rotation.rot_backup2', 'rotation.rot_backup3')
    if os.path.exists('rotation.rot_backup1'):
        shutil.copy('rotation.rot_backup1', 'rotation.rot_backup2')
    shutil.copy('rotation.rot', 'rotation.rot_backup1')
    prev = 0
    prev_line = ''
    for line in fileinput.input('rotation.rot', inplace=True):
        new_line = line
        i=0
        i1=0
        i2=0
        start=False
        space=False
        num = 0.0
        for c in line:
            if c == ' ' or c == '/t':
                if not space:
                    space = True
                if start:
                    i2 = i
                    age = line[i1:i2]
                    num = float(age)
                    if num == 0.0:
                        prev = 1
                        print(line, end = '')
                    elif num == min_age:
                        prev_line = line
                        prev = 2
                    else:
                        if num > max_age:
                            if prev == 2:
                                print(prev_line, end = '')
                            prev = 0
                        elif prev > 0:
                            new_line = line[:i1] + str(min_age) + line[i2:]
                            prev = 0
                            print(new_line, end = '')
                        print(line, end = '')
                    break
            elif space:
                start = True
                space = False
                i1 = i
            i+=1
    fileinput.close()

if Fix_flowlines:
    print('Fixing flowlines...')
    rotation_mod = pygplates.RotationModel('rotation.rot')
    mors = pygplates.FeatureCollection(mor_file)
    newflows = pygplates.FeatureCollection()

    for mor in mors:
        left = mor.get_left_plate()
        right = mor.get_right_plate()
        start, end = mor.get_valid_time()
        if end > max_age:
            continue
        mor_points = []
        pygplates.reconstruct(mor, rotation_mod, mor_points, 0)
        if len(mor_points) < 1:
            continue
        points1 = mor_points[0].get_reconstructed_geometry()
        points2 = points1.get_points()

        times = range(0, int(min(max_age, start)+10), 10)

        flow = pygplates.Feature.create_flowline(
            pygplates.MultiPointOnSphere(points2),
            times,
            valid_time = (min(max_age, start), end),
            left_plate = left,
            right_plate = right)
        newflows.add(flow)

    os.remove(flow_file)
    newflows.write(flow_file)

print('Done!')