[R-sig-Geo] "Geographical CRS given to non-conformant data" for CRS of epsg:4326

Roger Bivand Roger.Bivand at nhh.no
Sat Dec 17 16:27:57 CET 2011 

On Sat, 17 Dec 2011, Roger Bivand wrote:

> On Fri, 16 Dec 2011, Agustin Lobo wrote:
>
>> Hi!
>> 
>> I've uploaded the file here:
>> https://sites.google.com/site/filestemp2/home/world_adm0.zip
>> 
>> It was a popular file, but actually cannot find it in the diva server any 
>> more
>> According to qgis, there are no coordinates beyond +/- 90 or +/-180
>> 
>> Extents:
>> In layer spatial reference system units : xMin,yMin -180,-90 :
>> xMax,yMax 180,83.623
>> but readOGR() probably makes a more rigorous check.
>
> Right. Thanks for access to the file. If I change the name of the *.prj file, 
> I see:
>
>> wrld <- readOGR(".", "world_adm0")
>> print(bbox(wrld), digits=12)
>             min            max
> x -180.000183105 180.0000000000
> y  -90.000000000  83.6230316162
>
> so it is going sufficiently far below -180 for the validity check at 
> .Machine$double.eps ^ 0.25 (0.0001220703), but passes with a slightly slacker 
> tolerance. I suggest providing a tolerance in an option in sp to resolve 
> this. The error is primarily to stop users assigning a +proj=longlat to 
> projected coordinates (yes, it does happen). I'll look to add an options 
> mechanism to sp to handle both error/warning and tolerance for the power to 
> .Machine$double.eps.

This is now committed to the rspatial R-forge project, module sp, revision 
1179. Now it is possible either to throw a warning rather than an error, 
and/or to alter the tolerance; the error message now shows the exceeding 
value(s):

> library(rgdal)
Loading required package: sp
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.8.1, released 2011/07/09
Path to GDAL shared files: /usr/local/share/gdal
Loaded PROJ.4 runtime: Rel. 4.7.1, 23 September 2009, [PJ_VERSION: 470]
Path to PROJ.4 shared files: (autodetected)
> wrld1 <- readOGR(".", "world_adm0")
OGR data source with driver: ESRI Shapefile
Source: ".", layer: "world_adm0"
with 209 features and 3 fields
Feature type: wkbPolygon with 2 dimensions
Error in validityMethod(as(object, superClass)) :
   Geographical CRS given to non-conformant data: -180.000183105
> set_ll_warn(TRUE)
[1] TRUE
> wrld1 <- readOGR(".", "world_adm0")
OGR data source with driver: ESRI Shapefile
Source: ".", layer: "world_adm0"
with 209 features and 3 fields
Feature type: wkbPolygon with 2 dimensions
Warning message:
In validityMethod(as(object, superClass)) :
   Geographical CRS given to non-conformant data: -180.000183105
> set_ll_warn(FALSE)
[1] FALSE
> set_ll_TOL(0.2)
[1] 0.2
> wrld1 <- readOGR(".", "world_adm0")
OGR data source with driver: ESRI Shapefile
Source: ".", layer: "world_adm0"
with 209 features and 3 fields
Feature type: wkbPolygon with 2 dimensions
>

Please let me know if this helps,

Roger

