just glancing through the email chain. looks like you were able to get things plotted....but you are still having an issue with some points?
here is a link for information regarding the grid staggering used in WRF. http://www2.mmm.ucar.edu/wrf/users/docs/user_guide_V3/users_guide_chap3.htm > Some remarks: > (1) I just took the first pair of coordinates as derived from gdalinfo(" > geo_em.d01.nc") > you will find 4 different coordinate pairs (i did not proof which one is > right > The data is staggered (as outlined by Dominik) So some of the corner > coordinates belongs to the staggered data and the others coordinates to > the unstaggered ones. On Tue, Feb 23, 2016 at 4:30 PM, Agus Camacho <agus.cama...@gmail.com> wrote: > Thanks Michael, that was a good job and despite the map looks nice now, > the locations cant still be plotted adequately, either as raw lon lat or as > spatial points with different projections. > > > x=read.csv("C:/Users/Agus/Dropbox/Functional Biogeography - > Copy/scripts/class 4 maxent model/clean urosaurus records.csv") > x=x[,1:3] > colnames(x) > coordinates(x)=~lon+lat > projection(x)="+proj=lcc +lat_1=48 +lat_2=33 +lon_0=-100 +ellps=WGS84" > x=spTransform(x, CRS("+proj=lcc +lat_1=25.0 +lat_2=45.0 +lat_0=38.0 > +lon_0=-100.0 +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs")) > points(x) > > 2016-02-23 16:24 GMT-07:00 Michael Sumner <mdsum...@gmail.com>: > >> >> >> On Wed, 24 Feb 2016 at 06:07 Dominik Schneider < >> dominik.schnei...@colorado.edu> wrote: >> >>> This looks like WRF <http://www.wrf-model.org/index.php> data. I just >>> dealt with this. >>> The data is on a sphere as opposed to WGS84 so you need +ellps=sphere >>> +a=6370000 +b=6370000 +units=m >>> >>> +proj=lcc which is usually what wrf is run with. >>> The tricky part is: >>> +lat_1=25.0 +lat_2=45.0 +lat_0=38.0 +lon_0=-100.0 >>> because every WRF run is different (the WRF Preprocessing System >>> optimizes the projection for the domain). >>> and then there is probably no shift so you need(?) +x_0=0 +y_0=0 >>> >>> This gives: >>> +proj=lcc +lat_1=25.0 +lat_2=45.0 +lat_0=38.0 +lon_0=-100.0 >>> +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs >>> >>> >> This looks right to me: >> >> library(raster) >> library(rgdal) >> prj <- "+proj=lcc +lat_1=25.0 +lat_2=45.0 +lat_0=38.0 +lon_0=-100.0 >> +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs" >> >> r <- raster("results_us_future_output_none_0.nc", varname = "dctmx") >> ## use print(r) to see the nc header dump >> >> ## lons/lats from NetCDF >> lon <- raster("results_us_future_output_none_0.nc", varname = "lon") >> lat <- raster("results_us_future_output_none_0.nc", varname = "lat") >> >> ## the true projected coordinates (??) >> xy <- project(cbind(values(lon), values(lat)), prj) >> >> ## looks right >> plot(xy, pch = ".") >> >> ## this fails, so it's not exactly right >> ##grd <- rasterFromXYZ(cbind(xy, 0) >> >> ## could use sp::points2grid tools to control tolerances, or just fudge it >> >> ex <- extent(xy) >> resolution <- c(round(max(diff(sort(unique(xy[,1]))))), >> round(max(diff(sort(unique(xy[,2])))))) >> ex <- ex + c(-1, 1, -1, 1) * rep(resolution / 2, 2) >> >> ## finally >> mp <- setExtent(r, ex) >> projection(mp) <- prj >> >> library(maptools) >> data(wrld_simpl) >> namer <- spTransform(subset(wrld_simpl, NAME %in% c("United States", >> "Canada", "Mexico")), prj) >> >> plot(mp) >> plot(namer, add = TRUE) >> >> >> ("Why does this crucial metadata get thrown away?", he sighs . . .) >> >> HTH >> >> >> >> >>> But, wrf users like to give out lat and long so you need to assign it: >>> +proj=longlat +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs >>> >>> and then reproject the lat/long to lcc coordinates using this string: >>> +proj=lcc +lat_1=25.0 +lat_2=45.0 +lat_0=38.0 +lon_0=-100.0 >>> +ellps=sphere +a=6370000 +b=6370000 +units=m +no_defs >>> >>> One word of caution, make sure you received the correct coordinates. >>> Some variables are run cell center while some are run at cell edge. It >>> looks like from your .nc file it was made by your collaborator so I assume >>> they are right. >>> >>> That said, another word of caution, I found that the XLAT and XLONG >>> variables from WRF output aren't very precise. There is a "geogrid" file >>> from the preprocessing system that has the domain corners, resolution, nrow >>> and ncol from which you can make a better grid using the native projection >>> system (in my case it was a 4km grid). I suggest you try to get those. >>> >>> I hope this helps... I have to run but wanted to save people too much >>> head scratching. I can get you running with more help tonight if you need. >>> Dominik >>> >>> >>> On Tue, Feb 23, 2016 at 11:27 AM, Agus Camacho <agus.cama...@gmail.com> >>> 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 >>>> >>>> >>>> >>>> 2016-02-23 2:25 GMT-07:00 Michael Sumner <mdsum...@gmail.com>: >>>> >>>> > On Tue, 23 Feb 2016 at 20:09 Roger Bivand <roger.biv...@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 >>>> > > >> fr...@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-boun...@r-project.org] On >>>> Behalf Of >>>> > > Agus Camacho >>>> > > >> Sent: 22. februar 2016 19:20 >>>> > > >> To: t...@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_...@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@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.biv...@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@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@r-project.org >>>> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo >>>> > >>>> >>>> >>>> >>>> -- >>>> Agustín Camacho Guerrero. >>>> Doutor em Zoologia. >>>> Laboratório de Herpetologia, Departamento de Zoologia, Instituto de >>>> Biociências, USP. >>>> Rua do Matão, trav. 14, nº 321, Cidade Universitária, >>>> São Paulo - SP, CEP: 05508-090, Brasil. >>>> >>>> [[alternative HTML version deleted]] >>>> >>>> _______________________________________________ >>>> R-sig-Geo mailing list >>>> R-sig-Geo@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 >> >> > > > -- > Agustín Camacho Guerrero. > Doutor em Zoologia. > Laboratório de Herpetologia, Departamento de Zoologia, Instituto de > Biociências, USP. > Rua do Matão, trav. 14, nº 321, Cidade Universitária, > São Paulo - SP, CEP: 05508-090, Brasil. > [[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo