search

GPlates / gplately

GPlately is a Python package to interrogate tectonic plate reconstructions.
https://gplates.github.io/gplately/
GNU General Public License v2.0
60 stars 13 forks source link

consider feature type pygplates.FeatureType.gpml_transform when getting the transforms from topologies #227

Closed michaelchin closed 3 months ago

michaelchin commented 3 months ago

Currently PTT only consider pygplates.FeatureType.gpml_mid_ocean_ridge when getting ridges and transforms from topologies. see the code below.

https://github.com/GPlates/gplately/blob/a29e3aee09271b981a834ab95985f8a2a16baece/gplately/ptt/resolve_topologies.py#L224-L227

One of the clients mentioned that features with type pygplates.FeatureType.gpml_transform should also be considered.

I am not sure if other feature types should also be considered, such as pygplates.FeatureType.gpml_aseismic_ridge

jcannon-gplates commented 3 months ago

I just realised that pygplates.FeatureType.gpml_transform features are accounted for here - which gets extracted from self.other which comes from ptt.resolve_topologies.resolve_topologies_into_features() and is anything other that a gpml:SubductionZone or gpml:MidOceanRidge.

So I wonder if it's better/easier to have a user option in the plotting function(s) to include gpml:Transforms when plotting transform segments of gpml:MidOceanRidges (that is, whether to copy self.misc_transforms into self.transforms; and maybe also copy it into self.ridge_transforms). Though it's a bit icky having a user option just for that - most users won't know what that's for (I suppose could be an option for advanced users of some sort) - or could just decide one way or the other, and not have a user option.

@nickywright I'm guessing you wanted the option to plot gpml:Transforms the same as transform segments in gpml:MidOceanRidges? Maybe we don't bother with a user option and just change the code so that it plots them the same?

PS: I'll be away tomorrow (Weds).

nickywright commented 3 months ago

Sorry @jcannon-gplates - emails from GitHub get buried in a folder!

Currently plot_transforms / get_transforms gets the result from the ridges/transforms splitting (i.e. the result from a gpml: MidOceanRidge), and so it's a lot harder to get the result of an actual gpml: Transform type imo.

But intuitively, I would've thought plot/get transforms would have meant from actual gpml:Transforms, and the result of the ridges and transforms splitting would be called something completely different (or personally, not be used at all because clearly at 430 Ma - random time I picked - it doesn't work).

See two examples here (made with pygmt @michaelchin in case you need an example?). transform_test_0Ma.pdf transform_test_430Ma.pdf

It gets more noticeable as you go into the deep past with the more recent models (here I've used Müller et al. 2022) because I've noticed that the gpml: Transforms feature types are used more in the deeper past.

The 430 Ma plot also illustrates how the ridges and transform code isn't working very well right now at all (see plot on the left). I would imagine this would have major implications for if you were trying to calculate ridge vs transform plate boundary lengths, or perhaps crustal production at this time and it used these ridge segments/got rid of the transform segments from the gpml:MidOceanRidge, because now you're actually missing most of your actual ridge segments at this time period. I'm not sure what the codes do but I've mentioned this to Dietmar.

The code to reproduce this all (including Cartopy figures) is in this gist

jcannon-gplates commented 3 months ago

Thanks @nickywright, I agree that the ridge/transform separation doesn't appear to be working too well.

@michaelchin

In fact Dietmar just said:

I believe we should drop the functionality to separate ridges from transforms from pyGPlates and GPlately. It is dangerous because end users (like me for instance ...) may assume that it largely works, even though it does not, and if we end up using the resulting MOR lengths for calculations of MOR production, this can lead to estimates that are way off from what they ought to be.

Aside from the plot/get function changes I mention below, I think we can probably still leave the ridge/transform separation functionality within ptt.resolve_topologies but just not use it for now (eg, don't use in the above plotting). In the GPlates meeting we can discuss whether to completely remove the file ptt.separate_ridge_transform_segments.py or not. Actually I can see it's also used by ptt.ridge_spreading_rate.py so maybe we don't remove anything (but just don't expose ridge/transform separation in the plotting).

Which leads to Dietmar's next point...

For making maps that include plate boundaries, we should simply stick to always plotting MORs and their transforms together without differentiation, as their automatic separation will never work.

So we could, for now, just ignore the self.ridges and self.transforms (ie, the split ridge/transform segments) returned here and just use the self.ridge_transforms (ie, the full gpml:MidOceanRidge geometry).

That might mean a bunch of function renaming is needed. For example, plot_ridges_and_transforms() could be renamed plot_ridges() that would plot gpml:MidOceanRidge (effectively removing the existing plot_ridges() which currently only plots the split-out ridge segments of gpml:MidOceanRidge). And also having plot_transforms() plot gpml:Transform (instead of the split-out transform segments of gpml:MidOceanRidge). And same for their get function equivalents.

Regarding Dietmar's proposed way forward for crustal production along MORs:

In terms of computing crustal production rates along mid-ocean ridges, the way this used to be done many years ago was simply by masking age grids, say by masking out all crust older than 3 my, and then adding the total area of crust in km2 of age 0-3 Ma and dividing by 3 to obtain crustal production rate in km2/my (the 3 my threshold could be changed, but I am pretty sure that 1my will not be enough, as you might get no grid cells along slowly spreading ridges). We should return to doing it this way, as it will be much more reliable to estimate crustal production, despite artefacts in age grids. ... I’d like to see a new function in pyGPlates, and GPlately, to compute crustal production rates based on a set of paleo-age grids, and for our own information ...

So that would be more of a GPlately thing that pyGPlates (since age grids are involved).

Ultimately I can imagine functionality in pyGPlates that could seed points along MORs and move them over a time interval (similar to the first phase of the new age gridding workflow), Then join those moved points into a new line and calculate the intervening area (between original and moved lines). Or something along those lines (no pun intended). I'll have a think about that when I do the plate boundary divergence/convergence in pyGPlates - but that's for the longer term.

major implications for if you were trying to calculate ridge vs transform plate boundary lengths

So the computation of crustal production rates by itself probably wouldn't help here. In which case, in reference to the above-mentioned pyGPlates crustal production rate functionality, I should investigate whether that can also separate ridge/transform segments somehow (without running into the same stage rotation issues that's currently causing all our problems). If that works, it could ultimately replace the functionality of ptt.separate_ridge_transform_segments.py. But that's longer term too.

