In RPN files the latitudes and longitudes of all grid points are coded. Depending on the grid type they are coded in the parameters IG1-IG4 (for A, B, G, and L grids) or in the tictacs, '^^', '>>', (Z grids) or in the toctoc, '^>' (U grids, Yin-Yang).

Shell

Obtain longitudes and latitudes in shell

infile=...  # RPN input file with one record and tictac's (if needed for the grid type - if xrec can display the field you have everything you need)
outfile=... # RPN output file (including path or at least './')

# Create field with constant value 180 / PI = 57.29578
r.diag xlin   $infile  tmp  -a 0 -b 57.29578

r.diag ggtrig tmp    tmp2 -kind IDF -lon  # ==> outfile=longitude*infile
r.diag newnam tmp2   lon  -name LO
r.diag ggtrig tmp    tmp2 -kind IDF       # ==> outfile=latitude*infile
r.diag newnam tmp2   lat  -name LA

# Combine LO and LA in output file
r.diag joinup $outfile lon lat

Python 3

Obtain longitudes and latitudes in Python

#!/usr/bin/env python3

#####################################################################################
#
# Before opening Python you need to first load the Python version (in this case Python3) you want to use with:
#     module load python3/miniconda3
# and then load the module needed to read RPN files with Python with:
#     module load python3/python-rpn
#     module load python3/outils-divers
# and to get access to basic Python packages source:
#     source activate base_plus
#
#####################################################################################

# Importation of modules and library
import rpnpy.librmn.all as rmn              # Module to read RPN files
from rotated_lat_lon import RotatedLatLon     # Module to project field on native grid (created by Sasha Huziy)

# Read one record
filename = '...'                            # Name of RPN file to read
varname  = '...'                            # Name of variable to read
fid = rmn.fstopenall(filename,rmn.FST_RO)   # Open the file
rec = rmn.fstlir(fid,nomvar=varname)        # Read the full record of variable 'varname'
field = rec['d']                            # Assign 'field' to the data of 'varname'

# Read 2-D latitudes & longitudes - if needed
mygrid = rmn.readGrid(fid,rec)              # Get the grid information for the (LAM) Grid -- Reads the tictac's
latlondict = rmn.gdll(mygrid)               # Create 2-D lat and lon fields from the grid information
lat = latlondict['lat']                     # Assign 'lat' to 2-D latitude field
lon = latlondict['lon']                     # Assign 'lon' to 2-D longitude field

# Close RPN input file
rmn.fstcloseall(fid)                        # Close the RPN file








  • Aucune étiquette
Écrire un commentaire...