Note: This tutorial has been updated in October 2019 to reflect the fact that the GDS has been discontinued. The general principles of this tutorial have not changed, but now you need to download a small data file to your local computer to follow the tutorial.
Since GrADS cannot easily perform high level data analysis techniques, (barring a few complicated workarounds) it is useful to save data so it can be used with more advanced software. While it is often logical to skip GrADS altogether, it can occasionally be easier to use GrADS as an intermediary step since it's a very intuitive scripting language for gridded datasets.
This tutorial will teach you how to save data using GrADS using the "write" function by guiding you through an example that will save latitude/longitude and monthly average surface temperature data from NCEP reanalysis into a .txt file.
Okay, so as always we will need to start by opening at data file. For this tutorial, I am using the monthly mean 2 meter air-temperature reanalysis product. To follow this tutorial, you will need to download the "air.mon.mean.nc" file from the ESRL PSD NCEP 1 Reanalysis website. Once you download the file, and save it to your data directory, open the file in GrADS using "sdfopen"
'sdfopen /home/your_datapath/air.mon.mean.nc'
Again, you are welcome to query the file to look at some of the metadata, but it's not required that you do so. It is however good practice to plot the data and make sure it looks resonable before proceeding. This dataset spans a long range of time, and for this tutorial, I am going to randomly choose February 2000 as the time value. To specify the date instead of the timestep when picking a time, use the 'set time' (as opposed to 'set t') command. Note that the time format must be %hhZ%dd%MMM%yyyy as in the example below.
'set time 00Z01FEB2000'
'set gxout shaded'
'd air'
Now that we have an idea what the data looks like on a map, we can move on to saving it into a .txt file. Saving data from GrADS into an outside file is actually pretty straightforward. The hard part is figuring out how to organize each piece of data exactly how you want it.
The first thing to do is to determine what data we would like to save, and how to organize it. I find it a good practice to include the x/y grid indices when saving data into a text file. So, the final text output will have five data columns, the x/y grid point indicies, the lat/lon value at each grid point, and the actual temperature data. It's also good practice to include a header line at the top of the text file to specific which column is which. So, the first step is to open the file and write header. This is done by using the "write" function in GrADS.
write('NCEP.txt', 'X Y LON LAT TEMP')
Once this command is executed in GrADS the file NCEP.txt will be created in your folder with the text "X Y LON LAT TEMP" on top.
Now that the data file has be started, it's time to start filling it up with data. To do this, you need to set the GrADS output to write our data out to the screen instead of plot it on a map. This is done by setting 'gxout' to print. In addition to that, you need to specify how GrADS should format our printed data. This is done with the 'prnopts' command using C formating syntax.
'set gxout print'
'set prnopts %6.2f 1 1' ;*%6.2f: c format, 1: values to plot on each line, ;space between values
Now that GrADS is set to print to the screen, and the print formatting is set, you are now ready to save the data. The way to do this is to loop through all the grid-points and save the data to the file one grid-point at a time. In this example, the nested loop to do this will first run across the x (longitude) direction, then through y (latitude), saving the information into rows as we go. To begin, we will first need to set limits on our domain so that the loops know when to stop. Do this by the use of the 'query' function
'q dims'
xline=sublin(result,2) ;* 2nd line
yline=sublin(result,3) ;* 3rd line
xmax=subwrd(xline,13) ;*13th word on xline
ymax=subwrd(yline,13) ;*13th word on yline
say 'X grid-points: 'xmax
say 'Y grid-points: 'ymax
This code will print out the following to the screen:
X grid-points: 193
Y grid-points: 94
Now that domain dimensions are known, you can initialize the loops. There is very little code within the nested loop. First the 'display' command (which because of the 'gxout' and 'prnopts' configuration described above) is used to print out the temperature data at each grid-point to the "result" variable. To save the data to the file, the only thingsleft to do are a) pull the temperature from the results data, b) pull the lat/lon values using the 'q dims' command, and c) write the data to a new row in the NCEP.txt file. The code is as follows:
y=1
while(y<=ymax)
x=1
while(x<=xmax)
'set x 'x
'set y 'y
'd air'
* NOTE: It may be useful to test this to find out where the data is contained with in the result
* It just so happens that in this case, the data is the 1st word of the 2nd line, this is not always true
tmp=sublin(result,2)
tmp=subwrd(tmp,1)
*Get Lat/Lon Data
'q dims'
lons=sublin(result,2)
lats=sublin(result,3)
lon=subwrd(lons,6)
lat=subwrd(lats,6)
*Save data to file
*Note the "append", so to add to the file instead of overwriting it
write('NCEP.txt', x' 'y' 'lon' 'lat' 'tmp, append)
x=x+1
endwhile
y=y+1
endwhile
That's actually all there is to it. This script will take a few moments to loop through all of the grid-points, and as it does so you will likely see the words "unknown command 0" printed a bunch of times on your screen. I have often found it helpful to include a "say" command inside the loop, so to know where I am in the saving process. When this script finishes you will have a lot of data saved into a .txt file that can be read into any plotting or statistical program. For example, to confirm GrADS saved the data correctly to text, I created the image below by loading NCEP.txt into a Python array and plotted it on a map using Matplotlib and Cartopy.
As always, I included a downloadable script for this tutorial: Download this script
Note: The downloadable script has not been updated since it was originally created, and will therefore not work when trying to access the GDS.
Since GrADS cannot easily perform high level data analysis techniques, (barring a few complicated workarounds) it is useful to save data so it can be used with more advanced software. While it is often logical to skip GrADS altogether, it can occasionally be easier to use GrADS as an intermediary step since it's a very intuitive scripting language for gridded datasets.
This tutorial will teach you how to save data using GrADS using the "write" function by guiding you through an example that will save latitude/longitude and monthly average surface temperature data from NCEP reanalysis into a .txt file.
Okay, so as always we will need to start by opening at data file. For this tutorial, I am using the monthly mean 2 meter air-temperature reanalysis product. To follow this tutorial, you will need to download the "air.mon.mean.nc" file from the ESRL PSD NCEP 1 Reanalysis website. Once you download the file, and save it to your data directory, open the file in GrADS using "sdfopen"
'sdfopen /home/your_datapath/air.mon.mean.nc'
Again, you are welcome to query the file to look at some of the metadata, but it's not required that you do so. It is however good practice to plot the data and make sure it looks resonable before proceeding. This dataset spans a long range of time, and for this tutorial, I am going to randomly choose February 2000 as the time value. To specify the date instead of the timestep when picking a time, use the 'set time' (as opposed to 'set t') command. Note that the time format must be %hhZ%dd%MMM%yyyy as in the example below.
'set time 00Z01FEB2000'
'set gxout shaded'
'd air'
 |
