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

[Roger Bivand]

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.

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