I am trying to use the General Oblique Transformation to create a projection that is a rotated Oblique Mercator ... but I cannot for the life of me figure out what on earth that thing is doing.
I started by trying to reproduce a standard Transverse Mercator centered at 0 latitude and 0 longitude, and I'm already stuck. I actually can do it successfully using a "modified" pole, but when I try to use a center point and a rotation, it does unexpected things. I want to be able to set it by the latter, not have to deal with the former.
Here's what I've tried, using pyproj to interface with proj4:
from pyproj import CRS, Transformer
geod = CRS.from_proj4("+proj=latlon")
custom1 = CRS.from_proj4("+proj=ob_tran +o_proj=tmerc +o_lat_p=90 +o_lon_p=0")
cust12geo = Transformer.from_crs(custom1, geod)
print("Custom General Oblique to mimic Transverse, using pole")
print("Origin: ", cust12geo.transform(0, 0))
print("+X: ", cust12geo.transform(111e3, 0))
print("+Y: ", cust12geo.transform(0, 111e3))
custom2 = CRS.from_proj4("+proj=ob_tran +o_proj=tmerc +o_lat_c=0 +o_lon_c=0 +o_alpha=0")
cust22geo = Transformer.from_crs(custom2, geod)
print("Custom General Oblique to mimic Transverse, using center and rotation")
print("Origin: ", cust22geo.transform(0, 0))
print("+X: ", cust22geo.transform(111e3, 0))
print("+Y: ", cust22geo.transform(0, 111e3))
What I get is the following printed out:
Custom General Oblique to mimic Transverse, using pole
Origin: (0.0, 3.508354649267438e-15)
+X: (0.9970792962983757, 3.507823426518576e-15)
+Y: (0.0, 1.0038490881892714)
Custom General Oblique to mimic Transverse, using center and rotation
Origin: (90.0, 3.508354649267438e-15)
+X: (90.0, -0.9970792962983749)
+Y: (91.00384908818927, 3.5078161886111948e-15)
According to the documentation here, I should be able to do this - I'm attempting the second of the three ways mentioned:
In addition to specifying an oblique projection, how to rotate the projection should be specified. This is done in one of three ways: Define a new pole, rotate the projection about a given point or define a new "equator" spanned by two points on the sphere.
With the above I thought that it would take a Transverse Mercator projection (with its default origin of (0,0)), and then "rotate" it around (0, 0) with a rotation of 0, meaning it shouldn't move. Why is the origin suddenly at (90, 0), with +X going south and +Y going east, instead of east and north respectively? I am clearly not understanding the inputs to this.
p.s. In case it might be relevant, I'm actually working in Cartopy and looking for a way to get something similar to its RotatedPole, but with coordinates in meters and without needing to futz with the pole (I have a center lat/lon and rotation angle I want to use). I can do that once I figure out how to work the proj4 inputs correctly!
ETA: I actually was able to get what I wanted by making a custom Cartopy Projection sub-class that exposed the "+gamma" parameter in an Oblique Mercator projection. But, I'm still interested in seeing where my understanding of "+proj=ob_tran" is faulty!
ETA: Added an issue to PROJ #4045, in case this is really a bug.
