Re: BEDMAP with GMT

Hi

The original script works with bm2_draft_bed.txt as it is an ascii file with x, y, and z values.

However, bedmap2_bed.txt is an ascii grid file.

Try the following to convert the ascii grid file to a GMT format grid:
xyz2grd bedmap2_bed.txt -H6 -F -ZTLA -N-9999 -R-3333500/3333500/-3333500/3333500 -I1000 -GBEDMAP2_GMT.grd

This should make you a GMT grid.  You can then use grd2xyz to dump out an x, y and z ascii file just like bm2_draft_bed.txt, and run it through the original script, or follow Pauls suggestion below to avoid gridding over the pole.

Tom

From: Paul Wessel [mailto:[email protected]]
Sent: 25 June 2014 09:25
To: [email protected]
Subject: Re: [GMT-HELP] BEDMAP with GMT

I am surprised this ever worked.
Using surface to grid geographic data over a pole is a horrid idea since surface is best for Cartesian data or latitude-banded data.
If you absolutely need to convert this polar stereographic grid to degrees then I would do something like this:
1. Create a blank grid with desired resolution for -R-180/180/-90/-60
2. Dump via grd2xyz and project coordinates with the -Js projection
3. Sample the stereotraphic bedmap grid at the projected coordinates
4. paste the result with the grd2xyz output and send via xyz2grd to make the desired grid

Or something similar.
-p


On Jun 24, 2014, at 11:39 AM, Sumant Jha <[email protected]<mailto:[email protected]>> wrote:


Hi All,

Raising this topic again.
For some reason, the following solution does not works with the final release of bedmap2 which is available at https://secure.antarctica.ac.uk/data/bedmap2/
I have the old bm2_draft_bed.txt, and its output comes exactly as what Tom has attached. However, when using bedmap2_bed.txt (which I think, is the later version of bm2_draft_bed.txt), the output comes as attached.

Can anyone help me understand what has changed since the draft version to final version of BEDMAP2? Or is it that I am doing something wrong in my code? I am using same code that Tom has given in his email with minor changes like replacing the bm2_draft_bed.txt by bedmap2_bed.txt, using the bedmap2_bed.txt as input, and changing the -R to -180/180/-90/-60.

The code for bedmap2_bed and its output are attached.

#!/bin/sh

#  bedmap_ant.sh
#
#
#  Created by Sumant Jha on 12/27/12.
#  Based on TOM JORDAN's script, dated Nov 21, 2012 in [GMT-HELP] BEDMAP with GMT
#  Modified on 06/24/2014
#
# script to get degree grid of BEDMAP2 data
# Note: For some reason, this script is only working on bm2_draft_bed.txt

# clean and order xyz data
echo "***********************************************"
echo " Cleaning and Ordering X Y Z data in desired sequence"
awk ' NR>1 { print $1*1, $2*1, $3*1 } ' bedmap2_bed.txt | awk ' { print $2, $3, $1 } ' > bedmap


# projecting bedmap data to degrees (default bedmap comes in UTM) NOTE -C option
echo " "
echo " "
echo "***********************************************"
echo " Projecting BEDMAP2 to Degrees"
mapproject bedmap -I -R-180/180/-90/-60 -Js0/-90/-71/1:1 -C -F -V> temp1

# Make new 1 min grid
echo " "
echo " "
echo "***********************************************"
echo " Creating 1 min grid for BEDMAP2"
surface temp1 -R-180/180/-90/-60 -I1m -V -Gbedmap_1min.grd

# make colour pallet
echo " "
echo " "
echo "***********************************************"
echo " Creating color pallett for displaying data"
grd2cpt -Z -V bedmap_1min.grd > bedmap1.cpt

gmtset PAPER_MEDIA A4


# Draw Antarctica Bedmap grid
echo " "
echo " "
echo "***********************************************"
echo "Drawing Antarctica BEDMAP grid"
grdimage bedmap_1min.grd -Cbedmap1.cpt -R-180/180/-90/-60 -Js0/-90/-71/1:50000000 -B30g300/15g15wEsN -K -V -P> BEDMAP_ant.ps
# Draw coast
pscoast -W1 -R -Js -O >> BEDMAP_ant.ps

# Create contours spaced at 200 intervals (optional)
#grdcontour Bedmap_ant.grd -J -R -W1/80/80/80 -C200 -P -S4 -O -A1000+f3+k80/80/80+s8t -G2i/10 -Djunk -V -K >> BEDMAP_Ross.ps
echo " "
echo " "
echo "***********************************************"
echo "File Generated"

#View
gs BEDMAP_ant.ps

<BEDMAP_ant.png>

