Re: [R-sig-Geo] Coverting CRS in R to OSGB 1936

2022-11-29 Thread Rob Robinson
Or use
?stars::st_transform

The EPSG for OSGB36 is 27700 (https://epsg.io/27700), so something like...
newraster <- st_transform(raster, CRS=27700)
should work, but you might also want to look at st_warp() if you need a
regular grid in your new (OS) raster

Cheers
Rob

*** Learn about Britain's Birds at www.bto.org/birdfacts
  **

Dr Rob Robinson, Associate Director - Research (he/him)
Hon Reader: Univ East Anglia | Visiting Researcher: Swiss Ornithological
Institute
British Trust for Ornithology, The Nunnery, Thetford, Norfolk, IP24 2PU
Ph: +44 (0)1842 750050   T: @btorobrob
E: rob.robin...@bto.org  W: www.bto.org/about-bto/our-staff/rob-robinson

 "How can anyone be enlightened, when truth is so poorly lit"



On Tue, 29 Nov 2022 at 18:00, Alaba Boluwade 
wrote:

> Try:
>
> library(raster)
>
> new_raster<-spTransform(old_raster,CRS("+proj=longlat +ellps=airy
> +no_defs"))
>
> Good luck
> Alaba
>
> 
> From: R-sig-Geo  on behalf of Nick Wray <
> nickmw...@gmail.com>
> Sent: November 29, 2022 10:25 AM
> To: r-sig-geo@r-project.org 
> Subject: [R-sig-Geo] Coverting CRS in R to OSGB 1936
>
> Hello   I have shapefile data (of the river Tweed catchment in N England
> and Scotland) which I have been able to plot precipitation and temperature
> data onto without too many problems.  These data have the OSGB 1936/British
> National Grid  CRS projection.
>
> But I now have a raster of land cover (from the Centre for Ecology and
> Hydrology) covering the catchment area which I want to incorporate into my
> processing.  I have downloaded the raster into R and its properties are
>
> class  : RasterLayer
> band   : 1  (of  5  bands)
> dimensions : 4914, 8861, 43542954  (nrow, ncol, ncell)
> resolution : 25, 25  (x, y)
> extent : 261400, 482925, 576750, 699600  (xmin, xmax, ymin, ymax)
> crs: +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=40
> +y_0=-10 +ellps=airy +units=m +no_defs
> source : LCM.tif
> names  : LCM
> values : 0, 255  (min, max)
>
> So it has a different CRS - my question is - how (in R) can I reset the
> raster CRS to be the same as my earlier data sets ie OSGB 1936 British
> National Grid?  I haven't been able to find anything obvious on the net.
> Thanks Nick Wray
>
> [[alternative HTML version deleted]]
>
> ___
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
> [[alternative HTML version deleted]]
>
> ___
> R-sig-Geo mailing list
> R-sig-Geo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

[[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo


Re: [R-sig-Geo] Coverting CRS in R to OSGB 1936

2022-11-29 Thread Alaba Boluwade
Try:

library(raster)

new_raster<-spTransform(old_raster,CRS("+proj=longlat +ellps=airy +no_defs"))

Good luck
Alaba


From: R-sig-Geo  on behalf of Nick Wray 

Sent: November 29, 2022 10:25 AM
To: r-sig-geo@r-project.org 
Subject: [R-sig-Geo] Coverting CRS in R to OSGB 1936

Hello   I have shapefile data (of the river Tweed catchment in N England
and Scotland) which I have been able to plot precipitation and temperature
data onto without too many problems.  These data have the OSGB 1936/British
National Grid  CRS projection.

But I now have a raster of land cover (from the Centre for Ecology and
Hydrology) covering the catchment area which I want to incorporate into my
processing.  I have downloaded the raster into R and its properties are

class  : RasterLayer
band   : 1  (of  5  bands)
dimensions : 4914, 8861, 43542954  (nrow, ncol, ncell)
resolution : 25, 25  (x, y)
extent : 261400, 482925, 576750, 699600  (xmin, xmax, ymin, ymax)
crs: +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=40
+y_0=-10 +ellps=airy +units=m +no_defs
source : LCM.tif
names  : LCM
values : 0, 255  (min, max)

So it has a different CRS - my question is - how (in R) can I reset the
raster CRS to be the same as my earlier data sets ie OSGB 1936 British
National Grid?  I haven't been able to find anything obvious on the net.
Thanks Nick Wray

[[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

[[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo


[R-sig-Geo] Coverting CRS in R to OSGB 1936

2022-11-29 Thread Nick Wray
Hello   I have shapefile data (of the river Tweed catchment in N England
and Scotland) which I have been able to plot precipitation and temperature
data onto without too many problems.  These data have the OSGB 1936/British
National Grid  CRS projection.

But I now have a raster of land cover (from the Centre for Ecology and
Hydrology) covering the catchment area which I want to incorporate into my
processing.  I have downloaded the raster into R and its properties are

class  : RasterLayer
band   : 1  (of  5  bands)
dimensions : 4914, 8861, 43542954  (nrow, ncol, ncell)
resolution : 25, 25  (x, y)
extent : 261400, 482925, 576750, 699600  (xmin, xmax, ymin, ymax)
crs: +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=40
+y_0=-10 +ellps=airy +units=m +no_defs
source : LCM.tif
names  : LCM
values : 0, 255  (min, max)

So it has a different CRS - my question is - how (in R) can I reset the
raster CRS to be the same as my earlier data sets ie OSGB 1936 British
National Grid?  I haven't been able to find anything obvious on the net.
Thanks Nick Wray

[[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo


Re: [R-sig-Geo] Convert geojson file to R

2022-11-29 Thread Barry Rowlingson
Ah ha, its a file with no line endings. Count the lines, zero:

$ wc -l countrymasks.geojson
0 countrymasks.geojson

Although adding a line ending to the file still produces the error.
Hmm. Even running it through `json_pp` to get it on multiple lines
results in the error:


$ more countrymasks_pp.geojson
{
   "features" : [
  {
 "geometry" : {
"coordinates" : [
...


> cm = sf::st_read("./countrymasks_pp.geojson")
Reading layer `countrymasks_pp' from data source
  `/nobackup/rowlings/Downloads/SO/countrymasks_pp.geojson' using
driver `GeoJSON'
Error in CPL_read_ogr(dsn, layer, query, as.character(options), quiet,  :
  attempt to set index 210/210 in SET_STRING_ELT

Weirdly weird, `ogrinfo` thinks there's 210 features, but `st_read`
with my "query" fix gets 214...


$ ogrinfo -so -al countrymasks_pp.geojson
INFO: Open of `countrymasks_pp.geojson'
  using driver `GeoJSON' successful.

Layer name: countrymasks_pp
Geometry: Unknown (any)
Feature Count: 210   <--- 210 features


> cm = sf::st_read("./countrymasks_pp.geojson", query="select * from 
> countrymasks_pp where 1 = 1")
Reading query `select * from countrymasks_pp where 1 = 1'
from data source
`/nobackup/rowlings/Downloads/SO/countrymasks_pp.geojson' using driver
`GeoJSON'
Simple feature collection with 214 features and 15 fields
Geometry type: MULTIPOLYGON


I now notice some of the rows have incomplete property sets, for
example Paraguay has:

"properties" : {
"ISIPEDIA" : "PRY",
"NAME" : "Paraguay",
"an_crop" : "t",
"an_range" : "t",
"asap0_id" : 178,
"asap_cntry" : "f",
"g1_units" : 17,
"isocode" : "PY",
"km2_crop" : 80032,
"km2_rang2" : 20908,
"km2_tot" : 399367,
"name0" : "Paraguay",
"name0_shr" : "Paraguay"
 },

but Palestine has:

 "properties" : {
"ISIPEDIA" : "PSE",
"NAME" : "Palestine, State of",
"isocode" : "PS",
"km2_crop" : 84,
"km2_rang2" : 3423,
"km2_tot" : 6224,
"name0" : "Palestine, State of"
 }

Then there "Caribbean island small states" which has a vector in its properties:

   "properties" : {
"ISIPEDIA" : "CSID",
"NAME" : "Caribbean island small states",
"country_codes" : [
   "BRB",
   "DMA",
   "VIR",
   "BHS",
   "GRD",
   "ATG",
   "ANT",
   "LCA",
   "BLZ",
   "CYM",
   "VCT"
]

I wonder if these irregularities are causing odd problems with the
various gdal ways of parsing this...

Ugh. Let's all use geopackages

B

On Tue, Nov 29, 2022 at 11:34 AM Edzer Pebesma
 wrote:
>
> Interestingly, what seems to works is
>
> readLines('countrymasks.geojson') |> st_read() -> r
>
> with a warning:
>
> Warning message:
> In readLines("countrymasks.geojson") :
>incomplete final line found on 'countrymasks.geojson'
>
>
> On 29/11/2022 00:58, Miluji Sb wrote:
> > Thank you. I will report this bug (I did not have the confidence to call
> > this a bug before).
> >
> > Even using your code, I get the same output.
> >
> > structure(list(X = c(-67.3804401, -67.36091, -67.38058,
> > -67.339709998, -67.3780199, -67.32211), Y = c(-55.565569996,
> > -55.584009899, -55.600410001, -55.614969997, -55.63521,
> > -55.640089997), L1 = c(1, 1, 1, 1, 1, 1), L2 = c(1, 1, 1,
> > 1, 1, 1), L3 = c(1, 1, 1, 1, 1, 1)), row.names = c(NA, 6L), class =
> > "data.frame")
> >
> > Thank you again.
> >
> > On Mon, Nov 28, 2022 at 11:13 PM Barry Rowlingson 
> > wrote:
> >
> >> This seems to be a weird bug in `st_read`. If you read it with an SQL
> >> query that matches every row it works:
> >>
> >>> js = st_read("./countrymasks.geojson", query="select * from countrymasks
> >> where 1 = 1")
> >> Reading query `select * from countrymasks where 1 = 1' from data
> >> source `/home/rowlings/Downloads/countrymasks.geojson' using driver
> >> `GeoJSON'
> >> Simple feature collection with 214 features and 15 fields
> >> Geometry type: MULTIPOLYGON
> >> Dimension: XY
> >> Bounding box:  xmin: -180 ymin: -55.79439 xmax: 180 ymax: 83.62742
> >> Geodetic CRS:  WGS 84
> >>
> >> But leave out the query and you get that C code level error. Another
> >> equivalent query would be "select * from countrymasks" (without the
> >> "where" clause) but this
> >> triggers the error too. Very odd. Worth reporting as a bug?
> >>
> >> Barry
> >>
> >>
> >>
> >>
> >> On Mon, Nov 28, 2022 at 8:08 PM Miluji Sb  wrote:
> >>>
> >>> Thank you for reply. When I try using sf, I get the following error;
> >>>
> >>> Error in CPL_read_ogr(dsn, layer, query, as.character(options), quiet,  :
> >>>attempt to set index 210/210 in SET_STRING_ELT.
> >>>
> >>> Thanks again!
> >>>
> >>> 

Re: [R-sig-Geo] Convert geojson file to R

2022-11-29 Thread Edzer Pebesma

Interestingly, what seems to works is

readLines('countrymasks.geojson') |> st_read() -> r

with a warning:

Warning message:
In readLines("countrymasks.geojson") :
  incomplete final line found on 'countrymasks.geojson'


On 29/11/2022 00:58, Miluji Sb wrote:

Thank you. I will report this bug (I did not have the confidence to call
this a bug before).

Even using your code, I get the same output.

structure(list(X = c(-67.3804401, -67.36091, -67.38058,
-67.339709998, -67.3780199, -67.32211), Y = c(-55.565569996,
-55.584009899, -55.600410001, -55.614969997, -55.63521,
-55.640089997), L1 = c(1, 1, 1, 1, 1, 1), L2 = c(1, 1, 1,
1, 1, 1), L3 = c(1, 1, 1, 1, 1, 1)), row.names = c(NA, 6L), class =
"data.frame")

Thank you again.

On Mon, Nov 28, 2022 at 11:13 PM Barry Rowlingson 
wrote:


This seems to be a weird bug in `st_read`. If you read it with an SQL
query that matches every row it works:


js = st_read("./countrymasks.geojson", query="select * from countrymasks

where 1 = 1")
Reading query `select * from countrymasks where 1 = 1' from data
source `/home/rowlings/Downloads/countrymasks.geojson' using driver
`GeoJSON'
Simple feature collection with 214 features and 15 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box:  xmin: -180 ymin: -55.79439 xmax: 180 ymax: 83.62742
Geodetic CRS:  WGS 84

But leave out the query and you get that C code level error. Another
equivalent query would be "select * from countrymasks" (without the
"where" clause) but this
triggers the error too. Very odd. Worth reporting as a bug?

Barry




On Mon, Nov 28, 2022 at 8:08 PM Miluji Sb  wrote:


Thank you for reply. When I try using sf, I get the following error;

Error in CPL_read_ogr(dsn, layer, query, as.character(options), quiet,  :
   attempt to set index 210/210 in SET_STRING_ELT.

Thanks again!

On Mon, Nov 28, 2022 at 1:50 PM Josiah Parry 

wrote:



You're going to want to read the file with sf.

Try object <- sf::st_read("~countrymasks.geojson")

On Mon, Nov 28, 2022 at 7:09 AM Miluji Sb  wrote:


Greetings everyone,

I would like to convert the geojson file (



https://drive.google.com/file/d/18h3sOjZg5jp5euLTWRi5mC40Sja8TZDN/view?usp=sharing

)
to a dataframe - essentially obtain which has coordinates matched to a
country.

I have tried the following;

###
states <- geojsonsf::geojson_sf("~/countrymasks.geojson")
geo <- geojsonsf::sf_geojson(states)
sf <- sf::st_read(geo, quiet = T )
df <- as.data.frame(sf::st_coordinates(sf) )
##

But I get the following output, I am a bit lost. Any help will be

highly

appreciated.

Best,

   structure(list(X = c(-67.3804401, -67.36091, -67.38058,
-67.339709998, -67.3780199, -67.32211), Y =
c(-55.565569996,
-55.584009899, -55.600410001, -55.614969997, -55.63521,
-55.640089997), L1 = c(1, 1, 1, 1, 1, 1), L2 = c(1, 1, 1,
1, 1, 1), L3 = c(1, 1, 1, 1, 1, 1)), row.names = c(NA, 6L), class =
"data.frame")

 [[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo





 [[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo




[[alternative HTML version deleted]]

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo


--
Edzer Pebesma
Institute for Geoinformatics
Heisenbergstrasse 2, 48151 Muenster, Germany
Phone: +49 251 8333081

___
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo