ericpante/marmap

lc.dist() paths crossing continents

Opened this issue · 1 comments

Hi,
I'm having some issues to plot the least-cost paths between stations in a global scale.

I have sampling stations globally distributed with Long ranging from -179.52 to 179.14 and Lat from -40.55 to 35.39.
First, I downloaded the NOAA bathymetry data using getNOAA.bathy() as followed:

world<-getNOAA.bathy(-180,180,-90,90,res=15,keep=TRUE) #also tried with 'res=4' and got similar results as described bellow

Second, I created a transition object to avoid continents:

trans1<-trans.mat(world,min.depth=-10)

Then I calculated the least-cost distances and paths between stations using lc.dist():

dist1 <-lc.dist(trans1, geo_ref_stations, res ="dist")

out1 <-lc.dist(trans1, geo_ref_stations, res ="path")

The resulting distance matrix seems to be OK (i.e.: correct distance among anti-meridian sites), but when I plot the paths to a global map using this code, they look messy (crossing continents) like the map in example bellow:

plot(world,xlim=c(-180,180),ylim=c(-90,90),deep=c(-5000,-200,0),shallow=c(-200,0,0),col=c("grey","blue","black"),step=c(1000,200,1),lty=c(1,1,1),lwd=c(0.6,0.6,1.2),draw=c(FALSE,FALSE,FALSE))

list.apply(out1, lines,col="green")->dummy
points(geo_ref_stations,pch=21,col="navy",bg=col2alpha("red",.9),cex=1)

map_plot1_noaa

I tried to solve it by uploading the global GEBCO bathymetric dataset (downloaded from their website) using readGEBCO.bathy() as followed:

world<-readGEBCO.bathy("./GEBCO_2020.nc",resolution = 15)

Then I repeated the codes above to calculate the least.cost distances and plot the paths to a global map. Now, the paths on the map look nice (see bellow), but the resulting distances are NOT precise since the function calculates the distances as if the paths had to go all around the globe to connect antimeridian sites which are actually geographically closed.

map_plot2_GEBCO

Any solutions for these issues? I wonder if there is a way to either correct the paths using the NOAA data OR simply consider the antimeridian paths with the global GEBCO database.

Cheers, Pedro

Dear Pedro,

I'm really sorry I wasn't able to come back to you sooner.

I believe the problem comes from the way longitudes are encoded. I know that when using ETOPO1, longitudes range between -180 and 180. I suspect that when you used the GEBCO data, the longitudes are actually encoded as 0-360.

I can thus confirm that the distances you computed using ETOPO1 are correct. The strange paths crossing the continents are due to the way lines are drawn on a planisphere. To simplify, when the shortest path between 2 locations "A" and "B" crosses the antemeridian, the path should be split in two parts, one to the East of the antemeridian, and one to the West of the antemeridian (paths A-X and B-X respectively, see attached file). However, since the world map is represented as a planisphere centered on the Greenwich meridian (instead of a continuous sphere), the path is drawn on the wrong side of the globe, hence crossing the Greenwich meridian instead of the antemeridian (path A-B, again, see attached file).

Using GEBCO data produces wrong distances because of the way lc.dist() work: it expects longitudes in the range -180;180.

Unfortunately, changing these behaviours is far from trivial and I'm afraid I won't be able to deal with it for quite some time. The good news is that the distances computed by lc.dist() are indeed correct when using ETOPO1, and can thus be safely used for other computations or publication. The possibility to draw the paths was kind of gimmicky and tested mostly at small scales.

If marmap users want to help, they are more than welcome to contribute to speed things up!

Rplot01