Hi I am using GMT 6.0 version in my windows system. I tried to plot Median Stacked profile in bathymetry area as can be seen in the figure but I am trying to plot stations depth as well on stacked line but no success. Can anyone suggest how can I approach this? I tried to use gmt project as well but somewhere I am making mistake. Can anyone please help me? I am ready to share data files.
This below is the script I am using but no success so far! The same example is given on gmt website as example 33.
set int=C:\programs\gmt6\GMT_Course\gridfiles\GEBCO_2019.nc
::set grd=C:\programs\gmt6\GMT_Course\gridfiles\etopo2.nc
::set range=127.3/131.3/27/30.4
set range=126/132/26/31
gmt begin SFM_Profiles jpg A+m0.5c -E600
::gmt makecpt -Cwiki-2.0.cpt -T-6000/4000 -Z
gmt grdcut @earth_relief_15s -G%int% -R%range% -V
gmt grdimage %int% -R%range% -JM18c -I+a300 -Cglobe -V
gmt coast -Df -Bxa1f0.5 -Bya1f0.5 -BWSEN+t"Bathymetry Depth" -W0.8p --FONT_TITLE=35,31,black -V
::gmt basemap -Lx1c/19c+w50k+l+f+jBL --FONT_TITLE=30,25,navy -V
gmt colorbar -C%cpt% -Bxa2000f1000+l"elevation" -By+l"m" -Dx8c/-1.8c/9c/0.5ch -I0.8 -V
:: Select two points along the ridge
echo 128 27 > ridge.txt
echo 130.5 30 >> ridge.txt
gmt plot -Sc0.1i -Gblue ridge.txt
gmt plot ridge.txt -W1p,blue
gmt grdtrack ridge.txt -G%int% -C360k/2k/10k+v -Sm+sstack.txt > table.txt
gmt plot -W0.5p table.txt
gmt plot OT_Stations.txt -i1,0 -Sa17p -W0.5p -Gyellow
gawk "BEGIN {FS=\",\"} {print $2, $1, $3}" OT_loc.csv | gmt text -F+f17p,Times-Roman,red,thick+jBL -Dj0.4p -N
gmt convert stack.txt -o0,5 > env.txt
gmt convert stack.txt -o0,6 -I -T >> env.txt
::gmt plot env.txt -R-400/400/-3250/250 -JX6i/2i -W1p,red -Bxa100f10+ukm -Bya1500f250+lElevation+um -BWnSe+tProfiles -X8i -Y3i -V
gmt plot -R-230/230/-8500/1000 -Bxafg1000+l"Distance from middle section (km)" -Byaf+l"Elevation (m)" -BWSne+t"MEDIAN STACKED PROFILE" -JX9i/4i -Glightgray env.txt -X8.5i -Y2i
::gmt plot env.txt -R-200/200/-3500/-2000 -JX6i/3i -W1p,gold -V
gmt plot -W3p stack.txt
::gmt project OT_Stations_depth.txt -Cenv.txt -E$(awk '{print $1, $2}' ridge.txt) -Fpz -V | gmt plot -Sa0.15i -Gblack -W0.5p
::gmt project OT_Stations_depth.txt -Cenv.txt -Q -Fpz -V | gmt plot -Sa0.15i -Gblack -W0.5p
::gmt plot OT_Stations_depth.txt -i1,0 -Sa0.15i -Gblack -W0.5p
::gmt project stack.txt -C128/27 -Fpz | gmt plot -W3p,black
::gmt project OT_Stations_depth.txt -C128/27 -Fpz | gmt plot -Sa0.15i -Gblack -W0.5p
::gmt plot OT_Stations_depth.txt -i1,0,3 -Sa0.15i -Gblack -W0.5p
::gmt legend -DjBL+w2i/2.3i+l2+o0.25c/20c -C0.15i/0.1i -F+p+g255/228/225+p+r
gmt end show
del *.info *.conf *.int *.history *.grd
pause
Hi @lalit.seismo Welcome to the GMT forum!
I think that you should use project as you did here:
gmt project stack.txt -C128/27 -Fpz | gmt plot -W3p,black
Probably you have to select the columns (with -i) that you want to plot.
I suggest to take a look at the output of gmt project. You could share here if you want.
Hi @Esteban82
Thank you for your reply.
Actually stack.txt file has multiple columns. I am sharing here please take a look. I am confused which column will be suitable because mean is calculated using these. stack.txt (17.9 KB)
Sorry, I got the wrong line.I mean the command with OT_Stantions_depth.txt
::gmt project OT_Stations_depth.txt -C128/27 -Fpz | gmt plot -Sa0.15i -Gblack -W0.5p
Think you use the ridge coordinates to set -C and -E in project, then project your data pints and use -Fqz where q is the cross distance.
And don’t be so hard on your self that you are using 6.0. Surely you can use 6.4 and avoid a bunch of bugs and enjoy new features?
Hi I tried using -i1,0 as file has first column lat second long and then depth in 3rd:
gmt project OT_Stations_depth.txt -C128/27 -Fqz | gmt plot OT_Stations.txt -Sa0.15i -i1,0 -Gblack -W0.5p
but it has plotted star randomly. Please see the figure.
And also this file looks like thisOT_Stations_depth.txt (541 Bytes)
Please suggest something. I am also trying to plot 0 to -600m depth stacked (black line) make pink or any other color so that it can be seen clearly.
Hi, Thanks for your suggestion. Actually -C option I am using as -C128/27 but unable to understand how to use -E.
Please suggest where I am making mistake:
gmt project OT_Stations_depth.txt -C128/27 -Fqz | gmt plot OT_Stations.txt -Sa0.15i -i1,0 -Gblack -W0.5p
txt file also have attached previously which has stations lat,long and depth.
Not so easy. This gives the q distances
gmt project -E128/27 -C130.5/30 -Fqz -Q OT_Stations_depth.txt
but those distances are not zero on the ridge line (like your stack). So maybe run the same project command on the C or E point and then get its q value and subtract from the star q values and that will work?
1 Like
@pwessel is this output wrong? I was expecting the depths value first.
echo 28 127 1127 | gmt project -C128/27 -E130.5/30 -fg -: -Fzq
1.30227717218 1127
I add -; because the input is Lat Long. But, it seems that -; also flip the output.
echo 28.6239 127.7732 1127 | gmt project -C128/27 -E130.5/30 -Q -Fqz
-685.31053591 1127
@lalit.seismo does this map make sense to you?
I use -: in project becasue your data is Lat Long (and GMT expects Long Lat). In plot, I multiply the depth by -(+s-1) because in your data it was positive.
Here is the script:
gmt begin test png
gmt basemap -R-230/230/-8500/1000 -Bxafg1000+l"Distance from middle section (km)" -Byaf+l"Elevation (m)" -BWSne+t"MEDIAN STACKED PROFILE" -JX9i/4i
gmt project OT_Stations_depth.txt -C128/27 -E130.5/30 -Q -Fqz | gmt plot temp -Sa0.15i -Gblack -W0.5p -i1,0+s-1
gmt end
1 Like
z is always last since may contain trailing text. Probably not well documented…
But here I got z first (it is the same example but I manually put the long first)
echo 127.7732 28.6239 1127 | gmt project -C128/27 -E130.5/30 -Q -Fzq
1127 123.183378965
I think it is very clear:
If output format is ASCII then z also includes any trailing text (which is placed at the end of the record regardless of the order of z in flags ).
Yes, I spoke before looking. Just committed a minor PR so that the usage message also mentions this.
Hi @Esteban82, Thank you for explaining it. Yes, plot must look like this.
gmt project OT_Stations_depth.txt -C128/27 -E130.5/30 -Q -: -Fqz | gmt plot -Sa0.15i -Gblack -W0.5p -i1,0+s-1
Only thing i changed for right depth reversed the -C and -E long/lat and removed temp
Though there is one outlier station. I did not understand how come it’s going outside the range? Can we write stations name or number as well here?
One more question i’ve, can we colored a particular section of stacked line? like from the depth o to -600 m by any how? For this i also searched a lot but no success so far!
If you want to label the stars you have to supply what you want the label to be. Probably pipe your gmt project stations via text using the +r modifier to -F to select running record number as the text.
As for color. If you mean color the solid stacked line it depends if you just want it to have a constant color between two z-values. If so, maybe run the stacked curve through sample1d to get intermediate points, then through gmt select to only output the section you want and then plot that with your color.
