[GRASS-user] making fcell raster maps using climate data
Pierre Roudier
pierre.roudier at gmail.com
Mon May 30 02:33:47 EDT 2011
> Pierre,
> I cannot thank you enough for all these detailed, step-by-step explanations.
No worries ;)
> I am familiar with "v.surf.rst" as an
> interpolation method, so I will use that.
Then you can probably insert a system call to GRASS instead of my
kriging routine in make_maps().
Pierre
> Thank you again,
> Bulent
>
> On Sun, May 29, 2011 at 11:14 PM, Pierre Roudier <pierre.roudier at gmail.com>
> wrote:
>>
>> Bulent,
>>
>> > Thank you for taking time and helping out. I attached a short version of
>> > data in (.xls) format. The file contains climate model results between 0
>> > and
>> > 800 BP, for every century, at a single weather station (Van-Turkey).
>>
>> Thanks. This is much clearer :)
>>
>> > When you say "reformatting" do you actually mean organizing data for
>> > importing as
>> > a vector point?
>>
>> Yes. Consider the following R session. I just exported your XLS data
>> sheet into CSV in LibreOffice.
>>
>> > df <- read.csv('Van_Precipitation.csv')
>> > df
>> UTM_E UTM_N ELEV DATE Jan Feb Mar Apr May
>> 1 502500 4398240 1661 0 35.35559 34.10211 41.08517 54.96244 53.24204
>> 2 502500 4398240 1661 -100 35.35434 34.10950 41.11511 54.98155 51.01671
>> 3 502500 4398240 1661 -200 34.95282 34.74009 43.04097 56.12360 45.02433
>> 4 502500 4398240 1661 -300 35.35213 34.12072 41.16004 55.01019 48.52734
>> 5 502500 4398240 1661 -400 35.20788 34.33214 41.78447 55.40257 18.66510
>> 6 502500 4398240 1661 -500 35.10954 34.49606 42.30229 55.70860 49.07919
>> 7 502500 4398240 1661 -600 35.06167 34.58229 42.58510 55.86880 48.44339
>> 8 502500 4398240 1661 -700 35.00383 34.69651 42.97523 56.08147 45.40806
>> 9 502500 4398240 1661 -800 35.01142 34.68574 42.94451 56.06451 45.70229
>> Jun Jul Aug Sep Oct Nov Dec
>> Annual
>> 1 17.729643 22.593137 0.7567056 11.717915 54.03994 41.26782 31.98477
>> 398.8373
>> 2 18.018884 22.797248 2.0428841 11.654079 53.91107 41.37436 31.99786
>> 398.3736
>> 3 10.662323 9.111732 4.5243246 4.254332 44.57887 48.04743 33.02368
>> 368.0845
>> 4 18.508462 23.089403 10.5560756 11.613746 53.75689 41.53623 32.01759
>> 405.2488
>> 5 21.111523 21.222403 19.7168880 10.175261 50.41907 43.68269 32.32056
>> 384.0406
>> 6 19.838364 3.158531 0.5505311 10.543923 47.83534 45.58733 32.58478
>> 376.7945
>> 7 15.610073 17.096418 8.2630408 12.148564 46.61419 46.60779 32.73860
>> 395.6199
>> 8 9.915775 8.735210 5.0867794 1.947931 45.02865 47.92435 32.96591
>> 365.7697
>> 9 10.045168 9.128862 5.3451146 2.768823 45.19434 47.83308 32.94591
>> 367.6698
>>
>> Then I'll use the reshape 2 to change the structure of df:
>>
>> > library(reshape2)
>> > res <- melt(df, id.vars=c(1:4), variable.name='month',
>> > value.name='model_output')
>> > head(res, 15)
>> UTM_E UTM_N ELEV DATE month model_output
>> 1 502500 4398240 1661 0 Jan 35.35559
>> 2 502500 4398240 1661 -100 Jan 35.35434
>> 3 502500 4398240 1661 -200 Jan 34.95282
>> 4 502500 4398240 1661 -300 Jan 35.35213
>> 5 502500 4398240 1661 -400 Jan 35.20788
>> 6 502500 4398240 1661 -500 Jan 35.10954
>> 7 502500 4398240 1661 -600 Jan 35.06167
>> 8 502500 4398240 1661 -700 Jan 35.00383
>> 9 502500 4398240 1661 -800 Jan 35.01142
>> 10 502500 4398240 1661 0 Feb 34.10211
>> 11 502500 4398240 1661 -100 Feb 34.10950
>> 12 502500 4398240 1661 -200 Feb 34.74009
>> 13 502500 4398240 1661 -300 Feb 34.12072
>> 14 502500 4398240 1661 -400 Feb 34.33214
>> 15 502500 4398240 1661 -500 Feb 34.49606
>>
>> Here, I got the value of the result model for each month in each year,
>> for each station.
>>
>> The following step is interpolation. This is not clear what you want
>> to do, but I guess what you want is to generate one raster layer per
>> month and per year.
>>
>> I propose to use the plyr library, which allows you to split your data
>> in function of one or more variables (for us, month and year), and
>> apply a function to it.
>>
>> Disclaimer: the following script is untested (that's just the general
>> idea) and should be launched from a R session WITHIN GRASS. Also, I
>> did not include the projection of the variables "grid" and "data".
>>
>> # Generate an interpolation grid (the same for each raster to be
>> generated)
>> resolution <- 5 # your resolution
>> grid <- expand.grid(x = seq(min(df$UTM_E), max(df$UTM_E), by =
>> resolution), y = seq(min(df$UTM_N), max(df$UTM_N), by = resolution))
>> grid <- points2grid(grid)
>>
>> # This is a function that generate a map from each "slice" of data
>> make_map <- function(data, grid){
>> # dependencies
>> require(sp)
>> require(rgdal)
>> require(spgrass6)
>> require(gstat) # for kriging
>> require(automap) # to automate kriging
>>
>> # converting the data into a SpatialPointsDataFrame
>> coordinates(data) <- ~UTM_E + UTM_N
>>
>> # INSERT YOUR FAVORITE INTERPOLATION METHOD HERE
>> res <- autoKrige(model_output ~ 1, data, grid)
>>
>> # write it as a raster here
>> layer_name <- paste("model_output_", unique(data$month), "-",
>> unique(data$DATE), sep ="") # a name for the current layer
>> writeRAST6(res, layer_name, zcol="var1.pred", ignore.stderr=TRUE,
>> useGDAL=TRUE)
>>
>> # return the result
>> return(res)
>> }
>>
>> # Here is the plyr call - split the data, apply the function
>> "make_map", which interpolates the data, creates a layer name and
>> writes it in GRASS. Finally, it writes back everything in "rasters",
>> which is a list.
>> library(plyr)
>> rasters <- dlply(.data = res, .variables = .(month, DATE), .fun =
>> make_map, grid = grid)
>>
>> Again, untested and probably buggy. Probably others can do better and
>> quicker with a bit of awk magic.
>>
>> Hope this helps,
>>
>> Pierre
>>
>>
>>
>>
>> > Cheers,
>> > Bulent
>> >
>> > On Sun, May 29, 2011 at 8:19 PM, Pierre Roudier
>> > <pierre.roudier at gmail.com>
>> > wrote:
>> >>
>> >> No you probably don't need to do so. What R allows you to do, is to
>> >> reformat your data frame so you do not need 400 files!
>> >>
>> >> Could you post a (short) overview of the structure of your data? From
>> >> what I understand, I'd use R and the plyr package to reformat your
>> >> data before attempting an import.
>> >>
>> >> Cheers,
>> >>
>> >> Pierre
>> >>
>> >> 2011/5/30 Bulent Arikan <bulent.arikan at gmail.com>:
>> >> > Hi Pierre,
>> >> > Thank you for answering. I use R too, so your reply was definitely
>> >> > helpful.
>> >> > I am a bit confused –and, I do not mean that I am asking for an
>> >> > answer
>> >> > from
>> >> > you when writing this– about the steps involving GRASS. So far I
>> >> > tried
>> >> > "r.in.ascii" and "r.in.bin" without success.Simply, the CSV file I
>> >> > have
>> >> > has
>> >> > 400 rows (each row represents one year) and I am not sure if I need
>> >> > to
>> >> > make
>> >> > individual CSV files for each year (means 400 files!!!) before
>> >> > importing
>> >> > into GRASS.
>> >> > Thank you for writing.
>> >> > Bulent
>> >> > On Sun, May 29, 2011 at 5:40 PM, Pierre Roudier
>> >> > <pierre.roudier at gmail.com>
>> >> > wrote:
>> >> >>
>> >> >> Hi Bulent,
>> >> >>
>> >> >> My first idea - but I may be biased as I come from the R world -
>> >> >> would
>> >> >> be to make the imports steps of your data from R. R has a XLS
>> >> >> reader(but I usually find it safer to first export your spreadsheet
>> >> >> to
>> >> >> a CSV file).
>> >> >>
>> >> >> Hope this helps,
>> >> >>
>> >> >> Pierre
>> >> >>
>> >> >> 2011/5/30 Bulent Arikan <bulent.arikan at gmail.com>:
>> >> >> > Dear List,
>> >> >> > I have an Excel file that contains precipitation data in 13
>> >> >> > columns
>> >> >> > (mean
>> >> >> > values for each month and annual average) from a region between
>> >> >> > certain
>> >> >> > years (in rows): eventually showing temporal changes in
>> >> >> > precipitation
>> >> >> > values. The data are taken at a weather station whose coordinates
>> >> >> > are
>> >> >> > known.
>> >> >> > I want to use this file in order to make an fcell raster map in
>> >> >> > GRASS.
>> >> >> > After
>> >> >> > some research I found few leads and mostly on interpolation of
>> >> >> > data.
>> >> >> > However, I am confused about the actual process turning the data
>> >> >> > in
>> >> >> > Excel
>> >> >> > into a raster map. I will greatly appreciate any suggestions as to
>> >> >> > which
>> >> >> > module to use and how to arrange the data.
>> >> >> >
>> >> >> > Thank you,
>> >> >> > --
>> >> >> > BÜLENT
>> >> >> >
>> >> >> > _______________________________________________
>> >> >> > grass-user mailing list
>> >> >> > grass-user at lists.osgeo.org
>> >> >> > http://lists.osgeo.org/mailman/listinfo/grass-user
>> >> >> >
>> >> >> >
>> >> >>
>> >> >>
>> >> >>
>> >> >> --
>> >> >> Scientist
>> >> >> Landcare Research, New Zealand
>> >> >
>> >> >
>> >> >
>> >> > --
>> >> > BÜLENT ARIKAN, PhD
>> >> > Postdoctoral Scholar
>> >> > Center for Social Dynamics and Complexity &
>> >> > School of Human Evolution and Social Change
>> >> > Arizona State University
>> >> > Tempe - AZ
>> >> > 85287-2402
>> >> >
>> >>
>> >>
>> >>
>> >> --
>> >> Scientist
>> >> Landcare Research, New Zealand
>> >
>> >
>> >
>> > --
>> > BÜLENT ARIKAN, PhD
>> > Postdoctoral Scholar
>> > Center for Social Dynamics and Complexity &
>> > School of Human Evolution and Social Change
>> > Arizona State University
>> > Tempe - AZ
>> > 85287-2402
>> >
>>
>>
>>
>> --
>> Scientist
>> Landcare Research, New Zealand
>
>
>
> --
> BÜLENT ARIKAN, PhD
> Postdoctoral Scholar
> Center for Social Dynamics and Complexity &
> School of Human Evolution and Social Change
> Arizona State University
> Tempe - AZ
> 85287-2402
>
--
Scientist
Landcare Research, New Zealand
More information about the grass-user
mailing list