4

Since update of sp package i get a warning i did not get before:

 1: In showSRID(uprojargs, format = "PROJ", multiline = "NO") :
  Discarded ellps unknown in CRS definition: +proj=stere +lat_0=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +R=6370040 +units=km +no_defs
2: In showSRID(uprojargs, format = "PROJ", multiline = "NO") :
  Discarded datum unknown in CRS definition

This is a code i use to define CRS

 sp::CRS(
      "+proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +to_meter=1000 +no_defs"
    )

I try to create a projection and it was working before :) But since update i get this warning. What am i doing wrong? Thank you in advance,

Best regards, Andreas

Andreas
  • 397
  • 4
  • 18
  • 37

2 Answers2

2

In your definition, you have parameters a and b --- these refer to the radii of the ellipsoid representation of the earth. As they are the same, you are using a sphere. The new version of the PROJ library require that, if you use the proj4 notation, you must used the WGS84 (or closely related) datum. Otherwise, you need to provide the EPSG code, or WKT notation. So something like this:

p4 <- "+proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k=0.93301270189 +units=km +datum=WGS84"
sp::CRS(p4)
#CRS arguments:
# +proj=stere +lat_0=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +units=km  +datum=WGS84 +units=m +no_defs

Since there does not seem to be a EPSG code for this CRS, you would need to provide WKT description. I found this one

wkt <- 'PROJCS["unknown",
GEOGCS["unknown",
    DATUM["unknown",
        SPHEROID["unknown",6370040,0]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]]],
PROJECTION["Polar_Stereographic"],
PARAMETER["latitude_of_origin",90],
PARAMETER["central_meridian",10],
PARAMETER["scale_factor",0.93301270189],
PARAMETER["false_easting",0],
PARAMETER["false_northing",0],
UNIT["kilometre",1000,
    AUTHORITY["EPSG","9036"]],
AXIS["Easting",SOUTH],
AXIS["Northing",SOUTH]]'

But that still does not work because of the datum being unknown 

CRS(SRS_string=wkt)

#CRS arguments:
# +proj=stere +lat_0=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +R=6370040 +units=km +no_defs 
#Warning messages:
#1: In showSRID(SRS_string, format = "PROJ", multiline = "NO") :
#  Discarded ellps unknown in CRS definition: +proj=stere +lat_0=90 +lon_0=10 +k=0.93301270189 +x_0=0 +y_0=0 +R=6370040 +units=km +no_defs
#2: In showSRID(SRS_string, format = "PROJ", multiline = "NO") :
#  Discarded datum unknown in CRS definition
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • 1
    Hi Robert, thank you for your reply. Im reading Radolan data: https://stackoverflow.com/questions/45612285/make-a-gps-lat-long-lookup-table-from-radolan-data-in-r-convert-stereographic-p . to_meter=1000 is the size of a cell. a and b i think define overall area. If i use your suggestion i get the following waring : In rgdal::rawTransform(projto_int, projfrom, nrow(xy), xy[, 1], : 3376620 projected point(s) not finite – Andreas Aug 12 '20 at 15:52
  • perhaps because I omitted the `units=km`. I added that back in – Robert Hijmans Aug 12 '20 at 19:40
  • Yes this seem to remove the warning. Strangely i get different result regarding my plot now. Either i did it wrong before or something is wrong now. – Andreas Aug 12 '20 at 19:51
1

There were some big changes in sp and other spatial packages in the recent releases. See, for example:

https://www.r-spatial.org/r/2020/03/17/wkt.html

and

https://www.r-spatial.org/r/2020/06/17/s2.html

Relevant to your issues is there is a difference between crs and CRS. Perhaps things are getting stuck there.

Tim
  • 79
  • 8
  • Thank you for you response Tim. I also discovered that there was a major update i just can not grasp the problem that i have now. I thought maybe someone see the issue right away. – Andreas Aug 12 '20 at 13:31