How extend profile?

Q&A
1

I have a profile (see below) that lengths ca. 15 km. and I want extend the profile that same distance in both direction. I.e. I want to add two extra points (at the beginning and the end following its azimuth) at a distance of 15 km (so the new profile would be 45 km).

I think it should be easy but I can’t figure out how to do it.

#  Long Lat
-66.04971205	-55.14285297
-66.01518493	-55.11524701
-65.97886793	-55.08903702
-65.94340697	-55.06219498
-65.90790703	-55.03662696
2

Not sure about easy (i.e., scripting) but I would do this:

  1. Run gmt mapproject -Af (or -Ab, cannot recall) on your data file to get the azimuths you need at either end.
  2. Use gmt vector -Ttazim/15k to get your two points.
3

Ok, thanks.

And then I should make something like this to have a new profile, right?

gmt vector **-Tt** *azim*/15k > New_profile   # new start point
gmt convert data.txt -o0,1 >> New_profile
gmt vector **-Tt** *azim*/15k >> New_profile  # new end point
4

Yes, something like that if you need to script this for multiple use. So azim1 and azim2 then since they are likely different. And I am guessing -Ab for the first point (to get azimuth backward) and -Af for the last point. Maybe mapproject -Ab file | gmt convert -Ef to get first azim1 and map project -Af file | gmt convert -El to get the last azim.

5

And the solution is

# az2 is the forward azim at point2
julia> d, az1, az2 = invgeod([-66.04971205 -55.14285297], [-65.90790703 -55.03662696])
([14893.888506208705], [37.49803389341683], [37.381746701225865])

# Extend point 2 15 km along azim2
julia> geod([-65.90790703 -55.03662696], 37.3817, 15000)
([-65.76584878309455 -54.92947611687853], 37.26535673380841)


# Same story but revert point1 & point2
julia> d, az1, az2 = invgeod([-65.90790703 -55.03662696], [-66.04971205 -55.14285297])
([14893.888506208705], [-142.61825329877414], [-142.50196610658315])

julia> geod([-66.04971205 -55.14285297], -142.50197, 15000)
([-66.19329013066837 -55.24966629463536], -142.38407625074115)

or

-66.19329013066837 -55.24966629463536
-66.04971205 -55.14285297
-65.90790703 -55.03662696
-65.76584878309455 -54.92947611687853

extended=[-66.19329013 -55.249666294635; -66.04971205 -55.14285297; -65.90790703 -55.03662696; -65.765848783 -54.929476117];
imshow(extended, lc=:blue, marker=:circ, mc=:red)

GMTjl_tmp
GMTjl_tmp1809×1211 29.8 KB

6

OK, but if it is an important thing to do then I could implement gmt spatial -L15k and that would be even shorter and accessible to all wrappers.

7

Yes, thanks I will definitively used it.

8

For the record.
A new feature was added in GMT 6.4: gmt spatial -W.