Ken Mankoff wrote > Dear List, > > Most of the time when I work with data in GRASS it is a provided on a > standard GIS grid (EPSG code) and I can work with it easily. I create a > location using the "-c" option when I start GRASS, or I import it into an > existing mapset with r.import. > > I'm now working with data and all I know is this: > > +proj=stere +lat_0=90 +lat_ts=71 +lon_0=-39 +k=1 +x_0=0 +y_0=0 > +datum=WGS84 > +units=m +no_defs” > > In the past, when working with a few variables and this information, and > lat/lon grids in the NetCDF file, I would export the data (using Python) > to > lat,lon,data, and then r.in.ascii and project it that way. However, now > I'm > working with 14,000 high resolution daily records on the above projection. > > It seems the only way to do this efficiently is to work in the projection > of the data, and r.import the few other variables to this projection. Or > perhaps cdo or gdal or some other 3rd party tool can quickly re-project > everything for me to a known EPSG code. I've solved this before with > vectors, ogr2ogr, and arbitrary proj4 strings. > > Can anyone here suggest how to best get this data into GRASS? > > I've tried using gdalwarp: > > gdalwarp -of GTIFF -s_srs "+proj=stere +lat_0=90 +lat_ts=71 +lon_0=-39 > +k=1 > +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs" -t_srs EPSG:3413 > NETCDF:"${file}":RU out.tif > > But am concerned because when I then run gdalinfo on the out.tif file, it > reports strange bounds such as: > > Origin = (-909.965548780931840,1546.900601999830769) > Pixel Size = (7.477558376763198,-7.477558376763198) > Corner Coordinates: > Upper Left ( -909.966, 1546.901) (165d27'58.34"E, 89d59' 0.36"N) > Lower Left ( -909.966, -1264.661) ( 80d44'10.54"W, 89d59' 8.22"N) > Upper Right ( 794.918, 1546.901) (107d48' 8.56"E, 89d59' 2.20"N) > Lower Right ( 794.918, -1264.661) ( 12d50'53.17"W, 89d59'10.36"N) > Center ( -57.524, 141.120) (157d10'37.34"E, 89d59'54.94"N) > > Which seems to be very close to the pole (both lower and upper at >89 > deg). > This should be covering a much large portion of the Arctic.
e.g. you can use a proj4 string in the startup wizard to create a new location. just tested with your proj4 string, g.proj -j +proj=stere +lat_0=90 +lat_ts=71 +lon_0=-39 +k=1 +x_0=0 +y_0=0 +no_defs +a=6378137 +rf=298.257223563 +towgs84=0.000,0.000,0.000 +to_meter=1 HTH ----- best regards Helmut -- Sent from: http://osgeo-org.1560.x6.nabble.com/Grass-Users-f3884509.html _______________________________________________ grass-user mailing list grass-user@lists.osgeo.org https://lists.osgeo.org/mailman/listinfo/grass-user