>
> Had you looked at using the wrld_simpl dataset that is disributed with 
> maptools:
>
> data(wrld_simpl)
>
> and joining countries with the codes included? It is a bit out od date in not 
> having the Sudan split; it is about twice the object size of the file you are 
> using, plotting somewhat more slowly:
>
>> object.size(wrld_simpl)
> 6349008 bytes
>> object.size(wrld)
> 3250784 bytes
>> plot(wrld_simpl)
>> plot(wrld, border="red", add=TRUE)
>
> The file simplifies sea ice in Antarctica, and discards many small islands, 
> otherwise they match pretty well.
>
> Roger
>
>> 
>> Could it be possible that readOGR()  read the file with a warning
>> and print the offending coordinates?
>> 
>> readOGR() works fine with this alternative file:
>> http://www.gadm.org/data/gadm_v1_lev0_shp.zip
>> 
>> but it's huge and cannot actually work with it and have not found a
>> third alternative.
>> 
>> Thanks,
>> 
>> Agus
>> 
>> 2011/12/16 Roger Bivand <Roger.Bivand at nhh.no>:
>>> On Fri, 16 Dec 2011, Robert J. Hijmans wrote:
>>> 
>>>> Agus,
>>>> 
>>>> This is likely because there are some longitude coordinates that are  <
>>>> -180 or > 180.
>>>> OGR does not complain, but, as these are lon/lat data, sp runs
>>>> sp:::.ll_sanity  and throws an error.
>>>> 
>>>> I think it would be preferable if sp gave a warning rather than an error,
>>>> at least when longitudes are within a reasonable number (-360 -- 360 ?). 
>>>> I
>>>> think that a longitude of, e.g. 181 is perfectly interpretable, and
>>>> sometimes useful. For example to plot Russia, or Fiji in one piece in
>>>> lon/lat coordinates.
>>> 
>>> 
>>> As Agus knows, unless the file is made available, help will likely not be
>>> forthcoming. Without a traceback(), we don't know where the error is being
>>> thrown. The fact that an out of domain number is interpretable in one
>>> context doesn't necessarily make it acceptable in others, I believe. We
>>> could go to an option setting, which would slacken the error to a warning 
>>> if
>>> set. Opinions?
>>> 
>>> Roger
>>> 
>>> 
>>>> 
>>>> Robert
>>>> 
>>>> 
>>>> On Fri, Dec 16, 2011 at 4:24 AM, Agustin Lobo <alobolistas at gmail.com>
>>>> wrote:
>>>> 
>>>>> Hi!
>>>>> I' trying to read this shapefile that is just in lon,lat WGS84
>>>>>> 
>>>>>> ogrInfo(dsn=datad,layer="world_adm0v2")
>>>>> 
>>>>> Source:
>>>>> "/media/Iomega_HDD/JACOB/ExIE2011DATA/ExFireBorneoDATA/world_adm0",
>>>>> layer: "world_adm0v2"
>>>>> Driver: ESRI Shapefile number of rows 209
>>>>> Feature type: wkbPolygon with 2 dimensions
>>>>> +proj=longlat +datum=WGS84 +no_defs
>>>>> Number of fields: 3
>>>>>      name type length typeName
>>>>> 1      NAME    4     40   String
>>>>> 2 GMI_CNTRY    4      3   String
>>>>> 3    REGION    4     25   String
>>>>> 
>>>>> but I get:
>>>>>> 
>>>>>> world <-readOGR(dsn=datad,layer="world_adm0v2")
>>>>> 
>>>>> OGR data source with driver: ESRI Shapefile
>>>>> Source:
>>>>> "/media/Iomega_HDD/JACOB/ExIE2011DATA/ExFireBorneoDATA/world_adm0",
>>>>> layer: "world_adm0v2"
>>>>> with 209 features and 3 fields
>>>>> Feature type: wkbPolygon with 2 dimensions
>>>>> Error in validityMethod(as(object, superClass)) :
>>>>>  Geographical CRS given to non-conformant data
>>>>> I've tried
>>>>>> 
>>>>>> world
>>>>> 
>>>>> 
>>>>> <-readOGR(dsn=datad,layer="world_adm0v2",p4s=CRSargs(CRS("+init=epsg:4326")))
>>>>> but get the same result
>>>>> 
>>>>> can't this CRS (epsg:4326)
>>>>> +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0
>>>>> be read by readOGR() ?
>>>>> 
>>>>> Thanks
>>>>> Agus
>>>>> 
>>>>> _______________________________________________
>>>>> R-sig-Geo mailing list
>>>>> R-sig-Geo at r-project.org
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>> 
>>>> 
>>>>        [[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, NHH Norwegian School of Economics,
>>> Helleveien 30, N-5045 Bergen, Norway.
>>> voice: +47 55 95 93 55; fax +47 55 95 95 43
>>> e-mail: Roger.Bivand at nhh.no
>>> 
>>> 
>>> _______________________________________________
>>> 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, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no

More information about the R-sig-Geo mailing list