On Nov 21, 2012, at 6:46 AM, Jordan, Tom A.R.M. <[email protected]<mailto:[email protected]>> wrote:


Hi All

One solution to the problem is below.

First get the XYZ data for BEDMAP2 from the same ftp site Roi mentioned (it's in a different directory).
ftp://ftp.bas.ac.uk/ptf/bm2_ftp/text_version_of_bed/

You can then re-project the data to degrees with mapproject, the trick is the -C option.

Re grid the degree data and off you go.

Alternatively you could dump the coast line from pscoast (-m option).
Project the coast line to stereographic meters, same map project but no -I.
Then simply plot the projected coast on top of the bedmap grid extracted as Roi suggested using psxy -Jx.

Hope this is useful.

Tom

# script to get degree grid of BEDMAP2 data, its slow but works
# clean and order xyz data
awk ' NR>1 { print $1*1, $2*1, $3*1 } ' bm2_draft_bed.txt | awk ' { print $2, $3, $1 } ' > temp

# project data to degrees NOTE -C option
mapproject temp -I -R0/360/-90/0 -Js0/-90/-71/1:1 -C -F > temp1

# Make new 5 min grid
surface temp1 -R0/360/-90/-60 -I5m -Gbedmap_5min.grd

# make colour pallet
grd2cpt -Z bedmap_5min.grd > temp.cpt

# draw grid
grdimage bedmap_5min.grd -Ctemp.cpt -R0/360/-90/-60 -Js0/-90/-71/1:50000000 -K > temp.ps

# draw coast
pscoast -W1 -R -Js -O >> temp.ps

# view
gv temp.ps


From: World-Wide GMT Usage and Help Mailing List [mailto:[email protected]] On Behalf Of Roi Granot
Sent: 20 November 2012 17:49
To: [email protected]<mailto:[email protected]>
Subject: Re: [GMT-HELP] BEDMAP with GMT

There is a new Bedmap compilation (Bedmap2, doi:10.5194/tcd-6-4305-2012) which you can download here: ftp://ftp.nerc-bas.ac.uk/ptf/bm2_ftp/readme%20and%20data%20files/

With gdal you can convert it:
dal_translate -of GMT bm2_draft_bed.tif bedmap2.grd

and plot:
grdimage -Jx-.000001/-.000001 -Ibedmap2.grd.int -B5000000 -R-3333500/3333500/-3333500/3333500  bedmap2.grd -Cbed.cpt -V > bedmap2.ps

Roi

Dr. Roi Granot
Department of Geological & Environmental Sciences
Ben Gurion University of the Negev, Beer Sheva
P.O.Box 653 Israel 84105
Tel: 972-8-6477-509






On Nov 20, 2012, at 3:43 PM, Marcus Badger wrote:



Hi Sumant,

I downloaded the "Grounded ice-sheet bed elevation" in Ascii format from BAS from here:

http://www.antarctica.ac.uk//bas_research/data/access/bedmap/download/<http://www.antarctica.ac.uk/bas_research/data/access/bedmap/download/>

Then used xyz2grd to convert to grid. If you already have the grid then it might be fine. What do  you get when you run

grdinfo BEDMAP_groundbed.grd

For my grid I get:
groundbed.grd: Title: groundbed.grd
groundbed.grd: Command: xyz2grd groundbed.asc -E -Ggroundbed.grd
groundbed.grd: Remark:
groundbed.grd: Pixel node registration used
groundbed.grd: Grid file format: nf (# 18) GMT netCDF format (float)  (COARDS-compliant) [DEFAULT]
groundbed.grd: x_min: -2661600 x_max: 2678400 x_inc: 5000 name: x nx: 1068
groundbed.grd: y_min: -2149967 y_max: 2195033 y_inc: 5000 name: y ny: 869
groundbed.grd: z_min: -2495.36010742 z_max: 4364.85400391 name: z
groundbed.grd: scale_factor: 1 add_offset: 0

So you can see that the x and y are not in latlon but metres relative to the south pole so in my script I just set -R to the min and max. As this is already projected data you should then just need to overlay the appropriate latlon grid using psbasemap -Js with the appropriate parameter to match the projection of the grid. The docs give "we used the Polar Stereographic projection with 71°S as the latitude of true scale and 0°E as the central meridian" although I'm not sure what parameters you therefore need for -Js

Marcus.


----------------------------------------------------------
Dr Marcus Badger
Lecturer

School of Earth Sciences
University of Bristol
Wills Memorial Building
Queen's Road
BRISTOL BS8 1RJ
U.K.

Office: G40
Tel: +44 (0) 117 95 45491 or +44 (0) 117 33 16795
Mobile: +44 (0) 774 07 96811
Email: [email protected]<mailto:[email protected]>
----------------------------------------------------------

On 16/11/2012 17:20, Sumant Jha wrote:


Hi Markus,

Pardon my ignorance on this, as I am very new to GMT. Do we need to convert the BEDMAP grid to some other format? Where are you getting the xyz data for the same?

I was using the BEDMAP grid that came with Geoware.

You are right about the -R values. I noticed the same but was not sure, how to go about converting my -R values to that of BEDMAP values. Your version is something what I am trying to achieve, and then overlay the lat/longs. Not sure how to do that.

Sumant


On Nov 15, 2012, at 2:41 AM, Marcus Badger <[email protected]<mailto:[email protected]>> wrote:

Hi Sumant,

How far have you got? ie in what way is it not coming out right? I've managed to make this into a grid using xyz2grd and plot it using

grdimage groundbed.grd -V -Cpalette.cpt -B2000000NseW -R-2661600/2678400/-2149967/2195033 -JX5i -Yc -Xc -P -K >Antarctica.ps

but the grid is not in lat-lon but from the docs page "In this system, coordinates are expressed in metres with the origin at the South Pole." (which might be where you problems start) and so I am unsure how to convert this to something you can overlay coastline and gridlines. My version attached.


Marcus.

----------------------------------------------------------
Dr Marcus Badger
Lecturer

School of Earth Sciences
University of Bristol
Wills Memorial Building
Queen's Road
BRISTOL BS8 1RJ
U.K.

Office: G40
Tel: +44 (0) 117 95 45491 or +44 (0) 117 33 16795
Mobile: +44 (0) 774 07 96811
Email: [email protected]<mailto:[email protected]>
----------------------------------------------------------

On 13/11/2012 20:34, Sumant Jha wrote:
Guys - I fear this question has been asked a lot before - but, i cannot
seem to find an answer to it. I am trying to use the BEDMAP data in GMT.
The problem it seems to me is the projection. I tried using Jx and then
Js to overlay the lat-longs - but, somehow, it did not come out right. I
am hoping if someone out here have worked with GMT and BEDMAP data - and
could may be help me solve this issue.
I am using the following code to plot it:

grdimage
/Users/sumantjha/Documents/Geoware/DATA/misc/BEDMAP_groundbed.grd -V
-C/Users/sumantjha/Documents/Work/GMT/GMT4.5.8/share/cpt/GMT_gebco.cpt
-B20NseW -R155/180/-80/-75 -Jx10i -X1.2i -Y6i -K -V -P > Antarctica.ps

pscoast -R -Js0/-90/50i/1 -W -Df -Ir/.5/114/187/243 -G30/144/255 -V -O
-K -P >> Antarctica.ps

Any help will be appreciated.

Thanks
Sumant
To unsubscribe, send the message "signoff gmt-help" to
[email protected]<mailto:[email protected]>

To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]>
<Antarctica.jpg>

To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]>


