[R-sig-Geo] Help

[Micha Silver]

On 09/01/2020 21:02, Bakary Faty wrote:
> Thank you for appreciated reply,
>
> I explane you exactly what I want to do with this R code attached.
> I want to adapt this code to my data to build an isohyet map.
> But i have some difficulties to adapt it to my case.
> I will be very happy when you will help my to adapt this R code (attached)
> to my case.
> You can find attached the you R code and my data file.
> Best regards
>


I made two small changes in your code, and it works fine:

  * First I used the suggestion (earlier in this thread) to create your
    grid.
  * Then, you had an error in your call to autoKrige.
  * After getting that right, I created isohyetal lines with 
rasterToContour

Here's my version


library(automap)
library(raster)
library(rgdal)

## READ INPUT FILE
rain_data <- read.csv(file="DATAFILE_FOR_ISOHYET.csv")

point_coords <- rain_data[c("Lon","Lat")]
coordinates(rain_data) <- point_coords
p4str <- CRS("+init=epsg:4326")
proj4string(rain_data) <- p4str

## CONVERTION TO UTM
p4str_UTM <- CRS("+init=epsg:32633")
rain_data_UTM <- spTransform(rain_data, p4str_UTM)


bb <- bbox(rain_data_UTM)
minx <- bb[1,1]
maxx <- bb[1,2]
miny <- bb[2,1]
maxy <- bb[2,2]

## EACH PIXEL WILL BE 1000 METERS
pixel <- 1000
grd <- expand.grid(x=seq(minx, maxx, by=pixel), y=seq(miny, maxy, 
by=pixel))
coordinates(grd) <- ~x+y
gridded(grd) <- TRUE
proj4string(grd) <- p4str_UTM


## KRIGING, USING AUTOKRIGE WHICH CREATES A BEST GUESS VARIOGRAM
# OK_rain <- autoKrige(Rainfall_value ~ 1, rain_data_UTM, grd)
# There is no variable "Rainfall_value" in your data,
# It is called RAIN_DATA
OK_rain = autoKrige(formula = RAIN_DATA ~1,
                     input_data = rain_data_UTM,
                     new_data = grd)

## TRASFORM TO RASTER
rain_rast <- raster(OK_rain$krige_output)

summary(rain_rast)


# Minimumn is 540, max is 1735
# So create isohyetal lines about every 100 mm and plot

isohyets = rasterToContour(rain_rast, nlevels = 12)
plot(rain_rast)
lines(isohyets, add = TRUE)


> Le jeu. 9 janv. 2020 à 18:41, Bakary Faty <bakaryfaty using gmail.com 
> <mailto:bakaryfaty using gmail.com>> a écrit :
>
>     Thank you for appreciated reply,
>
>     I explane you exactly what I want to do with this R code attached.
>     I want to adapt this code to my data to build an isohyet map.
>     But i have some difficulties to adapt it to my case.
>     I will be very happy when you will help my to adapt this R code
>     (attached)
>     to my case.
>     You can find attached the you R code, my data file and my sahefile
>     of watershed.
>
>     Best regards
>
>
>     Le jeu. 9 janv. 2020 à 17:47, Ben Tupper <btupper using bigelow.org
>     <mailto:btupper using bigelow.org>> a écrit :
>
>         Welcome to r-sig-geo!
>
>         I don't think that you haven't provided us enough information
>         so that we can help.  On the other hand, does the example
>         below using expand.grid help?
>
>         minx <- 20
>         maxx <- 25
>         miny <- 31
>         maxy <- 36
>         pixel <- 1
>         grd <- expand.grid(x = seq(minx, maxx, by=pixel), y =
>         seq(miny, maxy, by=pixel))
>
>         Ben
>
>         On Thu, Jan 9, 2020 at 11:41 AM Bakary Faty
>         <bakaryfaty using gmail.com <mailto:bakaryfaty using gmail.com>> wrote:
>
>
>             Dear,
>
>             I'm writing to express my wish to join R-sig-geo list users.
>             I was doing a search on the net to know how to build an
>             isohyet map and I came across this R code.
>             However, I stumbled upon a problem from the line :
>             grd <- expand.grid(x=seq(minx, maxx, by=pixel),
>             y=seq(miny, maxy, by=pixel)),
>             I get the following error message:
>             default method not implemented for type 'S4'. I want to
>             know how resolve this error.
>
>             Also, I would like to ask you only at the line level:
>             minx <- rain_data_UTM at bbox[1,1]
>             maxx <- rain_data_UTM at bbox[1,2]
>             miny <- rain_data_UTM at bbox[2,1]
>             maxy <- rain_data_UTM at bbox[2,2],
>             if I should limit myself to "rain_data_UTM" or write
>             completely:
>             rain_data_UTM at bbox[,].
>              By the way, this is the pointfile I reconstructed.
>             You can find it attached to the mail.
>
>             Thanks in advance
>
>             Best regards
>
>
>
>             Bakary
>             _______________________________________________
>             R-sig-Geo mailing list
>             R-sig-Geo using r-project.org <mailto:R-sig-Geo using r-project.org>
>             https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
>
>         -- 
>         Ben Tupper
>         Bigelow Laboratory for Ocean Science
>         West Boothbay Harbor, Maine
>         http://www.bigelow.org/
>         https://eco.bigelow.org
>
>
>
>     -- 
>
>
>
>     Bakary
>
>
>
> -- 
>
>
>
> Bakary
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo using r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

-- 
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
https://orcid.org/0000-0002-1128-1325