So, for now, I guess it's the crustal production age-gridding approach (in GPlately), and just hide the ridge/transform separation functionality (without removing it completely). And the renaming of plot/get functions, which might have some backward compatibility issues (eg, plot_ridges() now unexpectedly including transform segments of mid-ocean ridge features).

michaelchin commented 3 months ago

The solution proposed by the client might also close https://github.com/GPlates/gplately/issues/228

michaelchin commented 3 months ago

Thanks @jcannon-gplates @nickywright

I summarise the things which need to be done below. Let me know if I missed anything or misunderstood anything.

nickywright commented 3 months ago

@michaelchin , that sounds correct to me. Thank you!

Perhaps you can also please clarify what get_misc_transforms does? (vs get_transforms)? I don't know what a 'misc transform' is? Maybe as part of this process, the docs can be updated to say that they are 'getting' (plotting) the gpml:Transform (same for get_ridges() and gpml:MidOceanRidge etc)

michaelchin commented 3 months ago

Perhaps you can also please clarify what get_misc_transforms does? (vs get_transforms)? I don't know what a 'misc transform' is?

I don't have the answer right now. I will ask around. I am not the original author of the code. I was told to look after GPlately a few months ago. Some parts of GPlately are still unknown to me.

@GPlates/gplately-dev Does anyone know the answer to @nickywright 's question above?

jcannon-gplates commented 3 months ago

I don't know what a 'misc transform' is?

It's the gpml:Transform features (according to this in plot.py).

I would suggest we remove get_misc_transforms/plot_misc_transforms now that get_transforms/plot_transforms will now take its place.

I also wonder if we should rename get_ridges/plot_ridges to get_mid_ocean_ridges/plot_mid_ocean_ridges since it will now get/plot gpml:MidOceanRidge features ? But get_ridges/plot_ridges is fine really. Just mentioning because it enables us (in the future) to use get_ridges for ridge sections of MOR features (if we ever get the splitting of MOR into ridges/transforms to work again) - although, if that happens, maybe we'd call it something else like get_mid_ocean_ridge_sections/get_mid_ocean_transform_sections instead (but nothing to worry about now).

michaelchin commented 3 months ago

According to client DM's latest email, it seemed that he has changed his mind. He said and I quote

...I have come back to the view that we actually need to make the current approach work better...

I think we should hold further actions until DM has made up his mind.

jcannon-gplates commented 3 months ago

I think we should hold further actions until DM has made up his mind.

I see two issues here:

  1. Whether to plot (and get) gpml:MidOceanRidge (ie, the proposed plot_ridges()/get_ridges()) and plot/get gpml:Transform (ie, the proposed plot_transforms()/get_transforms()?
  2. Now that we're keeping the ability to separate ridge and transform segments within a gpml:MidOceanRidge (because, as you say, Dietmar decided we should keep it), where should those segments go (ie, in the functions in point 1 or new functions)?

I suggest we should continue to do point 1. And point 2 should be done in new functions. I think Nicky agrees...

But intuitively, I would've thought plot/get transforms would have meant from actual gpml:Transforms, and the result of the ridges and transforms splitting would be called something completely different (or personally, not be used at all because clearly at 430 Ma - random time I picked - it doesn't work).

...and then the issue of the ridge/transform splitting not working very well can hopefully be improved in the future (I'll look into it when I do plate boundary convergence/divergence in pyGPlates).

So I think we can do point 1 without waiting, but we wait before doing point 2 until we're certain we still want to split MORs into ridge/transform segments. Maybe we can that discuss in the next meeting.

jcannon-gplates commented 3 months ago
  1. Now that we're keeping the ability to separate ridge and transform segments within a gpml:MidOceanRidge (because, as you say, Dietmar decided we should keep it), where should those segments go (ie, in the functions in point 1 or new functions)?

From the recent meeting, we decided we probably won't explicitly separate into ridge and transform segments. Instead each segment will get a value in the range [0,1] for how ridge-like it is. This'll be used to calculate things like ridge lengths, crustal production, ridge spreading rates.

So the idea of separating into ridge and transform segments no longer really applies. So, when the time comes (once the necessary functionality is in pyGPlates) we can take another look at restructuring the current ridge-transform separation code (and ridge spreading code).

So I think we can now just assume that we only need to plot/get gpml:MidOceanRidge and gpml:Transform types, ie, point 1 above...

  1. ... plot (and get) gpml:MidOceanRidge (ie, the proposed plot_ridges()/get_ridges()) and plot/get gpml:Transform (ie, the proposed plot_transforms()/get_transforms()
michaelchin commented 3 months ago

OK, I will