[GRASS-user] Issues when reprojecting from a GCS to a PCS

Hanlie Pretorius hanlie.pretorius at gmail.com
Thu Jul 15 09:35:11 EDT 2010


Hi Everyone,

I'm working in GRASS 6.4RC6 on Win XP and am trying to reproject TRMM
rainfall data from its native GCS (WGS84) to a local PCS. I'm using
the PCS because my other important data, specifically DEMs, are
provided in it. TRMM data come in 0.25°x0.25° squares.

My workfow is as follows:

1. Import the TRMM data into a WGS84 location using r.in.xyz. This
gives a raster with the following information:

 +----------------------------------------------------------------------------+
 | Layer:    3B42.000215.0 at PERMANENT        Date: Thu Jul 15 14:42:27 2010
 | Mapset:   PERMANENT                      Login of Creator: hanlie
 | Location: world_wgs84
 | DataBase: F:\Hanlie\grassdata
 | Title:    Raw x,y,z data binned into a raster grid by cell mean ( 3B42.000
 | Timestamp: none
 |---------------------------------------------------------------------------
 |
 |   Type of Map:  raster               Number of Categories: 255
 |   Data Type:    FCELL
 |   Rows:         6
 |   Columns:      3
 |   Total Cells:  18
 |        Projection: Latitude-Longitude
 |            N:     27:15S    S:     28:45S   Res:  0:15
 |            E:     28:45E    W:        28E   Res:  0:15
 |   Range of data:    min = 0.000000  max = 0.070614
 |
 |   Data Source:
 |    F:\Hanlie\UCT\M.Sc\Data\TRMM\2000\02_Februarie\3B42.000215.0.6.nc.lieb.
 |
 |
 |   Data Description:
 |    generated by r.in.xyz
 |
 |   Comments:
 |    r.in.xyz input="F:\Hanlie\UCT\M.Sc\Data\TRMM\2000\02_Februarie\3B42.\
 |    000215.0.6.nc.lieb.txt" output="3B42.000215.0" method="mean" type="F\
 |    CELL" x=2 y=1 z=3 zscale=1.0 percent=100
 |
 +----------------------------------------------------------------------------+

2. Set the region to this raster file and create a vector map from it
using v.in.region. This gives a vector with the following information:

 +----------------------------------------------------------------------------+
 | Layer:           lieb_trmm_region at PERMANENT
 | Mapset:          PERMANENT
 | Location:        world_wgs84
 | Database:        F:\Hanlie\grassdata
 | Title:
 | Map scale:       1:1
 | Map format:      native
 | Name of creator: hanlie
 | Organization:
 | Source date:     Thu Jul 15 14:55:13 2010
 |----------------------------------------------------------------------------|
 |   Type of Map:  vector (level: 2)
 |
 |   Number of points:       0               Number of areas:      1
 |   Number of lines:        0               Number of islands:    1
 |   Number of boundaries:   1               Number of faces:      0
 |   Number of centroids:    1               Number of kernels:    0
 |
 |   Map is 3D:              No
 |   Number of dblinks:      0
 |
 |         Projection: Lat/Lon
 |               N:            27:15S    S:            28:45S
 |               E:            28:45E    W:               28E
 |
 |   Digitization threshold: 0
 |   Comments:
 |
 +----------------------------------------------------------------------------+

3. Change to the PCS location and reproject the TRMM region using
v.proj. This gives a vector with the following information:

 +----------------------------------------------------------------------------+
 | Layer:           lieb_trmm_region at PERMANENT
 | Mapset:          PERMANENT
 | Location:        sa_tm_29deg_E
 | Database:        F:\Hanlie\grassdata
 | Title:
 | Map scale:       1:1
 | Map format:      native
 | Name of creator: hanlie
 | Organization:
 | Source date:     Thu Jul 15 14:55:13 2010
 |----------------------------------------------------------------------------|
 |   Type of Map:  vector (level: 2)
 |
 |   Number of points:       0               Number of areas:      1
 |   Number of lines:        0               Number of islands:    1
 |   Number of boundaries:   1               Number of faces:      0
 |   Number of centroids:    1               Number of kernels:    0
 |
 |   Map is 3D:              No
 |   Number of dblinks:      0
 |
 |         Projection: Transverse Mercator
 |               N: -3015356.39959204    S: -3181970.90748457
 |               E:    -24418.1519951    W:   -99037.39660484
 |
 |   Digitization threshold: 0
 |   Comments:
 |
 +----------------------------------------------------------------------------+

4. Set the region of the PCS to the TRMM region:

g.region -p vect=lieb_trmm_region at PERMANENT nsres=25 ewres=25

which gives:

projection: 99 (Transverse Mercator)
zone:       0
datum:      ** unknown (default: WGS84) **
ellipsoid:  wgs84
north:      -3015356.39959204
south:      -3181970.90748457
west:       -99037.39660484
east:       -24418.1519951
nsres:      24.99842579
ewres:      24.9980719
rows:       6665
cols:       2985
cells:      19895025


I use this resolution because I don't really know what else to use.
When I reproject the TRMM GRID into the PCS, I get polygons with
differing areas, so I can't use the TRMM grid size as resolution:

v.report map=lieb_trmm_squares at PERMANENT option=area units=kilometers
Displaying column types/names for database connection of layer 1:
cat|LabelX|LabelY|LblOffsetX|LblOffsetY|Label|gRow|gColumn|RowCol|area
1|28|-28.5|-20|-20|28||||679.062503849138
2|28|-28.25|-20|-20|28||||680.623828116383
3|28|-28|-20|-20|28||||682.172234597964
4|28|-27.75|-20|-20|28||||683.70769882544
5|28.25|-28.75|-20|-20|28.25||||677.428350142518
6|28.25|-28.5|-20|-20|28.25||||679.002141385312
7|28.25|-28.25|-20|-20|28.25||||680.563040996854
8|28.25|-28|-20|-20|28.25||||682.111024322052
9|28.25|-27.75|-20|-20|28.25||||683.646066944101
10|28.25|-27.5|-20|-20|28.25||||685.168144684034
11|28.5|-28.5|-20|-20|28.5||||678.961903132767
12|28.5|-28.25|-20|-20|28.5||||680.522519692921
13|28.5|-28|-20|-20|28.5||||682.070220966439
14|28.5|-27.75|-20|-20|28.5||||683.604982570918
15|28.5|-27.5|-20|-20|28.5||||685.126780361735


5. Reproject the TRMM rainfall from the GCS to the PCS, which gives no
data in the file. The extents also seem to be wrong:

 +----------------------------------------------------------------------------+
 | Layer:    3B42.000215.0 at PERMANENT        Date: Thu Jul 15 15:01:05 2010
 | Mapset:   PERMANENT                      Login of Creator: hanlie
 | Location: sa_tm_29deg_E
 | DataBase: F:\Hanlie\grassdata
 | Title:     ( 3B42.000215.0 )
 | Timestamp: none
 |----------------------------------------------------------------------------|
 |
 |   Type of Map:  raster               Number of Categories: 255
 |   Data Type:    FCELL
 |   Rows:         5552
 |   Columns:      1996
 |   Total Cells:  11081792
 |        Projection: Transverse Mercator
 |            N: -3029230.52590567    S: -3168021.78589357   Res: 24.99842579
 |            E: -36667.20722417    W: -86563.35872871   Res: 24.9980719
 |   Range of data:    min = nan  max = nan
 |
 |   Data Description:
 |    generated by r.proj
 |
 |   Comments:
 |    r.proj input="3B42.000215.0" location="world_wgs84" mapset="PERMANEN\
 |    T" method="cubic" resolution=25.0
 |
 +----------------------------------------------------------------------------+


Can someone perhaps see where I'm going wrong? I have hundreds of
these files to import, so I need to get the workflow right and then
automate the task.

Thanks
Hanlie


More information about the grass-user mailing list