To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]>

To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]>



--
This message (and any attachments) is for the recipient only. NERC
is subject to the Freedom of Information Act 2000 and the contents
of this email and any reply you make may be disclosed by NERC unless
it is exempt from release under the Act. Any material supplied to
NERC may be stored in an electronic records management system.
To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]> <temp.jpg>

Mailing list for GMT discussions of all kinds. If you are not sure you have found a bug, discuss it here first. To formally report bugs or request features, please register and add New Issue on gmt.soest.hawaii.edu<http://gmt.soest.hawaii.edu> To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]> Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on gmt.soest.hawaii.edu<http://gmt.soest.hawaii.edu> instead.

Mailing list for GMT discussions of all kinds. If you are not sure you have found a bug, discuss it here first. To formally report bugs or request features, please register and add New Issue on gmt.soest.hawaii.edu To unsubscribe, send the message "signoff gmt-help" to [email protected]<mailto:[email protected]> Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on gmt.soest.hawaii.edu instead.

________________________________
This message (and any attachments) is for the recipient only. NERC is subject to the Freedom of Information Act 2000 and the contents of this email and any reply you make may be disclosed by NERC unless it is exempt from release under the Act. Any material supplied to NERC may be stored in an electronic records management system.

Mailing list for GMT discussions of all kinds.  If you are not sure you have found a bug, discuss it here first.
To formally report bugs or request features, please register and add New Issue on gmt.soest.hawaii.edu
To unsubscribe, send the message "signoff gmt-help" to [email protected]
Note: gmt-help will become obsolete on Sept 1, 2014 - please use forum on gmt.soest.hawaii.edu instead.