| Mean 2m air temperature from Feb 2000 |
Now that we have an idea what the data looks like on a map, we can move on to saving it into a .txt file. Saving data from GrADS into an outside file is actually pretty straightforward. The hard part is figuring out how to organize each piece of data exactly how you want it.
The first thing to do is to determine what data we would like to save, and how to organize it. I find it a good practice to include the x/y grid indices when saving data into a text file. So, the final text output will have five data columns, the x/y grid point indicies, the lat/lon value at each grid point, and the actual temperature data. It's also good practice to include a header line at the top of the text file to specific which column is which. So, the first step is to open the file and write header. This is done by using the "write" function in GrADS.
write('NCEP.txt', 'X Y LON LAT TEMP')
Once this command is executed in GrADS the file NCEP.txt will be created in your folder with the text "X Y LON LAT TEMP" on top.
Now that the data file has be started, it's time to start filling it up with data. To do this, you need to set the GrADS output to write our data out to the screen instead of plot it on a map. This is done by setting 'gxout' to print. In addition to that, you need to specify how GrADS should format our printed data. This is done with the 'prnopts' command using C formating syntax.
'set gxout print'
'set prnopts %6.2f 1 1' ;*%6.2f: c format, 1: values to plot on each line, ;space between values
Now that GrADS is set to print to the screen, and the print formatting is set, you are now ready to save the data. The way to do this is to loop through all the grid-points and save the data to the file one grid-point at a time. In this example, the nested loop to do this will first run across the x (longitude) direction, then through y (latitude), saving the information into rows as we go. To begin, we will first need to set limits on our domain so that the loops know when to stop. Do this by the use of the 'query' function
'q dims'
xline=sublin(result,2) ;* 2nd line
yline=sublin(result,3) ;* 3rd line
xmax=subwrd(xline,13) ;*13th word on xline
ymax=subwrd(yline,13) ;*13th word on yline
say 'X grid-points: 'xmax
say 'Y grid-points: 'ymax
This code will print out the following to the screen:
X grid-points: 193
Y grid-points: 94
Now that domain dimensions are known, you can initialize the loops. There is very little code within the nested loop. First the 'display' command (which because of the 'gxout' and 'prnopts' configuration described above) is used to print out the temperature data at each grid-point to the "result" variable. To save the data to the file, the only thingsleft to do are a) pull the temperature from the results data, b) pull the lat/lon values using the 'q dims' command, and c) write the data to a new row in the NCEP.txt file. The code is as follows:
y=1
while(y<=ymax)
x=1
while(x<=xmax)
'set x 'x
'set y 'y
'd air'
* NOTE: It may be useful to test this to find out where the data is contained with in the result
* It just so happens that in this case, the data is the 1st word of the 2nd line, this is not always true
tmp=sublin(result,2)
tmp=subwrd(tmp,1)
*Get Lat/Lon Data
'q dims'
lons=sublin(result,2)
lats=sublin(result,3)
lon=subwrd(lons,6)
lat=subwrd(lats,6)
*Save data to file
*Note the "append", so to add to the file instead of overwriting it
write('NCEP.txt', x' 'y' 'lon' 'lat' 'tmp, append)
x=x+1
endwhile
y=y+1
endwhile
That's actually all there is to it. This script will take a few moments to loop through all of the grid-points, and as it does so you will likely see the words "unknown command 0" printed a bunch of times on your screen. I have often found it helpful to include a "say" command inside the loop, so to know where I am in the saving process. When this script finishes you will have a lot of data saved into a .txt file that can be read into any plotting or statistical program. For example, to confirm GrADS saved the data correctly to text, I created the image below by loading NCEP.txt into a Python array and plotted it on a map using Matplotlib and Cartopy.
 |
| NCEP.txt data plotted using Cartopy |
As always, I included a downloadable script for this tutorial: Download this script
Note: The downloadable script has not been updated since it was originally created, and will therefore not work when trying to access the GDS.
