# High resolution x,y dataset for Mollweide projection convert to latitude and longitude coordinates

I have a 2D array that's 144000 x 64000 that has some values that are irrelevant to the question. Let's say the array can be used like an X, Y coordinate system. I somehow need to get the Latitude and Longitude based from the Mollweide projection.

I took a look at the Wiki inverse transformation formula, but I believe I'm doing the math incorrectly. Let's say I have some X,Y coordinate (130000, 50000), how would I get the latitude and longitude from it?

Note: the actual array is composed of several arrays that are 4000 x 4000, with the number of arrays being 36 x 16.

My current math is based on this SO question:

R = (36 * 4000) / 2 / sqrt(2) = 50911.68825 radians
theta = arcsin(50000 / R / sqrt(2) = 43.9829 radians
latitude = arcsin[(2 * theta + sin(2 * theta)) / pi] = DOMAIN ERROR
=> since (2 * theta + sin(2 * theta)) / pi = 28.31853329
longitude = 130000 * pi / (2 * R * sqrt(2) * cos(theta)) = 3.941594775


So neither latitude nor longitude make ANY sense. Am I supposed to convert theta to something below 2 * Pi? What math am I doing incorrectly?

• Are you sure you want to do this the hard way? An open source lib like GDAL (implementing OSR) would be much faster and a whole lot more trustworthy. Commented Feb 25, 2020 at 4:51
• @MichaelStimson Please point me in the right direction :) I've been using rasterio for some population density .tif files and it returns these numpy arrays. Is there anything specific that GDAL can offer that helps? Rasterio also uses GDAL and this is my first time touching anything related to GIS. Commented Feb 25, 2020 at 4:54
• Isn't rasterio built on GDAL? Either way a point geometry can be generated, a spatial reference defined then projected, have a read of gdal.org/tutorials/osr_api_tut.html and see if that helps. This one too, for python pcjericks.github.io/py-gdalogr-cookbook/projection.html Commented Feb 25, 2020 at 4:56
• Also, the projection origin is at the center of your image, not a corner, so you need to take that into account as well in the calculation by subtracting half of your image width and height to your pixel x,y coordinates before passing your coordinates in the mollweide formulas. Ex: to evaluate pixel (130000,50000) use (130000-72000,50000-32000) -> (58000,18000) Commented Feb 25, 2020 at 5:31
• My apologies, sqrt(2) was correct Commented Feb 25, 2020 at 5:39

Also, since the projection origin is at the center of the image, not a corner, you have to take this into account as well before passing your x,y values in the Mollweide formulas. Ex: to evaluate pixel (130000,50000) use (130000-72000,50000-36000) -> (58000,18000)