[R-sig-Geo] adapting spatial points and wrld_smpl to a reference system implicit in a .nc file
Roger Bivand
Roger.Bivand at nhh.no
Tue Feb 23 19:56:34 CET 2016
On Tue, 23 Feb 2016, Agus Camacho wrote:
> Thanks heaps to all for your effort. If I go to another GEOSTAT ill bring
> more giant crab this time.
>
> The creator of the .nc file also looked at this webpage:
> http://www.pkrc.net/wrf-lambert.html
> It seemed like the right proj4 string might be this one:
>
> +proj=lcc +lat_1=25.0 +lat_2=45.0 +lat_0=38.0
> +lon_0=-100.0 +a=6370 +b=6370 +towgs84=0,0,0 +no_defs
>
> However this projection also does not allow me to adequately plot the
> locations on the raster.
>
> Here is the .nc file. it contains several layers.
>
> https://www.dropbox.com/s/yto3linsgom3zi7/results_us_future_output_none_0.nc?dl=0
>
Thanks for the link. The proj4 string is missing +units=km (the spherical
axes are certainly not in metres), but is wrong anyway.
This says:
> gI <- GDALinfo("results_us_future_output_none_0.nc")
Error in .local(.Object, ...) :
No UNIDATA NC_GLOBAL:Conventions attribute
The values of the coordinate variables are:
library(ncdf4)
r <- nc_open("results_us_future_output_none_0.nc")
pts <- cbind(c(ncvar_get(r, "lon")), c(ncvar_get(r, "lat")))
> summary(pts)
V1 V2
Min. :-151.30 Min. :12.36
1st Qu.:-120.44 1st Qu.:25.45
Median :-100.00 Median :35.72
Mean :-100.00 Mean :35.71
3rd Qu.: -79.56 3rd Qu.:45.91
Max. : -48.70 Max. :58.33
which are obviously not a grid:
pts1 <- SpatialPoints(pts)
plot(pts1)
# gridded(pts1) <- TRUE fails miserably
This is definitely not lcc, irrespective of what the originator says.
Mike - any forensics?
Roger
>
>
> 2016-02-23 2:25 GMT-07:00 Michael Sumner <mdsumner at gmail.com>:
>
>> On Tue, 23 Feb 2016 at 20:09 Roger Bivand <Roger.Bivand at nhh.no> wrote:
>>
>>> On Tue, 23 Feb 2016, Alex Mandel wrote:
>>>
>>>> I made an attempt at it too. Investigating the original data, I'm not
>>>> sure that the projection information supplied is correct for the data
>>>> linked. When I load up the data in a unprojected space, the coordinates
>>>> don't look at all similar to any Lambert projected data I have, they
>>>> actually look like Lat/Lon in some unprojected coordinate system,
>>>> perhaps a different spheroid than expected.
>>>
>>> Does anyone have a link to the original data? Is is possible that this is
>>> the General Oblique Transformation used by modellers - that is something
>>> that feels like longlat but is recentred and oblique? Example at the very
>>> end of my GEOSTAT talk last year (slides 81-83):
>>>
>>> http://geostat-course.org/system/files/geostat_talk_150817.pdf
>>>
>>> Roger
>>>
>>>
>> For what it is worth, the General Oblique Transformation is not the only
>> example - it's very common for modellers to have a mesh that has the
>> "mostly-properties" of a projection, but is not actually describable with
>> standard transform + affine parameters. The main cases that I've seen are
>> polar stereographic, equal area or oblique Mercator. Often they really are
>> simple transforms and you can reconstruct without loss, but it's not
>> usually possible to tell without exploration. It's an interesting
>> dis-connect to see code that builds a mesh with certain properties, then
>> only stores longitudes and latitudes - when it could be done with standard
>> tools and be stored and used much more efficiently.
>>
>> (I've seen Lambert Conformal Conic and Lambert Azimuthal Equal Area
>> terminology conflated in this context too. )
>>
>> I'm also interested to explore the original data.
>>
>> Cheers, Mike.
>>
>>
>>
>>>>
>>>> -Alex
>>>>
>>>> On 02/22/2016 10:17 PM, Frede Aakmann Tøgersen wrote:
>>>>> Hi
>>>>>
>>>>> I tried to make it work but I had to give up. I wanted to reproject
>> the
>>> Lamberth conformal conic coordinates to long-lat but it didn't work.
>>>>>
>>>>> Perhaps someone can see what I did wrong. Here is what I did (data in
>> R
>>> binary format and figure in png format both attached):
>>>>>
>>>>> library(raster)
>>>>> library(maptools)
>>>>> data(wrld_simpl)
>>>>>
>>>>> r <- raster("raster.grd")
>>>>> projection(r)
>>>>> ## > NA
>>>>>
>>>>> uro <- read.table("clean urosaurus records.csv", h = TRUE, sep = ",")
>>>>> coordinates(uro) <- ~lon+lat
>>>>>
>>>>> ## Set projections for the 3 data sets
>>>>>
>>>>> ## Lamberth's confocal conic projection with given parameters
>>>>> crs(r) <- "+proj=lcc +lat_0=38.0 +lon_0=-100 +lat_1=25.0 +lat_2=45.0
>>> +ellps=WGS84"
>>>>> projection(r)
>>>>>
>>>>> ## Assume that lon, lat are geographical coordinates (degrees decimal)
>>>>> proj4string(uro) <- CRS("+proj=longlat +datum=WGS84")
>>>>>
>>>>> ## wrld_simpl is in geographical coordinates
>>>>> proj4string(wrld_simpl)
>>>>>
>>>>> ## Make figure in png format
>>>>> ## Of course plotting data with 2 different projections will give
>>>>> ## some distortions
>>>>> pdf("uro.png")
>>>>>
>>>>> plot(r)
>>>>> points(uro)
>>>>> plot(wrld_simpl, add = TRUE) # World will be clipped to extent of 'r'
>>>>>
>>>>> dev.off()
>>>>>
>>>>> extent(r)
>>>>> ## class : Extent
>>>>> ## xmin : -131.4368
>>>>> ## xmax : -68.56323
>>>>> ## ymin : 12.35567
>>>>> ## ymax : 50.26619
>>>>>
>>>>> ## Reproject the raster to long-lat
>>>>> ## This doesn't work (collapsed domain)
>>>>> rp <- projectRaster(r, crs = "+proj=longlat +datum=WGS84 +no_defs
>>> +ellps=WGS84 +towgs84=0,0,0")
>>>>>
>>>>> ## Because
>>>>>> extent(rp)
>>>>> ## class : Extent
>>>>> ## xmin : -100.0015
>>>>> ## xmax : -99.68557
>>>>> ## ymin : 37.70658
>>>>> ## ymax : 38.00046
>>>>>
>>>>> ## Save data in R binary format
>>>>> save(list = c("r", "uro", "wrld_simpl"), file = "uro.RData")
>>>>>
>>>>>
>>>>> Yours sincerely / Med venlig hilsen
>>>>>
>>>>> Frede Aakmann Tøgersen
>>>>> Specialist, M.Sc., Ph.D.
>>>>> Plant Performance & Modeling
>>>>>
>>>>> Technology & Service Solutions
>>>>> T +45 9730 5135
>>>>> M +45 2547 6050
>>>>> frtog at vestas.com
>>>>> http://www.vestas.com
>>>>>
>>>>> Company reg. name: Vestas Wind Systems A/S
>>>>> This e-mail is subject to our e-mail disclaimer statement.
>>>>> Please refer to www.vestas.com/legal/notice
>>>>> If you have received this e-mail in error please contact the sender.
>>>>>
>>>>>
>>>>>
>>>>> -----Original Message-----
>>>>> From: R-sig-Geo [mailto:r-sig-geo-bounces at r-project.org] On Behalf Of
>>> Agus Camacho
>>>>> Sent: 22. februar 2016 19:20
>>>>> To: tech at wildintellect.com
>>>>> Cc: r-sig-geo
>>>>> Subject: Re: [R-sig-Geo] adapting spatial points and wrld_smpl to a
>>> reference system implicit in a .nc file
>>>>>
>>>>> Thanks Alex, but the locations still fall in the sea when i plot them
>>> using
>>>>> your recommended Solution. I looked at the sites you proposed and they
>>> have
>>>>> other values for lat_1, lat_0, etc..
>>>>>
>>>>> 2016-02-22 11:04 GMT-07:00 Alex M <tech_dev at wildintellect.com>:
>>>>>
>>>>>> On 02/22/2016 09:50 AM, Agus Camacho wrote:
>>>>>>> Dear all,
>>>>>>>
>>>>>>> Im trying to overlap these points:
>>>>>>>
>>>>>>
>>>
>> https://www.dropbox.com/s/awdclg4cvsdngej/clean%20urosaurus%20records.csv?dl=0
>>>>>>>
>>>>>>> and a wrld_simpl object:
>>>>>>> library(maptools)
>>>>>>> data(wrld_simpl)
>>>>>>>
>>>>>>> Over this raster layer
>>>>>>>
>>>>>>
>>>
>> https://www.dropbox.com/sh/qcw174tgogpnz7s/AAByDc3TeyFe3W4nEqTFix6Oa?dl=0
>>>>>>>
>>>>>>> This rastr comes from a .nc file without a reference system. The
>>> author
>>>>>> of
>>>>>>> that .nc file gave me the following data about the .nc.
>>>>>>>
>>>>>>> The projection is *Lambert conformal conic* projection
>>>>>>> CEN_LAT = 38.0
>>>>>>> CEN_LON = -100.0
>>>>>>> TRUELAT1 = 25.
>>>>>>> TRUELAT2 = 45.
>>>>>>>
>>>>>>> However, despite i have gone through many sites in the internet, i
>>> cant
>>>>>>> figure it out:
>>>>>>>
>>>>>>> a) if that is all the data i need to set a reference system for my
>>> points
>>>>>>> and the wrld_simp object.
>>>>>>>
>>>>>>> b) how to change a typical CRS object with such data
>>>>>>>
>>>>>>> Ex.CRS ("+proj=lcc+lat_0=38.0+lon0_2=-100+ellps=WGS84")
>>>>>>>
>>>>>>> Where do i enter the TRUELAT and CENLAT values?
>>>>>>> Are there any site that explains easily what the fields in the CRS
>>> mean
>>>>>> and
>>>>>>> how to change them?
>>>>>>>
>>>>>>> Thanks in advance.
>>>>>>>
>>>>>>
>>>>>> https://github.com/OSGeo/proj.4/wiki/GenParms
>>>>>> https://trac.osgeo.org/proj/wiki/GenParms
>>>>>>
>>>>>> I believe:
>>>>>> +lat_0 = CEN_LAT Latitude of origin
>>>>>> +lat_1 = TRUELAT1 Latitude of first standard parallel
>>>>>> +lat_2 = TRUELAT2 Latitude of second standard parallel
>>>>>> +lon_0 = CEN_LON Central meridian
>>>>>>
>>>>>> proj strings are defined by the proj4 libary. It's website listed
>> above
>>>>>> and the associated mailing lists or gis stackexchange would be the
>>>>>> places to get help on it.
>>>>>> https://lists.osgeo.org/mailman/listinfo/metacrs
>>>>>>
>>>>>> It often helps to browse similar projections on
>>>>>> http://spatialreference.org/
>>>>>> http://epsg.io/
>>>>>>
>>>>>> Enjoy,
>>>>>> Alex
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>
>>>> _______________________________________________
>>>> R-sig-Geo mailing list
>>>> R-sig-Geo at r-project.org
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>> --
>>> Roger Bivand
>>> Department of Economics, Norwegian School of Economics,
>>> Helleveien 30, N-5045 Bergen, Norway.
>>> voice: +47 55 95 93 55; fax +47 55 95 91 00
>>> e-mail: Roger.Bivand at nhh.no
>>> http://orcid.org/0000-0003-2392-6140
>>> https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
>>> http://depsy.org/person/434412
>>> _______________________________________________
>>> R-sig-Geo mailing list
>>> R-sig-Geo at r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>> --
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> Kingston Tasmania 7050 Australia
>>
>> [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>
>
>
--
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 91 00
e-mail: Roger.Bivand at nhh.no
http://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en
http://depsy.org/person/434412
More information about the R-sig-Geo
mailing list