Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I asked a question before that nobody answered...but I kind of figured it out myself. But I still remain with some things that I just don't seem to get right. Can anybody please find the time to help me? It really seems like something trivial... Thanks!

-

I want to transform a geographicly gridded latlon data-set into a polar stereographic projection with gdalwarp. For this I need to add the correct metadata to my gridded data-set, and I can't seem to get this right.

It is regularly spaced with 0.04 degrees from 0:360 longitude and -86:-60 Latitude (Antarctica). Only thing I have is the grid, so it comes without a header. I however have added the current header to the dataset:

ncols        9001
nrows        651
xllcorner    0
yllcorner    -86.000000000000
dx           0.04
dy           0.04

The first problem I get is that when I gdalinfo it I get:

Corner Coordinates:
Upper Left  (   0.0000000, -59.9600000) 
Lower Left  (   0.0000000, -86.0000000) 
Upper Right (     360.040,     -59.960) 
Lower Right (     360.040,     -86.000) 
Center      ( 180.0200000, -72.9800000) 

Which seems wrong (as it is 0.04 off). Possibly I have to use a gridsize of 9000/650 instead of 9001/651, but also this is confusing to me as I then omit a line. Do I set the header correctly? Other than defining the grid and the coordinates like this, what do I need to do with the projection data? Need I add -s_srs WGS84 (this is what I do) or something like -s_srs EPSG:3031 (South Polar Stereographic)? This gives me the following:

Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
    SPHEROID["WGS 84",6378137,298.257223563,
        AUTHORITY["EPSG","7030"]],
    AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]]
Origin = (0.000000000000000,-59.960000000000001)
Pixel Size = (0.040000000000000,-0.040000000000000)

I want to finally project it onto the following system:

PROJCS["unnamed",
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]],
PROJECTION["Polar_Stereographic"],
PARAMETER["latitude_of_origin",-71],
PARAMETER["central_meridian",0],
PARAMETER["scale_factor",1],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["metre",1,
    AUTHORITY["EPSG","9001"]]]
 Origin = (-2800500.000000000000000,2800500.000000000000000)
Pixel Size = (1000.000000000000000,-1000.000000000000000)
Metadata:
AREA_OR_POINT=Area
Image Structure Metadata:
 INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (-2800500.000, 2800500.000) ( 45d 0' 0.00"W, 54d40' 3.36"S)
Lower Left  (-2800500.000,-2800500.000) (135d 0' 0.00"W, 54d40' 3.36"S)
Upper Right ( 2800500.000, 2800500.000) ( 45d 0' 0.00"E, 54d40' 3.36"S)
Lower Right ( 2800500.000,-2800500.000) (135d 0' 0.00"E, 54d40' 3.36"S)
Center      (   0.0000000,   0.0000000) (  0d 0' 0.01"E, 90d 0' 0.00"S)

Now there are 2 things occuring with my results when I warp it like this:

gdalwarp -s_srs WGS84 -r cubic data1.txt data2.tif

(i use either an existing data-set (I want to remap so...) or set the -te and -tr parameters, and get the same results....)

  1. I get a warped file which fails to do the entire file (there is probably something going wrong at the 0/360 longitude transition. As expected it omits the value that falls at -360.04
  2. I get only half of the file!! Only the right side of the projection is done...

So my questions are:

  1. How do i correctly set the correct projection parameters and how do I deal with polar stereographic data? Why can't I properly set the latitude and longitude range?
  2. What in earth is going wrong when it only projects the right side?

Thanks for the help!

share|improve this question
    
Nobody help me? I am totally stuck :( –  jmvanwessem May 8 '12 at 13:37

1 Answer 1

I can't answer all your questions, but I have noticed a couple of things that might help. First, your apparent off-by-0.04 is because by default, GDAL assumes the centre of a pixel represents the coordinate, not the corner. So by shifting your xllcorner and yllcorner by -0.02 will end up with the right coverage.

Secondly, most CRSs that I'm aware of have a range between -180 and +180 degrees longitude, whereas your bounds are 0 - 360. Try setting your xllcorner value to -180.02 (assuming you're using Greenwich as your prime meridian).

share|improve this answer

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.