# Inconsistent result when reconstructing a feature with rotaion angle and the function "reconstruct"

Hi all,

I have met a problem when reconstructing a feature using the Muller et al. 2019 model.

I did the reconstruction using the “reconstruct” function as below:
reconstructed_geometry = []
pygplates.reconstruct(feature_to_be_reconstructed, rotation_model, reconstructed_geometry,
reconstruction_time=time)
print reconstructed_geometry[0].get_reconstructed_geometry().to_lat_lon()

I also tried a different way:
plateid = feature_to_be_reconstructed.get_reconstruction_plate_id()
rotation_angle = rotation_model.get_rotation(time, plateid, 0, 0)
print (rotation_angle*(assigned_point_features_10[i].get_geometry())).to_lat_lon()

I got two different results. This is confusing. I thought the two methods are equivalent. Am I missing something?

Let me know if anyone can help me with this. Thanks!

How is `assigned_point_features_10` related to `feature_to_be_reconstructed`? Are they the same?

Yes. They are the same.

The initial purpose is to partition a point feature into plates to get the plate id. The point was at (4.322765776364536, -172.00856066713223) at 130 Ma. And the partitioning function gives a plate id of 926 and its present-day location at (45.601579512218066, 30.569660725405726).

When I reconstruct the point back to 130 Ma using the pygplates.reconstruct function, I get the same location as the initial point. But we I compute the rotation angle using the time 130 Ma and the plate id 926 and then rotate it back to 130 Ma, I get a different result.

I suspect they might have different geometries (even though they represent the same feature).

It sounds like you’re using a non-zero (ie, 130Ma) reconstruction time in the partitioning function. In which case I’m assuming `feature_to_be_reconstructed` contains that point at 130Ma (`4.322765776364536, -172.00856066713223`) and `assigned_point_features_10` is the output of the partitioning function which contains the present day point location (`45.601579512218066, 30.569660725405726`).

Is that right? If so, then they contain different geometries. If not, then can you post some more code (such as how you partition)?

And maybe confirm whether the feature geometries are the same with:

``````print feature_to_be_reconstructed[0].get_geometry().to_lat_lon()
print assigned_point_features_10[i].get_geometry().to_lat_lon()
``````

Blockquote
It sounds like you’re using a non-zero (ie, 130Ma) reconstruction time in the partitioning function. In which case I’m assuming `feature_to_be_reconstructed` contains that point at 130Ma (`4.322765776364536, -172.00856066713223` ) and `assigned_point_features_10` is the output of the partitioning function which contains the present day point location (`45.601579512218066, 30.569660725405726` ).
Blockquote
Yes, exactly. But I could not understand why they contain different geometries. Below is a clearer demo.

Thanks for the clarification, I see what you mean now.

So the problem, as I now understand it, is you’re reconstructing `assigned_point_feature[0]` two different ways and getting two different results. I had previously thought you were comparing reconstructions of `assigned_point_feature[0]` (containing present day geometry) with reconstructions of `point_feature` (containing geometry already reconstructed to 130Ma). But that’s not the case here.

In that case, I think the problem is due to plate `926` having a non-zero rotation at present day (relative to anchor plate 0). And since you are explicitly specifying a `from_time` of `0` in:

``````rotation_angle = rotation_model.get_rotation(time, pid, 0, 0)
``````

…(ie, the first `0`) then the non-zero present day rotation is excluded from the returned rotation. In other words, it’s returning the rotation of plate `926` at 130Ma relative to 0Ma (which excludes the fact that the geometry stored in your feature also undergoes a non-zero rotation just to get to 0Ma).

``````rotation_angle = rotation_model.get_rotation(time, pid)
``````

…and see if you now get the same reconstructed results.

The reason for this distinction is explained at the bottom of RotationModel.get_rotation() in the following note:

Explicitly setting `from_time` to zero can give a different result than not specifying `from_time` at all if the moving plate (or fixed plate) has a non-zero finite rotation at present day (relative to the anchor plate). However all present-day finite rotations should ideally be zero (identity), so typically there should not be a difference.

…and so it’s probably better to avoid non-zero present day rotations whenever possible as it usually leads to confusing results.

The above note assumes you’re using the latest pyGPlates public release (revision 28) which introduced that distinction to avoid inconsistencies that started appearing whenever there were non-zero present day rotations.

Thank you, John.

I can understand the issue now. I just made a test to delete the zeros in function call. The results are consistent now.