Best map projection for Europe - geostationary satellite?

[lucalelli]

Dear GMT users,

I’d like to brainstorm about the right map projection for a geostationary satellite swath. I am experimenting with various geographic projections and I haven’t found a satisfactory solution yet. The region is roughly “-55/75/25/80”. The data will always have this shape.

I would like to achieve two results: the first is not to distort the data and the satellite swath, the greater the curvature at high latitudes. The second is to have a final image that is pleasing to the eye without taking up too much white space between the data and the frame.

In the first example, the Winkel Tripel, as well as Robinson or similar, still distorts the rows and the satellite swath is strecthed, while keeping the white background to an acceptable size

A slightly better result for the strecthing is achieved with the Orthographic projection

or with a Lambert conic conformal

But with the latter two projections the frame is not rectangular and the resulting figure is not really usable. Experimenting outside GMT, one solution with cartopy allows the following result heading in the right direction.

Long story short:

1. What would be the best projection not to distort the data and the swath?
2. Can GMT mix a certain projection (say Orthgraphic) with a rectangular frame?

Thanks

[Joaquim]

Hi Luca,

Some comments and a suggestion.

• geostationary orbits do not have swaths (the sensor is always looking to the same place on Earth). From your figures I would say that the data is from a near-polar sol-synchronous orbit.
• It simply not possible to have no map distortions when one represent spherical data on a plane. What we can do is minimize them.

As for the suggestion, try the oblique Mercator projection with the Y axis aligned with the orbit inclination (around 80 degrees from equatorial plane).

[Andreas]

You could also use -R...+r for some of your projections? Then you’ll always get a rectangular map.

From The -R option:

3. -Rxlleft/ylleft/xuright/yuright +r. This form is useful for map projections that are oblique, making meridians and parallels poor choices for map boundaries. Here, we instead specify the lower left corner and upper right corner geographic coordinates, followed by the modifier +r. This form guarantees a rectangular map even though lines of equal longitude and latitude are not straight lines.

[mkononets]

@lucalelli

Can GMT mix a certain projection (say Orthgraphic) with a rectangular frame?

6. GMT Map Projections — GMT 6.6.0 documentation

there are several examples demonstrating how one can get rectangular frames for projections with non-rectangular coordinate grids, using the -R...+r syntax mentioned by @Andreas using lower left and upper right corner coordinates.

[lucalelli]

Here I implemented the suggestions of @Joaquim (left) and @Andreas @mkononets (right). The title of the figures shows the settings of basemap.

As Joaquim correctly pointed out, there is no way to eliminate the curvature. For this reason I labeled upfront in my mind any cylindric projection as evil. But I was wrong. The azimuthal rotation in oblique Mercator offers also flexibility and delivers quite a decent result. The output of both appraches is very similar or at least I cannot say what is the best.

Thanks for having a look at that!

[Joaquim]

Actually, I was expecting to see a ~perfect rectangle with the oblique Mercator. I think you don’t use the +r in this case.
Also, what data is that? Why do we see that curvature pattern pointing to a center not far outside the northern edge of the figure?

[lucalelli]

Hi @Joaquim , I added +r to the oblique Mercator but the curvature is not only still present, it has been exacerbated.

Data itself are my own retrievals. The image is a composite stitched from 12 consecutive hourly scans as below

The geolocation is from the MTG-S/S4 platform.

I agree with you that the data curvature seems odd, also given that these simple lines with cartopy have produced the image above

projection=ccrs.Orthographic(central_longitude=0.0, central_latitude=40.0))
plt.scatter(... , transform=ccrs.PlateCarree())
op.set_extent([-20, 40, 20, 85],crs=ccrs.PlateCarree() )

but I have not found a solution yet.

[Joaquim]

@lucalelli

.

Luca, I didn’t mean to say that +r would make the curvature goes away. In fact I’m puzzled of where it’s comming from. Specially since the mapplotlib figure does not show it. … but on the other hand plt.scatter(... makes think that it’s makking a scatter plot, whilst the GMT one looks like plotting a grid.

Any direct link where we can download a data sample? ESA sites are a pain to find data.