A smarter selection of the region when using -R...r with non-cartesian projections

KristofKoch commented 4 years ago

When using gridded data in non-cartesian projections with rectangular borders I run into some problems when selecting the region of interest. The current mechanism is to specify lower left and upper right corner of the region of interest when using -R...r. This creates some difficulties in higher latitudes when using it with projections like Lambert conic conformal -JL as illustrated in this crude mockup:

region_R_r

Explanation:

the green stars are the lower left and upper richt corner as given in -R...r
the blue outline is the (crude) representation of the projection -JL with the blue dashed lines representing the center of the projection
the black frame is the plotted area
the red fill shows areas without gridded data if same -R...r is used for select/grdimage/etc

Problems:

it is hard to get a plot centered on the projection center
it usually requires several iterations to get a region definition large enough to fill the plot area with gridded data. Not very convenient with large files.

Feature request, open for debate:

add a behavior similar to -R without r where one specifies the northern and southern most latitudes and can use one of the following longitude specifications:
- longitudes at mid-latitude taken from projection center in -J, proposed -R...rm for mid (see example 1 below)
- max longitudes on plot, only a point in the corner, proposed -R...rx for max (see example 2 below)
- min longitudes visible visible on longest latitude on plot, there are some more but these must be on the plot, proposed -R...rn for min (see example 3 below)
have GMT calculate the parameters for the "expanded" projection in the background for use with select/grdimage/etc to cover the red areas.

Example 1: region_middle The longitudes are referenced to the latitude of the projection center.

Example 2: region_max The longitudes are the min longitudes to be present on the plot.

Example 3: region_min The longitudes are the max longitudes displayed on the longest line of latitude.

In all cases the parameters for the "expanded" projection are calculated in the background.

While making this feature request it came to my mind that this could be seen as two different requests. One for the different ways to define the region and one for the calculation of the expanded projection parameters in the background. Please advise if I should split them.

anbj commented 4 years ago

This is a great suggestion, and I will vouch for it. My experience is the same; -R...r is often painful to do and requires several iterations, as said above. I think I asked about something similar some time ago, but I cant find the link.

PaulWessel commented 4 years ago

I find it way easier to specify the region in, say, km relative to the desired projection center, e.g.,

-R-1500/1500/-1000/1000+uk

Isn't that a simpler mechanism?

gmt basemap -R-1500/1500/-1000/1000+uk -JL45/35/25/45/6i -Bafg -pdf map

KristofKoch commented 4 years ago

@PaulWessel I must admit I haven't thought of this possibility. Let me toy around with it a bit and I'll report back.

anbj commented 4 years ago

Paul, if I would like to make a map created with

pscoast --MAP_FRAME_TYPE=plain -RNO -JS10/90/100/10c -B0 -png r_no -W

to be rectangular, what would be the best way? Is it possible with the method you mentioned above?

What about having an option to say to GMT, "make the map with given -R and put it in a rectangular box" (see figure, sorry for the quality):

r_no_f

Sorry to hijack the thread @KristofKoch (hope it's OK). I think we want the same thing.

KristofKoch commented 4 years ago

I found the recommendation by @PaulWessel to specify the distances from projection center quite useful but it doesn't always fit my bill if I have to fit some fixes in the extreme corners onto the plot. Then it is back to trial and error. There the idea of @anbj comes in very handy (thank you for your input) and is even simpler than my ideas presented in the original post.

Is there a way of tweaking the projection parameters to get a plot with specified dimensions when using Pauls idea? Width is no problem but height is. Say I want to fit a map on a page of DIN A4 paper (210 mm x 297 mm). I was able to achieve this with the current -R...r behavior by calculating the upper right coordinate. Am I missing something?

Additionally I found the @earth_relief_ grids to behave differently when comparing @earth_relief_01s and @earth_relief_03s to the rest. 01s and 03s behave like -JL... -R...r with missing data on western and northern edges while all other detail levels plot the data right up to the edge.

gmt begin relief_selection
  gmt grdimage @earth_relief_03s \
    -JL6/57/56/58/20c \
    -R-150/150/-150/150+uk \
    -Cdem2 -I+d -Bag
gmt end show

relief_selection

anbj commented 4 years ago

Any news about this?

stale[bot] commented 4 years ago

This issue has been automatically marked as stale because it has not had activity in the last 90 days. It will be closed if no further activity occurs within 7 days. Thank you for your contributions.

anbj commented 3 years ago

Ping?

stale[bot] commented 3 years ago

This issue has been automatically marked as stale because it has not had activity in the last 90 days. It will be closed if no further activity occurs within 7 days. Thank you for your contributions.

anbj commented 3 years ago

This is a better example of what I meant, this time made with professional software.

Enter normal -R and get your pie slice. Append e.g. f and make gmt fill in the blanks so you get a square/rectangular map.

ps=r-selection.ps

gmt set MAP_FRAME_TYPE plain

gmt pscoast -R-20/70/60/90 -JS-20/90/6c -W -B0 -K -B+t"-R-20/70/60/90" > $ps

gmt pscoast -X8c -R-20/60/70/60r -J -W -O -K -Slightblue -Ggray -B0 -B+t"-R-20/70/60/90f" >> $ps
gmt pscoast -R-20/70/60/90 -J -W -B0 -O >> $ps

gmt psconvert -Tg -A -Z $ps

rm gmt.*

r-selection

anbj commented 2 years ago

Ping

maxrjones commented 2 years ago

Following up from https://github.com/GenericMappingTools/gmt/pull/6474#issuecomment-1077506566, do you think an example of that mapproject use-case in the new gmt-examples site would be sufficient for this issue?

anbj commented 2 years ago

For me, yes, to a very large degree. But I can't answer on behalf of @KristofKoch!

anbj commented 2 years ago

Just for reference; what I wanted was implemented in #6669.

PaulWessel commented 1 year ago

@KristofKoch and @anbj , please revisit and see if #7657 provides a workable solution to some of these.

GenericMappingTools / gmt

A smarter selection of the region when using -R...r with non-cartesian projections #3099