[GRASSLIST:270] Fwd: Developing a forest fragmentation index

[maning sambale]

Dear list,

I hope somebody can help me on this.  I am a MS student studying
forest loss and fragmentation here in the Philippines. I am interested
in implementing Mr. Kurt Ritters forest fragmentation model
(http://www.ecologyandsociety.org/vol4/iss2/art3/) and James Hurd's
adaptation to watersehd units
(http://resac.uconn.edu/publications/tech_papers/pdf_paper/Forest_Fragmentation_ASPRS2002.pdf)
in my study area and relate the fragmentation to socio-economic
drivers.

Mr. Ritters gave me the C source code for his model, however, I am not
programmer so I don't know how to implement these.    Right now, I am
doing my image analysis (i.* modules) and classification using GRASS
and import the output to Arcview were EPA has an Attila extension
(http://www.epa.gov/esd/land-sci/attila/index.htm) for classifying
them into the forest fragmenation index.  What I want to do is create
a GRASS script for the forest fragmentation index.

I have been using GRASS for just a couple of months and I am not very
familiar with most of the commands.  I am willing to collaborate with
anybody from this list doing related research on forest fragmentation.
 As I have mentioned above, I am not a programmer thus my limitation
in understanding C code.  However, I am willing to do GRASS related
contributions related to my study.

The method is described below:

The basis for the forest fragmentation index was the fragmentation
model developed by Ritters et al. (2000) which was used and expanded
by Hurd et. al (2002) on Landsat TM data. The objective of the model
is to provide additional information beyond the commonly used measures
of forest proportions (i. e. percent of forest cover) by providing
information on the spatial configuration of fragmentation and
connectivity of forest within the study
region.

The amount of forest and its occurrence is measured as adjacent forest
pixels within a 5 x 5 windows surrounding each forest pixel. This
information will be used to classify the window by the type of
fragmentation. The result was stored at the location of the center
pixel. Thus, a pixel value in the derived map refers to
"between-pixel" fragmentation around the corresponding forest location
(Riitters et al. 2000). The model generates two values that
characterize a forest pixel located at the center of a sliding window
of fixed size, Pf and Pff:
Pf is the proportion of pixels in the window that are forested. Pff
(strictly) is the proportion of all adjacent (cardinal directions
only) pixel pairs that include at least one forest pixel, for which
both pixels are forested. Pff (roughly) estimates the conditional
probability that, given a pixel of forest, its neighbor is also
forest.

Pf = ( # of forest pixels) / (# of all non water pixels)
Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
with at least one pixel forest)

The classification model identifies six fragmentation categories: (1)
interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3) transitional,
0.4 < Pf < 0.6; (4) edge, Pf > 0.6 and Pf - Pff > 0; (5) perforated,
Pf > 0.6 and Pf – Pff < 0, and (6) undetermined, Pf > 0.6 and Pf = Pff

Any tips in implementing the generation of forest fragmentation index
either as a GRASS script or a step-by-step GRASS commands would be
very helpful (r.mapcalc syntax?).

I have also attached Ritters's C code as well as my previous
correspondence to him.  Also attached are images from Hurd's study.

Thank you very and looking forward for any response.

Cheers,

Maning Sambale


---------- Forwarded message ----------
From: Kurt H Riitters <[hidden email]>
Date: Oct 23, 2004 2:46 AM
Subject: Re: Developing a forest fragmentation index
To: emmanuel sambale <[hidden email]>
Cc: Kurt H Riitters <[hidden email]>


Emanuel,
  Attached is the main C program.  As I must say to everyone, I am not able
to provide technical support for using the program, and I do not guarantee
that it is bug-free.  You can share this with anyone you want to, but if
you share it, then you must answer any questions about it.

  Probably the interesting part will be the algorithm for the moving
window, and the definition of the functions (like Pf and Pff).  If it were
me, I'd think about adapting those parts to new programs of your design.

   The input / output format is assumed to be an old format that I adapted
from somewhere else...and it is not standard and is not used by anyone
else... so you will want to re-write the input / output routines to read
files in other formats.

   I do not have any written documentation, but there are a lot of comments
in the code that will help understand the program flow.

  Lastly, since this code is a research tool, there are many places
containing references to code or subroutines that are no longer included.

    Oh, I should say... this program uses a moving window to calculate the
X and Y values that are in the model that appears in that 2000 paper.  One
run of the program makes a map of Pf; a second run makes a map of Pff.
Once you have the output maps of the Pf and Pff values, it is necessary to
put them together to perform the classification into categories like
'perforated' 'edge' etc.  Arc or envi would be convenient for that.

  Best regards,

Kurt Riitters
RWU-4803
Forest Health Monitoring
US Forest Service
3041 Cornwallis Road
Research Triangle Park NC 27709
919-549-4015
[hidden email]

> > Emmanuel,
> >  Thank you for the message.
> >   Of course you may use the methods; they are public knowledge.
> >
> >   I wrote C programs to process the landcover maps.  I am happy to send
> > them to you, however, if you are not a programmer then it will be very
> > difficult to adapt the programs to your problems.  Even another
programmer
> > may find it difficult since I am not a very good programmer either.
> >
> >  I do not have any arc or envi scripts since I use the C programs...
But
> > it is a very simple moving window (arc: focal function) approach that
> > should be easy to implement in another language...  For example, arc
should
> > be able to do the "Pf" part of the model with the Grid: focalsum
function.

[Kari Salovaara]

maning sambale wrote:
> Dear list,
>
> I hope somebody can help me on this.  I am a MS student studying
> forest loss and fragmentation here in the Philippines.

Hi,

have You looked r.le programs. You can find more info from
- GRASS GIS 6.0.2 Reference Manual raster commands
- GRASS: Documentation - Special topics -  Landscape structure analysis
at http://grass.itc.it/

An example You can find from Open Source Free Software GIS - GRASS users
conference 2002
at
http://www.ing.unitn.it/~grass/conferences/GRASS2002/proceedings/proceedings/papers_by_title.html
- In Search Of Habitats ~ Testing A Gnu Approach (Nieminen Juhana): Full
Paper, Abstract

Regards,
Kari


--
Kari Salovaara
Hanko, Finland

[Martin Wegmann]

On Tuesday 21 March 2006 10:02, Kari Salovaara wrote:
[...]
> have You looked r.le programs. You can find more info from
> - GRASS GIS 6.0.2 Reference Manual raster commands
> - GRASS: Documentation - Special topics -  Landscape structure analysis
> at http://grass.itc.it/

currently r.li, a full rewrite of r.le, is under development. Indices are
currently implemented and hopefully soon r.li will enter cvs.
Please see previous announcement
http://grass.itc.it/pipermail/grass5/2006-January/020875.html.

If you are interested to contribute by providing code/algorithms or do some
testing etc. you are welcome.

regards, Martin

>
> An example You can find from Open Source Free Software GIS - GRASS users
> conference 2002
> at
> http://www.ing.unitn.it/~grass/conferences/GRASS2002/proceedings/proceeding
>s/papers_by_title.html - In Search Of Habitats ~ Testing A Gnu Approach
> (Nieminen Juhana): Full Paper, Abstract
>
> Regards,
> Kari

--
Martin Wegmann

DLR - German Aerospace Center
German Remote Sensing Data Center
@
Dept.of Geography
Remote Sensing and Biodiversity Unit
&&
Dept. of Animal Ecology and Tropical Biology
University of Wuerzburg
Am Hubland
97074 Würzburg

phone: +49-(0)931 - 888 4797
mobile: [phone]
fax:   +49-(0)931 - 888 4961
http://www.biota-africa.org
http://www.biogis.de

[Jáchym Čepický-2]

Hallo,
I have some experience with neighborhood operations on rastermaps, so I
could help you to make the module together.

Perhaps, you could be able make some prototype with r.mapcalc module?
(see http://grass.itc.it/gdp/raster/mapcalc-algebra.pdf and
http://grass.itc.it/gdp/raster/mapcalc.pdf)

Jachym

maning sambale wrote:

>Dear list,
>
>I hope somebody can help me on this.  I am a MS student studying
>forest loss and fragmentation here in the Philippines. I am interested
>in implementing Mr. Kurt Ritters forest fragmentation model
>(http://www.ecologyandsociety.org/vol4/iss2/art3/) and James Hurd's
>adaptation to watersehd units
>(http://resac.uconn.edu/publications/tech_papers/pdf_paper/Forest_Fragmentation_ASPRS2002.pdf)
>in my study area and relate the fragmentation to socio-economic
>drivers.
>
>Mr. Ritters gave me the C source code for his model, however, I am not
>programmer so I don't know how to implement these.    Right now, I am
>doing my image analysis (i.* modules) and classification using GRASS
>and import the output to Arcview were EPA has an Attila extension
>(http://www.epa.gov/esd/land-sci/attila/index.htm) for classifying
>them into the forest fragmenation index.  What I want to do is create
>a GRASS script for the forest fragmentation index.
>
>I have been using GRASS for just a couple of months and I am not very
>familiar with most of the commands.  I am willing to collaborate with
>anybody from this list doing related research on forest fragmentation.
> As I have mentioned above, I am not a programmer thus my limitation
>in understanding C code.  However, I am willing to do GRASS related
>contributions related to my study.
>
>The method is described below:
>
>The basis for the forest fragmentation index was the fragmentation
>model developed by Ritters et al. (2000) which was used and expanded
>by Hurd et. al (2002) on Landsat TM data. The objective of the model
>is to provide additional information beyond the commonly used measures
>of forest proportions (i. e. percent of forest cover) by providing
>information on the spatial configuration of fragmentation and
>connectivity of forest within the study
>region.
>
>The amount of forest and its occurrence is measured as adjacent forest
>pixels within a 5 x 5 windows surrounding each forest pixel. This
>information will be used to classify the window by the type of
>fragmentation. The result was stored at the location of the center
>pixel. Thus, a pixel value in the derived map refers to
>"between-pixel" fragmentation around the corresponding forest location
>(Riitters et al. 2000). The model generates two values that
>characterize a forest pixel located at the center of a sliding window
>of fixed size, Pf and Pff:
>Pf is the proportion of pixels in the window that are forested. Pff
>(strictly) is the proportion of all adjacent (cardinal directions
>only) pixel pairs that include at least one forest pixel, for which
>both pixels are forested. Pff (roughly) estimates the conditional
>probability that, given a pixel of forest, its neighbor is also
>forest.
>
>Pf = ( # of forest pixels) / (# of all non water pixels)
>Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
>with at least one pixel forest)
>
>The classification model identifies six fragmentation categories: (1)
>interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3) transitional,
>0.4 < Pf < 0.6; (4) edge, Pf > 0.6 and Pf - Pff > 0; (5) perforated,
>Pf > 0.6 and Pf � Pff < 0, and (6) undetermined, Pf > 0.6 and Pf = Pff
>
>Any tips in implementing the generation of forest fragmentation index
>either as a GRASS script or a step-by-step GRASS commands would be
>very helpful (r.mapcalc syntax?).
>
>I have also attached Ritters's C code as well as my previous
>correspondence to him.  Also attached are images from Hurd's study.
>
>Thank you very and looking forward for any response.
>
>Cheers,
>
>Maning Sambale
>
>
>---------- Forwarded message ----------
>From: Kurt H Riitters <[hidden email]>
>Date: Oct 23, 2004 2:46 AM
>Subject: Re: Developing a forest fragmentation index
>To: emmanuel sambale <[hidden email]>
>Cc: Kurt H Riitters <[hidden email]>
>
>
>Emanuel,
>  Attached is the main C program.  As I must say to everyone, I am not able
>to provide technical support for using the program, and I do not guarantee
>that it is bug-free.  You can share this with anyone you want to, but if
>you share it, then you must answer any questions about it.
>
>  Probably the interesting part will be the algorithm for the moving
>window, and the definition of the functions (like Pf and Pff).  If it were
>me, I'd think about adapting those parts to new programs of your design.
>
>   The input / output format is assumed to be an old format that I adapted
>from somewhere else...and it is not standard and is not used by anyone
>else... so you will want to re-write the input / output routines to read
>files in other formats.
>
>   I do not have any written documentation, but there are a lot of comments
>in the code that will help understand the program flow.
>
>  Lastly, since this code is a research tool, there are many places
>containing references to code or subroutines that are no longer included.
>
>    Oh, I should say... this program uses a moving window to calculate the
>X and Y values that are in the model that appears in that 2000 paper.  One
>run of the program makes a map of Pf; a second run makes a map of Pff.
>Once you have the output maps of the Pf and Pff values, it is necessary to
>put them together to perform the classification into categories like
>'perforated' 'edge' etc.  Arc or envi would be convenient for that.
>
>  Best regards,
>
>Kurt Riitters
>RWU-4803
>Forest Health Monitoring
>US Forest Service
>3041 Cornwallis Road
>Research Triangle Park NC 27709
>919-549-4015
>[hidden email]
>
>  
>
>>>Emmanuel,
>>> Thank you for the message.
>>>  Of course you may use the methods; they are public knowledge.
>>>
>>>  I wrote C programs to process the landcover maps.  I am happy to send
>>>them to you, however, if you are not a programmer then it will be very
>>>difficult to adapt the programs to your problems.  Even another
>>>      
>>>
>programmer
>  
>
>>>may find it difficult since I am not a very good programmer either.
>>>
>>> I do not have any arc or envi scripts since I use the C programs...
>>>      
>>>
>But
>  
>
>>>it is a very simple moving window (arc: focal function) approach that
>>>should be easy to implement in another language...  For example, arc
>>>      
>>>
>should
>  
>
>>>be able to do the "Pf" part of the model with the Grid: focalsum
>>>

[Kari Salovaara]

Martin Wegmann wrote:

>
> currently r.li, a full rewrite of r.le, is under development. Indices are
> currently implemented and hopefully soon r.li will enter cvs.
> Please see previous announcement
> http://grass.itc.it/pipermail/grass5/2006-January/020875.html.
>
> If you are interested to contribute by providing code/algorithms or do some
> testing etc. you are welcome.
>
> regards, Martin
>

Thanks. My mistake that I didn't know about that, reason might be wrong
name of the list GRASS5.
I'm using only versions 6.x at the moment and was thinking it relates
only to versions 5.0. I used whole
night reading archives. Thanks again.
I shall check or start testing in few days as the idea and
implementation looks promising from the first sight.

See also a later message from Paolo Cavallini, [GRASSLIST:280] r.li
(landscape indices)
They have also some interesting articles under www.faunalia.com (few
even in english ;-) )

Best regards,
Kari

--
Kari Salovaara
Hanko, Finland

[maning sambale]

Thank you for the responses.  Will dig into to the manuals you recommended.
Will post later when I made some progress.  Right now i'm into the
generating the Pf which is simply r.neighbors and r.mapcalc commands .

What is bugging me is generating Pff, which counts from the moving
window the forest pixel pairs in 4 cardinal directions.  (see
attached)

Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
with at least one pixel forest)

Cheers,
Maning

PS.  how do I install r.li tarball in grass in cygwin?


On 3/22/06, Kari Salovaara <[hidden email]> wrote:

> Martin Wegmann wrote:
> >
> > currently r.li, a full rewrite of r.le, is under development. Indices are
> > currently implemented and hopefully soon r.li will enter cvs.
> > Please see previous announcement
> > http://grass.itc.it/pipermail/grass5/2006-January/020875.html.
> >
> > If you are interested to contribute by providing code/algorithms or do some
> > testing etc. you are welcome.
> >
> > regards, Martin
> >
>
> Thanks. My mistake that I didn't know about that, reason might be wrong
> name of the list GRASS5.
> I'm using only versions 6.x at the moment and was thinking it relates
> only to versions 5.0. I used whole
> night reading archives. Thanks again.
> I shall check or start testing in few days as the idea and
> implementation looks promising from the first sight.
>
> See also a later message from Paolo Cavallini, [GRASSLIST:280] r.li
> (landscape indices)
> They have also some interesting articles under www.faunalia.com (few
> even in english ;-) )
>
> Best regards,
> Kari
>
> --
> Kari Salovaara
> Hanko, Finland
>
>

--
|---------|----------------------------------------------------------|
| __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
|  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
|  |  /T  |Philippine Biodiversity Web Map                           |
| _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
|---------|----------------------------------------------------------|

[maning sambale]

---------- Forwarded message ----------
From: maning sambale <[hidden email]>
Date: Mar 23, 2006 11:48 AM
Subject: Re: [GRASSLIST:303] Re: Fwd: Developing a forest fragmentation index
To: GRASS user list <[hidden email]>


Thank you for the responses.  Will dig into to the manuals you recommended.
Will post later when I made some progress.  Right now i'm into the
generating the Pf which is simply r.neighbors and r.mapcalc commands .

What is bugging me is generating Pff, which counts from the moving
window the forest pixel pairs in 4 cardinal directions.  (see
attached)

Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
with at least one pixel forest)

Cheers,
Maning

PS.  how do I install r.li tarball in grass in cygwin?


On 3/22/06, Kari Salovaara <[hidden email]> wrote:

> Martin Wegmann wrote:
> >
> > currently r.li, a full rewrite of r.le, is under development. Indices are
> > currently implemented and hopefully soon r.li will enter cvs.
> > Please see previous announcement
> > http://grass.itc.it/pipermail/grass5/2006-January/020875.html.
> >
> > If you are interested to contribute by providing code/algorithms or do some
> > testing etc. you are welcome.
> >
> > regards, Martin
> >
>
> Thanks. My mistake that I didn't know about that, reason might be wrong
> name of the list GRASS5.
> I'm using only versions 6.x at the moment and was thinking it relates
> only to versions 5.0. I used whole
> night reading archives. Thanks again.
> I shall check or start testing in few days as the idea and
> implementation looks promising from the first sight.
>
> See also a later message from Paolo Cavallini, [GRASSLIST:280] r.li
> (landscape indices)
> They have also some interesting articles under www.faunalia.com (few
> even in english ;-) )
>
> Best regards,
> Kari
>
> --
> Kari Salovaara
> Hanko, Finland
>
>

--
|---------|----------------------------------------------------------|
| __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
|  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
|  |  /T  |Philippine Biodiversity Web Map                           |
| _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
|---------|----------------------------------------------------------|


--
|---------|----------------------------------------------------------|
| __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
|  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
|  |  /T  |Philippine Biodiversity Web Map                           |
| _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
|---------|----------------------------------------------------------|

[HamishB]

> What is bugging me is generating Pff, which counts from the moving
> window the forest pixel pairs in 4 cardinal directions.  (see
> attached)

Read about the NEIGHBORHOOD MODIFIER in the r.mapcalc help page.
That gives you access to the cell above, below, to the left, and to the
right of the cell being processed.

"..., map[0,1] refers to the cell one column to the right of the current
cell. This syntax permits the development of neighborhood-type filters
within a single map or across multiple maps."


Or if you feel really keen you can create a new filter for r.neighbors
by modifying one of the old ones?


r.mapcalc is very powerful though and should let you win.



Hamish

[maning sambale]

Dear Grass users,

Hi!

I've been trying to develop the fragmentation index using GRASS
following the suggestions from the GRASS list.  The script is a
step-by-step process how I came up with the finalindex map.  Probably
not the best way, since this is my first effort in GRASS scripting and
r.mapcalc.

you can find the script, filters and sample images here:
http://esambale.wikispaces.com/fragindex?f=print

Some concerns are:
does r.mfilter understand mask and null values?
what about the edges of the images?

#!/bin/sh

# Create pf map; Pf = (# forest pixels) / (# non-forest pixels)

# A=forest pixels
# B=non-forest pixels
# category vegation code
#1    Agriculture
#2    Brushland
#3    Clouds
#4    Forest
#5    Grassland
#6    River

# get map
r.mapcalc 1990veg=1990fragtest@PERMANENT

# set null values
r.null map=1990veg setnull=3,6

#recode data
# count number of forest pixels value=1
r.mapcalc "A=(if(1990veg == 4,1))+(if(1990veg == 1||1990veg ==
2||1990veg == 5,0))"

# count number of non-forest pixels value=1
r.mapcalc "B=if(1990veg == 1||1990veg == 2||1990veg == 4||1990veg == 5,1)"

# C number of forest pixels in a 5x5 window
# D number of non-forest pixels in a 5x5 window

r.neighbors input=A output=C method=sum size=5
r.neighbors input=B output=D method=sum size=5

# count number of forest pixels value=1

r.mapcalc "E = 1.0 * C"
r.mapcalc "F = 1.0 * D"
r.mapcalc "pf = (E/F)"

## pixels with both forest from the cardinal directions
#create filters
r.mfilter input=A output=Nf filter=filtrN.txt
r.mfilter input=A output=Sf filter=filtrS.txt
r.mfilter input=A output=Wf filter=filtrW.txt
r.mfilter input=A output=Ef filter=filtrE.txt

# recode filtrd to 2=1 1=0
r.mapcalc "Nf2=(if(Nf == 1||Nf == 0,0))+(if(Nf == 2,1))"
r.mapcalc "Sf2=(if(Sf == 1||Sf == 0,0))+(if(Sf == 2,1))"
r.mapcalc "Ef2=(if(Ef == 1||Ef == 0,0))+(if(Ef == 2,1))"
r.mapcalc "Wf2=(if(Wf == 1||Wf == 0,0))+(if(Wf == 2,1))"


#add maps
r.mapcalc "F3 = (Nf2 + Sf2 + Ef2 + Wf2)"

# recode
r.mapcalc "F4=(if(F3 == 4||Wf == 3||Wf == 2,1))+(if(Wf == 1,0))"

#count both forest pixels
r.neighbors input=F4 output=F5 method=sum size=5

r.mapcalc "F6 = 1.0 * F5"
r.mapcalc "pff = (F6/E)"

#frag model
#code    Category        Pf                    Pff
# 1     Patch            < 40%
# 2     Transitional    >= 40% but < 60%
# 3     Edge            > 60%                Pf - Pff > 0    1
# 4     Perforated        > 60%                Pf – Pff < 0    2
# 5     Interior        1
# 6        Undetermined    > 0.6                Pf = Pff        3

r.mapcalc "Pff2 = pf - pff"

r.mapcalc "Pff3 = (if (pff2 > 0,1)) + (if (pff2 < 0,2)) + (if (pff2 == 0,3))"

#0-.4        1
#.4 - .6    2
#.6 - 0.99    3
#1            4

r.mapcalc "Pff4 = (if (pff < 0.4,1)) + (if (pff >= 0.4 && pff < .6,2))
+ (if (pff >= 0.6 && pff < .99,3)) /
+ (if (pff == 1,4))"

r.mapcalc "F11 = pff4 == 1"
r.mapcalc "F21 = pff4 == 2"
r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"
r.mapcalc "F51 = pf == 1"
r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

r.mapcalc "index = if (if (F11==1,1)) + (if (F21==1,2)) + (if
(F31==1,3)) + (if (F41==1,4)) + (if (F51==1,5)) + /
(if (F61==1,6))"

r.mapcalc "indexfin = if (index == 0 ||index == 7 ||index == 8,0) + if
(index == 1,1) + if /
(index == 2,2) + if (index == 3,3) + if (index == 4,4) + if (index ==
5,5) + if (index == 6,6)"


g.remove rast=A,B,C,D,E,Ef,Ef2,F,F3,F4,F5,F6,Nf,Nf2,Sf,Wf,Wf2,F11,F21,F31,F41,F51,F61,Pff2,Pff3,Pff4,Sf2,index

Cheers,
Maning

On 3/22/06, maning sambale <[hidden email]> wrote:

> Thank you for the responses.  Will dig into to the manuals you recommended.
> Will post later when I made some progress.  Right now i'm into the
> generating the Pf which is simply r.neighbors and r.mapcalc commands .
>
> What is bugging me is generating Pff, which counts from the moving
> window the forest pixel pairs in 4 cardinal directions.  (see
> attached)
>
> Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
> with at least one pixel forest)
>
> Cheers,
> Maning
>
> PS.  how do I install r.li tarball in grass in cygwin?
>
> On 3/21/06, Jachym Cepicky <[hidden email]> wrote:
> > Hallo,
> > I have some experience with neighborhood operations on rastermaps, so I
> > could help you to make the module together.
> >
> > Perhaps, you could be able make some prototype with r.mapcalc module?
> > (see http://grass.itc.it/gdp/raster/mapcalc-algebra.pdf and
> > http://grass.itc.it/gdp/raster/mapcalc.pdf)
> >
> > Jachym
> >
> > maning sambale wrote:
> >
> > >Dear list,
> > >
> > >I hope somebody can help me on this.  I am a MS student studying
> > >forest loss and fragmentation here in the Philippines. I am interested
> > >in implementing Mr. Kurt Ritters forest fragmentation model
> > >(http://www.ecologyandsociety.org/vol4/iss2/art3/) and James Hurd's
> > >adaptation to watersehd units
> > >(http://resac.uconn.edu/publications/tech_papers/pdf_paper/Forest_Fragmentation_ASPRS2002.pdf)
> > >in my study area and relate the fragmentation to socio-economic
> > >drivers.
> > >
> > >Mr. Ritters gave me the C source code for his model, however, I am not
> > >programmer so I don't know how to implement these.    Right now, I am
> > >doing my image analysis (i.* modules) and classification using GRASS
> > >and import the output to Arcview were EPA has an Attila extension
> > >(http://www.epa.gov/esd/land-sci/attila/index.htm) for classifying
> > >them into the forest fragmenation index.  What I want to do is create
> > >a GRASS script for the forest fragmentation index.
> > >
> > >I have been using GRASS for just a couple of months and I am not very
> > >familiar with most of the commands.  I am willing to collaborate with
> > >anybody from this list doing related research on forest fragmentation.
> > > As I have mentioned above, I am not a programmer thus my limitation
> > >in understanding C code.  However, I am willing to do GRASS related
> > >contributions related to my study.
> > >
> > >The method is described below:
> > >
> > >The basis for the forest fragmentation index was the fragmentation
> > >model developed by Ritters et al. (2000) which was used and expanded
> > >by Hurd et. al (2002) on Landsat TM data. The objective of the model
> > >is to provide additional information beyond the commonly used measures
> > >of forest proportions (i. e. percent of forest cover) by providing
> > >information on the spatial configuration of fragmentation and
> > >connectivity of forest within the study
> > >region.
> > >
> > >The amount of forest and its occurrence is measured as adjacent forest
> > >pixels within a 5 x 5 windows surrounding each forest pixel. This
> > >information will be used to classify the window by the type of
> > >fragmentation. The result was stored at the location of the center
> > >pixel. Thus, a pixel value in the derived map refers to
> > >"between-pixel" fragmentation around the corresponding forest location
> > >(Riitters et al. 2000). The model generates two values that
> > >characterize a forest pixel located at the center of a sliding window
> > >of fixed size, Pf and Pff:
> > >Pf is the proportion of pixels in the window that are forested. Pff
> > >(strictly) is the proportion of all adjacent (cardinal directions
> > >only) pixel pairs that include at least one forest pixel, for which
> > >both pixels are forested. Pff (roughly) estimates the conditional
> > >probability that, given a pixel of forest, its neighbor is also
> > >forest.
> > >
> > >Pf = ( # of forest pixels) / (# of all non water pixels)
> > >Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
> > >with at least one pixel forest)
> > >
> > >The classification model identifies six fragmentation categories: (1)
> > >interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3) transitional,
> > >0.4 < Pf < 0.6; (4) edge, Pf > 0.6 and Pf - Pff > 0; (5) perforated,
> > >Pf > 0.6 and Pf � Pff < 0, and (6) undetermined, Pf > 0.6 and Pf = Pff
> > >
> > >Any tips in implementing the generation of forest fragmentation index
> > >either as a GRASS script or a step-by-step GRASS commands would be
> > >very helpful (r.mapcalc syntax?).
> > >
> > >I have also attached Ritters's C code as well as my previous
> > >correspondence to him.  Also attached are images from Hurd's study.
> > >
> > >Thank you very and looking forward for any response.
> > >
> > >Cheers,
> > >
> > >Maning Sambale
> > >
> > >
> > >---------- Forwarded message ----------
> > >From: Kurt H Riitters <[hidden email]>
> > >Date: Oct 23, 2004 2:46 AM
> > >Subject: Re: Developing a forest fragmentation index
> > >To: emmanuel sambale <[hidden email]>
> > >Cc: Kurt H Riitters <[hidden email]>
> > >
> > >
> > >Emanuel,
> > >  Attached is the main C program.  As I must say to everyone, I am not able
> > >to provide technical support for using the program, and I do not guarantee
> > >that it is bug-free.  You can share this with anyone you want to, but if
> > >you share it, then you must answer any questions about it.
> > >
> > >  Probably the interesting part will be the algorithm for the moving
> > >window, and the definition of the functions (like Pf and Pff).  If it were
> > >me, I'd think about adapting those parts to new programs of your design.
> > >
> > >   The input / output format is assumed to be an old format that I adapted
> > >from somewhere else...and it is not standard and is not used by anyone
> > >else... so you will want to re-write the input / output routines to read
> > >files in other formats.
> > >
> > >   I do not have any written documentation, but there are a lot of comments
> > >in the code that will help understand the program flow.
> > >
> > >  Lastly, since this code is a research tool, there are many places
> > >containing references to code or subroutines that are no longer included.
> > >
> > >    Oh, I should say... this program uses a moving window to calculate the
> > >X and Y values that are in the model that appears in that 2000 paper.  One
> > >run of the program makes a map of Pf; a second run makes a map of Pff.
> > >Once you have the output maps of the Pf and Pff values, it is necessary to
> > >put them together to perform the classification into categories like
> > >'perforated' 'edge' etc.  Arc or envi would be convenient for that.
> > >
> > >  Best regards,
> > >
> > >Kurt Riitters
> > >RWU-4803
> > >Forest Health Monitoring
> > >US Forest Service
> > >3041 Cornwallis Road
> > >Research Triangle Park NC 27709
> > >919-549-4015
> > >[hidden email]
> > >
> > >
> > >
> > >>>Emmanuel,
> > >>> Thank you for the message.
> > >>>  Of course you may use the methods; they are public knowledge.
> > >>>
> > >>>  I wrote C programs to process the landcover maps.  I am happy to send
> > >>>them to you, however, if you are not a programmer then it will be very
> > >>>difficult to adapt the programs to your problems.  Even another
> > >>>
> > >>>
> > >programmer
> > >
> > >
> > >>>may find it difficult since I am not a very good programmer either.
> > >>>
> > >>> I do not have any arc or envi scripts since I use the C programs...
> > >>>
> > >>>
> > >But
> > >
> > >
> > >>>it is a very simple moving window (arc: focal function) approach that
> > >>>should be easy to implement in another language...  For example, arc
> > >>>
> > >>>
> > >should
> > >
> > >
> > >>>be able to do the "Pf" part of the model with the Grid: focalsum
> > >>>
> >
> >
>
>
> --
> |---------|----------------------------------------------------------|
> | __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
> | '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
> |  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
> |  |  /T  |Philippine Biodiversity Web Map                           |
> | _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
> |---------|----------------------------------------------------------|
>
>

--
|---------|----------------------------------------------------------|
| __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
|  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
|  |  /T  |Philippine Biodiversity Web Map                           |
| _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
|---------|----------------------------------------------------------|


--
|---------|----------------------------------------------------------|
| __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
| '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
|  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
|  |  /T  |Philippine Biodiversity Web Map                           |
| _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
|---------|----------------------------------------------------------|

[HamishB]

> I've been trying to develop the fragmentation index using GRASS
> following the suggestions from the GRASS list.  The script is a
> step-by-step process how I came up with the finalindex map.

Great to see such quick progress. This flavour of raster processing has
historically been GRASS's strongest (or at least oldest) feature.

> Probably not the best way, since this is my first effort in GRASS
> scripting and r.mapcalc.

hey, if you get to the correct answer that's what's important. Clarity
is second. Efficiency is a distant third (unless you need real-time or
it requiers weeks of time to run).

> Some concerns are:
> does r.mfilter understand mask and null values?

Don't know. NULL and MASKed cells look the same to a module AFAIK.

> what about the edges of the images?

Cells beyond the edge of the should be treated as NULL. But, again don't
know exactly what happens without testing.



hints:

1) Use r.colors to set a color for each cat number.
e.g.
r.colors 1990veg color=rules << EOF
1 yellow
2 brown
3 grey
4 0:127:0
5 green
6 blue
EOF


2) Use r.support (or cat << EOF > $MAPSET/cats/$MAP ; text ; EOF) to
set cat names for each feature in a raster map. Then they show up in the
legend, queries, etc. Try 'r.param.scale param=feature' or r.support to
see what a cats file should look like.


3) r.mapcalc efficiency

You can combine some of your mapcalc commands into one step to save on
intermediary maps,

> # count number of forest pixels value=1
> r.mapcalc "E = 1.0 * C"
> r.mapcalc "F = 1.0 * D"
> r.mapcalc "pf = (E/F)"

r.mapcalc << EOF
E = 1.0 * C
F = 1.0 * D
pf = (E/F)
EOF

not sure if it saves any time. for above, why not just
r.mapcalc "pf = (C*1.0/D)"
 ?


Hamish

[Jáchym Čepický-2]

Hallo,

I have only few comments to your script (look after "!!!" string in the
script):

#!/bin/sh

# Create pf map; Pf = (# forest pixels) / (# non-forest pixels)

# A=forest pixels
# B=non-forest pixels
# category vegation code
#1    Agriculture
#2    Brushland
#3    Clouds
#4    Forest
#5    Grassland
#6    River

# !!!
#%Module
#%  description: r.forest
#%End
#%option
#% key: input
#% type: string
#% description: raster input map
#% required : yes
#%end
#%option
#% key: output
#% type: string
#% description: raster output map
#% required : yes
#%end

if [ "$1" != "@ARGS_PARSED@" ] ; then
    exec $GISBASE/etc/bin/cmd/g.parser "$0" "$@"
fi

# !!! get map
g.copy rast=$GIS_OPT_input,$GIS_OPT_output

# !!! set null values
r.null map=$GIS_OPT_output setnull=3,6

#recode data
# !!! count number of forest pixels value=1
r.mapcalc "A=(if(1990veg == 4,1))+(if(1990veg == 1||1990veg ==
2||1990veg == 5,0))"

# count number of non-forest pixels value=1
r.mapcalc "B=if(1990veg == 1||1990veg == 2||1990veg == 4||1990veg == 5,1)"

# C number of forest pixels in a 5x5 window
# D number of non-forest pixels in a 5x5 window

r.neighbors input=A output=C method=sum size=5
r.neighbors input=B output=D method=sum size=5

## pixels with both forest from the cardinal directions
#create filters
# !!! could this not be done with help of r.mapcalc neighbord
# operations?
r.mfilter input=A output=Nf filter=filtrN.txt
r.mfilter input=A output=Sf filter=filtrS.txt
r.mfilter input=A output=Wf filter=filtrW.txt
r.mfilter input=A output=Ef filter=filtrE.txt

# !!! recode filtrd to 2=1 1=0 and recode
r.mapcalc "F3 = eval ( \
        Nf2 = if(Nf == 1||Nf == 0,0))+(if(Nf == 2,1)), \
        Sf2=(if(Sf == 1||Sf == 0,0))+(if(Sf == 2,1)), \
        Ef2=(if(Ef == 1||Ef == 0,0))+(if(Ef == 2,1)), \
        Wf2=(if(Wf == 1||Wf == 0,0))+(if(Wf == 2,1)), \
        Nf2 + Sf2 + Ef2 + Wf2))"

# !!! count both forest pixels
r.neighbors input=F4 output=F5 method=sum size=5

#frag model
#code    Category        Pf                    Pff
# 1     Patch            < 40%
# 2     Transitional    >= 40% but < 60%
# 3     Edge            > 60%                Pf - Pff > 0    1
# 4     Perforated        > 60%                Pf ??? Pff < 0    2
# 5     Interior        1
# 6        Undetermined    > 0.6                Pf = Pff        3

r.mapcalc "Pff3 = eval ( \
           Pff2 = (1.0 * C) /(1.0 * D)  - (1.0 * F5)/E, \
          (if (pff2 > 0,1)) + (if (pff2 < 0,2)) + (if (pff2 == 0,3)) \
          )"

#0-.4        1
#.4 - .6    2
#.6 - 0.99    3
#1            4

# !!!
r.mapcalc "Pff4 = (if ((F6/E) < 0.4,1)) + (if ((F6/E) >= 0.4 && (F6/E) < .6,2))
+ (if ((F6/E) >= 0.6 && (F6/E) < .99,3)) / (if ((F6/E) == 1,4))"

r.mapcalc "F11 = pff4 == 1"
r.mapcalc "F21 = pff4 == 2"
r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"
r.mapcalc "F51 = pf == 1"
r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

r.mapcalc "index = if (if (F11==1,1)) + (if (F21==1,2)) + (if
(F31==1,3)) + (if (F41==1,4)) + (if (F51==1,5)) + /
(if (F61==1,6))"

r.mapcalc "indexfin = if (index == 0 ||index == 7 ||index == 8,0) + if
(index == 1,1) + if /
(index == 2,2) + if (index == 3,3) + if (index == 4,4) + if (index ==
5,5) + if (index == 6,6)"


g.remove rast=A,B,C,D,E,Ef,Ef2,F,F3,F4,F5,F6,Nf,Nf2,Sf,Wf,Wf2,F11,F21,F31,F41,F51,F61,Pff2,Pff3,Pff4,Sf2,index

I made only few changes, so you could have an idea, how to rewrite your
script a bit, so it does not perform so many calculations - or so it
does not to create so many temporary maps. I did not modify the end of
the script.

Hamish allready wrote you, that the order of priorities is

1. It works
2. Clarity
--
3. Efficiency

Your script worked and now you can try to make it work better. When you
finish these optimalisations, try to answer following question: Is this
script fast enought for my purposes? If our answer is "no", then you
have to rewrite it (or at least some parts of it) to C. Before you do
so, you should get really working prototype.

>
> Cheers,
> Maning
>
> On 3/22/06, maning sambale <[hidden email]> wrote:
> > Thank you for the responses.  Will dig into to the manuals you recommended.
> > Will post later when I made some progress.  Right now i'm into the
> > generating the Pf which is simply r.neighbors and r.mapcalc commands .
> >
> > What is bugging me is generating Pff, which counts from the moving
> > window the forest pixel pairs in 4 cardinal directions.  (see
> > attached)
> >
> > Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
> > with at least one pixel forest)
> >
> > Cheers,
> > Maning
> >
> > PS.  how do I install r.li tarball in grass in cygwin?
> >
> > On 3/21/06, Jachym Cepicky <[hidden email]> wrote:
> > > Hallo,
> > > I have some experience with neighborhood operations on rastermaps, so I
> > > could help you to make the module together.
> > >
> > > Perhaps, you could be able make some prototype with r.mapcalc module?
> > > (see http://grass.itc.it/gdp/raster/mapcalc-algebra.pdf and
> > > http://grass.itc.it/gdp/raster/mapcalc.pdf)
> > >
> > > Jachym
> > >
> > > maning sambale wrote:
> > >
> > > >Dear list,
> > > >
> > > >I hope somebody can help me on this.  I am a MS student studying
> > > >forest loss and fragmentation here in the Philippines. I am interested
> > > >in implementing Mr. Kurt Ritters forest fragmentation model
> > > >(http://www.ecologyandsociety.org/vol4/iss2/art3/) and James Hurd's
> > > >adaptation to watersehd units
> > > >(http://resac.uconn.edu/publications/tech_papers/pdf_paper/Forest_Fragmentation_ASPRS2002.pdf)
> > > >in my study area and relate the fragmentation to socio-economic
> > > >drivers.
> > > >
> > > >Mr. Ritters gave me the C source code for his model, however, I am not
> > > >programmer so I don't know how to implement these.    Right now, I am
> > > >doing my image analysis (i.* modules) and classification using GRASS
> > > >and import the output to Arcview were EPA has an Attila extension
> > > >(http://www.epa.gov/esd/land-sci/attila/index.htm) for classifying
> > > >them into the forest fragmenation index.  What I want to do is create
> > > >a GRASS script for the forest fragmentation index.
> > > >
> > > >I have been using GRASS for just a couple of months and I am not very
> > > >familiar with most of the commands.  I am willing to collaborate with
> > > >anybody from this list doing related research on forest fragmentation.
> > > > As I have mentioned above, I am not a programmer thus my limitation
> > > >in understanding C code.  However, I am willing to do GRASS related
> > > >contributions related to my study.
> > > >
> > > >The method is described below:
> > > >
> > > >The basis for the forest fragmentation index was the fragmentation
> > > >model developed by Ritters et al. (2000) which was used and expanded
> > > >by Hurd et. al (2002) on Landsat TM data. The objective of the model
> > > >is to provide additional information beyond the commonly used measures
> > > >of forest proportions (i. e. percent of forest cover) by providing
> > > >information on the spatial configuration of fragmentation and
> > > >connectivity of forest within the study
> > > >region.
> > > >
> > > >The amount of forest and its occurrence is measured as adjacent forest
> > > >pixels within a 5 x 5 windows surrounding each forest pixel. This
> > > >information will be used to classify the window by the type of
> > > >fragmentation. The result was stored at the location of the center
> > > >pixel. Thus, a pixel value in the derived map refers to
> > > >"between-pixel" fragmentation around the corresponding forest location
> > > >(Riitters et al. 2000). The model generates two values that
> > > >characterize a forest pixel located at the center of a sliding window
> > > >of fixed size, Pf and Pff:
> > > >Pf is the proportion of pixels in the window that are forested. Pff
> > > >(strictly) is the proportion of all adjacent (cardinal directions
> > > >only) pixel pairs that include at least one forest pixel, for which
> > > >both pixels are forested. Pff (roughly) estimates the conditional
> > > >probability that, given a pixel of forest, its neighbor is also
> > > >forest.
> > > >
> > > >Pf = ( # of forest pixels) / (# of all non water pixels)
> > > >Pff = (# of pixel pairs with both pixel forest) / (# of pixel pairs
> > > >with at least one pixel forest)
> > > >
> > > >The classification model identifies six fragmentation categories: (1)
> > > >interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3) transitional,
> > > >0.4 < Pf < 0.6; (4) edge, Pf > 0.6 and Pf - Pff > 0; (5) perforated,
> > > >Pf > 0.6 and Pf â?????? Pff < 0, and (6) undetermined, Pf > 0.6 and Pf = Pff
> > > >
> > > >Any tips in implementing the generation of forest fragmentation index
> > > >either as a GRASS script or a step-by-step GRASS commands would be
> > > >very helpful (r.mapcalc syntax?).
> > > >
> > > >I have also attached Ritters's C code as well as my previous
> > > >correspondence to him.  Also attached are images from Hurd's study.
> > > >
> > > >Thank you very and looking forward for any response.
> > > >
> > > >Cheers,
> > > >
> > > >Maning Sambale
> > > >
> > > >
> > > >---------- Forwarded message ----------
> > > >From: Kurt H Riitters <[hidden email]>
> > > >Date: Oct 23, 2004 2:46 AM
> > > >Subject: Re: Developing a forest fragmentation index
> > > >To: emmanuel sambale <[hidden email]>
> > > >Cc: Kurt H Riitters <[hidden email]>
> > > >
> > > >
> > > >Emanuel,
> > > >  Attached is the main C program.  As I must say to everyone, I am not able
> > > >to provide technical support for using the program, and I do not guarantee
> > > >that it is bug-free.  You can share this with anyone you want to, but if
> > > >you share it, then you must answer any questions about it.
> > > >
> > > >  Probably the interesting part will be the algorithm for the moving
> > > >window, and the definition of the functions (like Pf and Pff).  If it were
> > > >me, I'd think about adapting those parts to new programs of your design.
> > > >
> > > >   The input / output format is assumed to be an old format that I adapted
> > > >from somewhere else...and it is not standard and is not used by anyone
> > > >else... so you will want to re-write the input / output routines to read
> > > >files in other formats.
> > > >
> > > >   I do not have any written documentation, but there are a lot of comments
> > > >in the code that will help understand the program flow.
> > > >
> > > >  Lastly, since this code is a research tool, there are many places
> > > >containing references to code or subroutines that are no longer included.
> > > >
> > > >    Oh, I should say... this program uses a moving window to calculate the
> > > >X and Y values that are in the model that appears in that 2000 paper.  One
> > > >run of the program makes a map of Pf; a second run makes a map of Pff.
> > > >Once you have the output maps of the Pf and Pff values, it is necessary to
> > > >put them together to perform the classification into categories like
> > > >'perforated' 'edge' etc.  Arc or envi would be convenient for that.
> > > >
> > > >  Best regards,
> > > >
> > > >Kurt Riitters
> > > >RWU-4803
> > > >Forest Health Monitoring
> > > >US Forest Service
> > > >3041 Cornwallis Road
> > > >Research Triangle Park NC 27709
> > > >919-549-4015
> > > >[hidden email]
> > > >
> > > >
> > > >
> > > >>>Emmanuel,
> > > >>> Thank you for the message.
> > > >>>  Of course you may use the methods; they are public knowledge.
> > > >>>
> > > >>>  I wrote C programs to process the landcover maps.  I am happy to send
> > > >>>them to you, however, if you are not a programmer then it will be very
> > > >>>difficult to adapt the programs to your problems.  Even another
> > > >>>
> > > >>>
> > > >programmer
> > > >
> > > >
> > > >>>may find it difficult since I am not a very good programmer either.
> > > >>>
> > > >>> I do not have any arc or envi scripts since I use the C programs...
> > > >>>
> > > >>>
> > > >But
> > > >
> > > >
> > > >>>it is a very simple moving window (arc: focal function) approach that
> > > >>>should be easy to implement in another language...  For example, arc
> > > >>>
> > > >>>
> > > >should
> > > >
> > > >
> > > >>>be able to do the "Pf" part of the model with the Grid: focalsum
> > > >>>
> > >
> > >
> >
> >
> > --
> > |---------|----------------------------------------------------------|
> > | __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
> > | '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
> > |  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
> > |  |  /T  |Philippine Biodiversity Web Map                           |
> > | _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
> > |---------|----------------------------------------------------------|
> >
> >
>
>
> --
> |---------|----------------------------------------------------------|
> | __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
> | '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
> |  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
> |  |  /T  |Philippine Biodiversity Web Map                           |
> | _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
> |---------|----------------------------------------------------------|
>
>
> --
> |---------|----------------------------------------------------------|
> | __.-._  |"Ohhh. Great warrior. Wars not make one great." -Yoda     |
> | '-._"7' |"Freedom is still the most radical idea of all" -N.Branden|
> |  /'.-c  |Linux registered user #402901, http://counter.li.org/     |
> |  |  /T  |Philippine Biodiversity Web Map                           |
> | _)_/LI  |www.geocities.com/esambale/philbiodivmap/philbirds.html   |
> |---------|----------------------------------------------------------|
>

--
Jachym Cepicky
e-mail: [hidden email]
URL: http://les-ejk.cz

[maning sambale]

Dear grasslist,

I have made some revisions from my forest fragmentation script.  Thank
you for your contributions.  I have incorporated some of your
suggestions (those that I understood) as well as adding options for
specifying window size, removal of temporary maps and displaying
summary class statistics.

One of the problems I encounter is the processing of very large maps.
My P4 CPU clock shoots to 100% and then crashes.  I was thinking maybe
of cutting the images to manageable chunks and then process them by
batches.

Next is to overlay watershed basins and create summary statistics for
each class for each basin and import the class statistics to R
readable format.  For further statistical analysis.

Thanks for your ideas!

Maning

r.forest.frag script

#!/bin/sh
############################################################################
#
# MODULE:       r.forestfrag
#
# AUTHOR(S):    Emmanuel Sambale [hidden email] [hidden email]
#
# PURPOSE:      Creates forest fragmentation index map from a raster
vegetation file
#                The index map was based on Riitters, K., J. Wickham,
R. O'Neill, B. Jones,
#                and E. Smith. 2000. Global-scale patterns of forest
fragmentation. Conservation
#                Ecology 4(2): 3. [online] URL:
http://www.consecol.org/vol4/iss2/art3/
#
# COPYRIGHT:    (C) 2005 by the GRASS Development Team
#
#               This program is free software under the GNU General Public
#               License (>=v2). Read the file COPYING that comes with GRASS
#               for details.
#
#############################################################################

#%Module
#% description: creates forest fragmentation index from a GRASS raster map
#%End
#%flag
#%  key: r
#%  description: remove temporary maps
#%END
#%option
#% key: input
#% type: string
#% gisprompt: old,cell,raster
#% description: Name of vegetation map
#% required : yes
#%End
#%option
#% key: window
#% type: integer
#% description: window size default is 5
#% answer : 5
#% required : yes
#%END

if  [ -z "$GISBASE" ] ; then
        echo "You must be in GRASS GIS to run this program."
        exit 1
fi

if [ "$1" != "@ARGS_PARSED@" ] ; then
   exec $GISBASE/bin/g.parser "$0" "$@"
fi

unset LC_ALL
export LC_NUMERIC=C

#set to current input map region
g.region rast=$GIS_OPT_input -p

# get map
r.mapcalc 1990veg2=$GIS_OPT_input

# set null values
r.null map=1990veg2 setnull=3,6

#recode data where nonforest = 0 & forest = 1
r.mapcalc "A=(if(1990veg2 == 4,1))+(if(1990veg2 == 1||1990veg2 ==
2||1990veg2 == 5,0))"

echo "computing pf values ..."
# count number of non-forest pixels value=1
r.mapcalc "B=if(1990veg2 == 1||1990veg2 == 2||1990veg2 == 4||1990veg2 == 5,1)"

# C number of forest pixels in a 5x5 window
# D number of non-forest pixels in a 5x5 window

r.neighbors input=A output=C method=sum size="$GIS_OPT_window"
r.neighbors input=B output=D method=sum size="$GIS_OPT_window"

# count number of forest pixels value=1

# this one is better create pf map
r.mapcalc << EOF
E = 1.0 * C
F = 1.0 * D
pf = (E/F)
EOF

echo "computing pff values ..."
## pixels with both forest from the cardinal directions
# 0--x--0
# |  |  |
# x--x--x
# |  |  |
# 0--x--0

r.mapcalc "F4 = (A * A [1,0]) * (A * A [-1,0]) * (A * A [0,-1]) * (A * A [0,1])"

#count both forest pixels
r.neighbors input=F4 output=F5 method=sum size="$GIS_OPT_window"

# create pff map
r.mapcalc << EOF
F6 = 1.0 * F5
pff = (F6/E)
EOF

echo "computing fragmentation index ..."
# (1) interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3)
transitional, 0.4 < Pf < 0.6; (4) edge,
# Pf > 0.6 and Pf - Pff > 0; (5) perforated, Pf > 0.6 and Pf – Pff <
0, and (6) undetermined,
# Pf > 0.6 and Pf = Pffr.mapcalc "Pff2 = pf - pff"

r.mapcalc "Pff3 = (if (pff2 > 0,1)) + (if (pff2 < 0,2)) + (if (pff2 == 0,3))"

# recode pff to 0-.4, .4-.6, .6-.99, 1
r.mapcalc "Pff4 = (if (pff < 0.4,1)) + (if (pff >= 0.4 && pff < .6,2))
+ (if (pff >= 0.6 && pff < .99,3)) + (if (pff == 1,4))"


r.mapcalc "F11 = pff4 == 1"
r.mapcalc "F21 = pff4 == 2"
r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"
r.mapcalc "F51 = pf == 1"
r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

r.mapcalc "index = if (if (F11==1,1)) + (if (F21==1,2)) + (if
(F31==1,3)) + (if (F41==1,4)) + (if (F51==1,5)) + (if (F61==1,6))"

r.mapcalc "indexfin2= ( A * index )"

#create color codes

echo "creating color codes and categories ..."
r.colors indexfin2 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
7 145:207:97
8 145:207:98
EOF

#create categories
cat recl.txt
cat recl.txt | r.reclass indexfin2 out=indexfin3 title="frag index"

r.mapcalc indexfin4=indexfin3

if [ $GIS_FLAG_r -eq 1 ] ; then
  echo "Temporary files deleted ...."
  g.remove rast=A,B,C,D,E,F,F11,F21,F31,F4,F41,F5,F51,F6,F61,Pff2,Pff3,Pff4
  g.remove rast=indexfin3
  g.remove rast=indexfin2
  g.remove rast=index
fi

#create color codes
r.colors indexfin4 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
7 145:207:97
8 145:207:98
EOF

echo "generate map reports ..."
r.report map=indexfin4 units=h,p null=* nsteps=255 -e -i

reclass table
recl.txt
0 = 0 nonforest
1 = 1 patch
2 = 2 transitional
3 = 3 edge
4 = 4 perforated
5 = 5 interior
6 thru 8 = 6 undet

[Jáchym Čepický-2]

Hallo,

On Mon, Apr 10, 2006 at 06:10:33PM +0800, maning sambale wrote:

> Dear grasslist,
>
> I have made some revisions from my forest fragmentation script.  Thank
> you for your contributions.  I have incorporated some of your
> suggestions (those that I understood) as well as adding options for
> specifying window size, removal of temporary maps and displaying
> summary class statistics.
>
> One of the problems I encounter is the processing of very large maps.
> My P4 CPU clock shoots to 100% and then crashes.  I was thinking maybe
> of cutting the images to manageable chunks and then process them by
> batches.
>

If you are using mostly r.mapcalc, it should not crashe - in which step
does this happen?

some operations should be performed faster, if you would use g.copy and
r.reclass instead of r.mapcalc -- look after "!!!" string below

> Next is to overlay watershed basins and create summary statistics for
> each class for each basin and import the class statistics to R
> readable format.  For further statistical analysis.
>
> Thanks for your ideas!
>
> Maning
>
> r.forest.frag script

Jachym

>
#!/bin/sh
############################################################################
#
# MODULE:       r.forestfrag
#
# AUTHOR(S):    Emmanuel Sambale [hidden email] [hidden email]
#
# PURPOSE:      Creates forest fragmentation index map from a raster
vegetation file
#                The index map was based on Riitters, K., J. Wickham,
R. O'Neill, B. Jones,
#                and E. Smith. 2000. Global-scale patterns of forest
fragmentation. Conservation
#                Ecology 4(2): 3. [online] URL:
http://www.consecol.org/vol4/iss2/art3/
#
# COPYRIGHT:    (C) 2005 by the GRASS Development Team
#
#               This program is free software under the GNU General Public
#               License (>=v2). Read the file COPYING that comes with GRASS
#               for details.
#
#############################################################################

#%Module
#% description: creates forest fragmentation index from a GRASS raster map
#%End
#%flag
#%  key: r
#%  description: remove temporary maps
#%END
#%option
#% key: input
#% type: string
#% gisprompt: old,cell,raster
#% description: Name of vegetation map
#% required : yes
#%End
#%option
#% key: window
#% type: integer
#% description: window size default is 5
#% answer : 5
#% required : yes
#%END

if  [ -z "$GISBASE" ] ; then
        echo "You must be in GRASS GIS to run this program."
        exit 1
fi

if [ "$1" != "@ARGS_PARSED@" ] ; then
   exec $GISBASE/bin/g.parser "$0" "$@"
fi

unset LC_ALL
export LC_NUMERIC=C

#set to current input map region
g.region rast=$GIS_OPT_input -p

# get map !!! r.mapcalc -> g.copy
g.copy  rast=$GIS_OPT_input,1990veg2

# set null values
r.null map=1990veg2 setnull=3,6

#recode data where nonforest = 0 & forest = 1
# !!! r.mapcalc -> r.reclass
#r.mapcalc "A=(if(1990veg2 == 4,1))+(if(1990veg2 == 1||1990veg2 ==
2||1990veg2 == 5,0))"
r.reclass in=1990veg2 out=A << EOF
4 = 1
2 = 5
2 = 5
end
EOF

echo "computing pf values ..."
# count number of non-forest pixels value=1
# r.mapcalc -> r.reclass !!!
#r.mapcalc "B=if(1990veg2 == 1||1990veg2 == 2||1990veg2 == 4||1990veg2 == 5,1)"
r.reclass in=1990veg2 out=B << EOF
1 = 1
2 = 1
4 = 1
5 = 1
end
EOF


# C number of forest pixels in a 5x5 window
# D number of non-forest pixels in a 5x5 window

r.neighbors input=A output=C method=sum size="$GIS_OPT_window"
r.neighbors input=B output=D method=sum size="$GIS_OPT_window"

# count number of forest pixels value=1

# this one is better create pf map
r.mapcalc << EOF
E = 1.0 * C
F = 1.0 * D
pf = (E/F)
EOF

echo "computing pff values ..."
## pixels with both forest from the cardinal directions
# 0--x--0
# |  |  |
# x--x--x
# |  |  |
# 0--x--0

r.mapcalc "F4 = (A * A [1,0]) * (A * A [-1,0]) * (A * A [0,-1]) * (A * A [0,1])"

#count both forest pixels
r.neighbors input=F4 output=F5 method=sum size="$GIS_OPT_window"

# create pff map
r.mapcalc << EOF
F6 = 1.0 * F5
pff = (F6/E)
EOF

echo "computing fragmentation index ..."
# (1) interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3) transitional, 0.4 < Pf < 0.6; (4) edge,
# Pf > 0.6 and Pf - Pff > 0; (5) perforated, Pf > 0.6 and Pf ??? Pff < 0, and (6) undetermined,
# Pf > 0.6 and Pf = Pffr.mapcalc "Pff2 = pf - pff"

r.mapcalc "Pff3 = (if (pff2 > 0,1)) + (if (pff2 < 0,2)) + (if (pff2 == 0,3))"

# recode pff to 0-.4, .4-.6, .6-.99, 1
r.mapcalc "Pff4 = (if (pff < 0.4,1)) + (if (pff >= 0.4 && pff < .6,2))
+ (if (pff >= 0.6 && pff < .99,3)) + (if (pff == 1,4))"


# !!! r.mapcalc -> r.reclass
# r.mapcalc "F11 = pff4 == 1"
# r.mapcalc "F21 = pff4 == 2"
# r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
# r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"
# r.mapcalc "F51 = pf == 1"
# r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

r.reclass in=pff4 out=F11 << EOF
1 = 1
end
EOF

r.reclass in=pff4 out=F21 << EOF
2 = 1
end
EOF

r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"

r.reclass in=pf out=F51 << EOF
1 = 1
end
EOF

r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

# !!! would this work too?
# r.mapcalc "index = if (if (F11==1,1)) + (if (F21==1,2)) + (if (F31==1,3)) + (if (F41==1,4)) + (if (F51==1,5)) + (if (F61==1,6))"
r.mapcalc "index = F11 + F21*2 + F31*3 + F41*4 + F51*5 + F61*6"

r.mapcalc "indexfin2= ( A * index )"

#create color codes

echo "creating color codes and categories ..."
r.colors indexfin2 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
7 145:207:97
8 145:207:98
EOF

#create categories
cat recl.txt
cat recl.txt | r.reclass indexfin2 out=indexfin3 title="frag index"

# r.mapcalc -> g.copy !!!
g.copy indexfin3,indexfin4

if [ $GIS_FLAG_r -eq 1 ] ; then
  echo "Temporary files deleted ...."
  g.remove rast=A,B,C,D,E,F,F11,F21,F31,F4,F41,F5,F51,F6,F61,Pff2,Pff3,Pff4
  g.remove rast=indexfin3
  g.remove rast=indexfin2
  g.remove rast=index
fi

#create color codes
r.colors indexfin4 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
7 145:207:97
8 145:207:98
EOF

echo "generate map reports ..."
r.report map=indexfin4 units=h,p null=* nsteps=255 -e -i

reclass table
recl.txt
0 = 0 nonforest
1 = 1 patch
2 = 2 transitional
3 = 3 edge
4 = 4 perforated
5 = 5 interior
6 thru 8 = 6 undet

--
Jachym Cepicky                                  
e-mail: [hidden email]          
URL: http://les-ejk.cz                     
GPG: http://les-ejk.cz/gnupg_public_key/   
-----------------------------------------  
OFFICE:                                    
Department of Geoinformation Technologies  
LDF MZLU v Brně                            
Zemědělská 3                                
613 00 Brno                                
e-mail: [hidden email]            
URL:    http://mapserver.mendelu.cz
Tel.:   [phone]

[HamishB]

Maning:
> > One of the problems I encounter is the processing of very large
> > maps.  My P4 CPU clock shoots to 100% and then crashes.  I was
> > thinking maybe of cutting the images to manageable chunks and then
> > process them by batches.
Jachym:
> If you are using mostly r.mapcalc, it should not crashe - in which
> step does this happen?


To help debug change first line to

!/bin/sh -x

to see shell script command info as it runs.
Then you can see when it breaks.

CPU to 100% is normal. "Crashing" (how?) usually happens when the system
runs out of memory and the kernel kills the process that is causing all
the trouble. r.mapcalc & r.neighbors only process little chunks
(row/rows) of the data at a time so it is unusual that they'd be
responsible.

Another possibility is a hardware problem where your processor starts to
overheat once it has been running full bore for a while. More random
system crashes rather than a single module crash in the same place for
that case.



Hamish

[Stefan Sylla]

Hi,

has this script ever been implemented into current grass releases?

However, I've been using this script right now and it helped me a lot to get the algorithm from Riitters et al. running, so thanks a lot for your work, Manning.

I made some changes to the script, because I did not understand some parts of it. If there is anybody who would be interested to discuss I'd be happy.

For example, I replaced the part

...
echo "computing fragmentation index ..."
# (1) interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3)
# transitional, 0.4 < Pf < 0.6; (4) edge,
# Pf > 0.6 and Pf - Pff > 0; (5) perforated, Pf > 0.6 and Pf – Pff < 0, and (6) undetermined, Pf > 0.6 and Pf = Pff

r.mapcalc "pff2 = pf - pff"

r.mapcalc "pff3 = (if (pff2 > 0,1)) + (if (pff2 < 0,2)) + (if (pff2 == 0,3))"

# recode pff to 0-.4, .4-.6, .6-.99, 1
r.mapcalc "pff4 = (if (pff < 0.4,1)) + (if (pff >= 0.4 && pff < .6,2)) + (if (pff >= 0.6 && pff < .99,3)) + (if (pff == 1,4))"


r.mapcalc "F11 = pff4 == 1"
r.mapcalc "F21 = pff4 == 2"
r.mapcalc "F31 = (pff4 == 3) && (pff3 == 1)"
r.mapcalc "F41 = (pff4 == 3) && (pff3 == 2)"
r.mapcalc "F51 = pf == 1"
r.mapcalc "F61 = (pff4 == 3) && (pff3 == 3)"

r.mapcalc "index = (if (F11==1,1)) + (if (F21==1,2)) + (if (F31==1,3)) + (if (F41==1,4)) + (if (F51==1,5)) + (if (F61==1,6))"
...

with

...
echo "computing fragmentation index ..."
# (1) interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3)
# transitional, 0.4 < Pf < 0.6; (4) edge,
# Pf > 0.6 and Pf - Pff < 0; (5) perforated, Pf > 0.6 and Pf – Pff > 0, and (6) undetermined, Pf > 0.6 and Pf = Pff

r.mapcalc "pf2 = pf - pff"


r.mapcalc "f1=if(pf<0.4,1,0)" # patch
r.mapcalc "f2=if(pf>=0.4 && pf<0.6,2,0)" # transitional
r.mapcalc "f3=if(pf>=0.6 && pf2<0,3,0)" # edge
r.mapcalc "f4=if(pf>0.6 && pf2>0,4,0)" # perforate
r.mapcalc "f5=if(pf==1 && pff==1,5,0)" # interior
r.mapcalc "f6=if(pf>0.6 && pf<1 && pf2==0,6,0)" # undetermined

r.mapcalc "index=f1+f2+f3+f4+f5+f6"
...

... because I could not figure out how the formulas for each fraction type (interior, patch, etc.) where implemented in the code; so I tried to get it a little clearer with the code above (please comment, because the results are different!). Furthermore, I discovered that there where some errata in the original writings from Riitters et al. for edge and perforated forest, so I also corrected for that.

If there's anybody out there to further discuss the script this would be great!

Stefan

----------------------------------------------------------------------------------------------------------------------
P.S. Here's the complete script that I am currently using and which is running locally on my machine:

#!/bin/sh
############################################################################
#
# MODULE:       r.forestfrag
#
# AUTHOR(S):    Emmanuel Sambale [hidden email] [hidden email]
#
# PURPOSE:      Creates forest fragmentation index map from a raster vegetation file
#                The index map was based on Riitters, K., J. Wickham, R. O'Neill, B. Jones,
#                and E. Smith. 2000. Global-scale patterns of forest fragmentation. Conservation
#                Ecology 4(2): 3. [online] URL: http://www.consecol.org/vol4/iss2/art3/
#
# COPYRIGHT:    (C) 2005 by the GRASS Development Team
#
#               This program is free software under the GNU General Public
#               License (>=v2). Read the file COPYING that comes with GRASS
#               for details.
#
#############################################################################

#%Module
#% description: creates forest fragmentation index from a GRASS raster map
#%End
#%flag
#%  key: r
#%  description: remove temporary maps
#%END
#%option
#% key: input
#% type: string
#% gisprompt: old,cell,raster
#% description: Name of vegetation map
#% required : yes
#%End
#%option
#% key: window
#% type: integer
#% description: window size default is 5
#% answer : 5
#% required : yes
#%END

if  [ -z "$GISBASE" ] ; then
        echo "You must be in GRASS GIS to run this program."
        exit 1
fi

#if [ "$1" != "@ARGS_PARSED@" ] ; then
#   exec $GISBASE/bin/g.parser "$0" "$@"
#fi

unset LC_ALL
export LC_NUMERIC=C

#set to current input map region
g.region rast=$GIS_OPT_input -p

# get map (assuming 'landclass96 is a raster with forest/non-forest categories)
r.mapcalc 1990veg2=$GIS_OPT_input

# set null values (water areas will not be included)
#r.null map=1990veg2 setnull=6

#recode data: nonforest = 0, forest = 1
r.mapcalc "A=(if(1990veg2==1,1,0))"

echo "computing pf values ..."
# create grid with all values = 1 (needed to calculate pf):
r.mapcalc "B=if(1990veg2,1,1)"

# C: number of forest pixels within a 3x3 moving window
# D: create grid with value=number of pixels within 3x3 moving window

r.neighbors input=A output=C method=sum size=3 --o
r.neighbors input=B output=D method=sum size=3 --o

# create pf map
r.mapcalc << EOF
E = 1.0 * C
F = 1.0 * D
pf = (E/F)
EOF

## generate grid where value=1 if forest pixel is surrounded by forest pixels in cardinal directions:
# 0--x--0
# |  |  |
# x--x--x
# |  |  |
# 0--x--0

# 3x3 moving window
r.mapcalc "F4 = (A * A [1,0]) * (A * A [-1,0]) * (A * A [0,-1]) * (A * A [0,1])"

# create grid with pixel-value = number of pixels completely surrounded by forest within 3x3 moving window:
r.neighbors input=F4 output=F5 method=sum size=3 --o

# create pff map
r.mapcalc << EOF
F6 = 1.0 * F5
pff = (F6/E)
EOF

echo "computing fragmentation index ..."
# (1) interior, for which Pf = 1.0; (2) patch, Pf < 0.4; (3)
# transitional, 0.4 < Pf < 0.6; (4) edge,
# Pf > 0.6 and Pf - Pff < 0; (5) perforated, Pf > 0.6 and Pf – Pff > 0, and (6) undetermined, Pf > 0.6 and Pf = Pff

r.mapcalc "pf2 = pf - pff"

r.mapcalc "f1=if(pf<0.4,1,0)" # patch
r.mapcalc "f2=if(pf>=0.4 && pf<0.6,2,0)" # transitional
r.mapcalc "f3=if(pf>=0.6 && pf2<0,3,0)" # edge
r.mapcalc "f4=if(pf>0.6 && pf2>0,4,0)" # perforate
r.mapcalc "f5=if(pf==1 && pff==1,5,0)" # interior
r.mapcalc "f6=if(pf>0.6 && pf<1 && pf2==0,6,0)" # undetermined

r.mapcalc "index=f1+f2+f3+f4+f5+f6"

r.mapcalc "indexfin2= ( A * index )"

r.mapcalc "indexfin2= ( A * index )"

#create color codes

echo "creating color codes and categories ..."
r.colors indexfin2 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
EOF

#create categories
cat recl.txt | r.reclass indexfin2 out=indexfin3 title="frag index"

r.mapcalc indexfin4=indexfin3

echo "Temporary files deleted ...."
  g.remove rast=A,B,C,D,E,F,f1,f2,f3,f4,f5,f6,pf2
  g.remove rast=indexfin3
  g.remove rast=indexfin2
  g.remove rast=index

#create color codes
r.colors indexfin4 color=rules << EOF
0 yellow
1 215:48:39
2 252:141:89
3 254:224:139
4 217:239:139
5 26:152:80
6 145:207:96
EOF

echo "generate map reports ..."
r.report map=indexfin4 units=h,p null=* nsteps=255 -e -i -n
d.rast.leg indexfin4

## create reclass table "recl.txt":
#cat > recl2.txt
#0 = 0 nonforest
#1 = 1 patch
#2 = 2 transitional
#3 = 3 edge
#4 = 4 perforated
#5 = 5 interior
# 6 thru 8 = 6 undef

[Markus Neteler]

On Sun, Feb 5, 2012 at 2:26 PM, Stefan Sylla <[hidden email]> wrote:

> Hi,
>
> has this script ever been implemented into current grass releases?
>
> However, I've been using this script right now and it helped me a lot to get
> the algorithm from Riitters et al. running, so thanks a lot for your work,
> Manning.
>
> I made some changes to the script, because I did not understand some parts
> of it. If there is anybody who would be interested to discuss I'd be happy.

...

In any case, a good place to maintain the script would be the
GRASS-Addons SVN, see

http://grass.osgeo.org/wiki/Addons

Markus
_______________________________________________
grass-user mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-user

[anusheema]

CONTENTS DELETED

The author has deleted this message.

[Markus Neteler]

On Sat, Dec 27, 2014 at 12:53 PM, anusheema <[hidden email]> wrote:
> Dear list members,
>
> I am also interested in implementing Dr. Kurt Ritters forest fragmentation
> model. I was, however, attempting this in R Statistical software, although I
> am pretty new to both R and GRASS.

Which GRASS version do you use?

> While reading forums on this model, I came across this post. I see that Mr.
> Maning Sambale and Mr. Stefan Sylla have developed an AddOn for GRASS:
> r.forestfrag (http://grasswiki.osgeo.org/wiki/GRASS_AddOns#r.forestfrag).
> But while trying to attempt this model on GRASS, it keeps showing an error
> (AddOn not found).

If you are using GRASS 6: the script has been written for GRASS 7. You
can install both 6 and 7 in parallel on the same computer. However, I
see a link to a GRASS 6 version of the script in the Wiki page.

If you are using GRASS 7: The reason is that it has not been
implemented as Python script but as Shell script. To my knowledge
g.extention (the extension manager) is not capable to deal with shell
scripts.

> If you could please share the script? I was trying to implement the script
> mentioned in this post, but without much luck. In R, with the help of this
> script, replicating Pf is easier, however, for Pff, implementing the code in
> cardinal direction only, is more complex.

You can simply download the script from here:
http://trac.osgeo.org/grass/browser/grass-addons/grass7/raster/r.forestfrag/r.forestfrag?format=txt
and put it into the search path in your GRASS session. Be sure to set
the "executable" flag.

Please let us know how it goes.

Markus
_______________________________________________
grass-user mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-user

[pvanbosgeo]

Hi Anusheema

As I replied on your question on my blogpost, there are two versions of the script; one for GRASS 6  and one for GRASS 7. In both cases you need to install the script manually (I am planning to rewrite the script as a python script, but that will take some time).

The script for GRASS 7 did indeed not work anymore for me due to syntax changes in some of the functions used. An updated version should be available on http://trac.osgeo.org/grass/browser/grass-addons/grass7/raster/r.forestfrag. If you try out and it doesn't work, let me know (and if so, please indicate whether you are running this on Linux or Windows).

Paulo

On Sat, Dec 27, 2014 at 2:38 PM, Markus Neteler <[hidden email]> wrote:

On Sat, Dec 27, 2014 at 12:53 PM, anusheema <[hidden email]> wrote:
> Dear list members,
>
> I am also interested in implementing Dr. Kurt Ritters forest fragmentation
> model. I was, however, attempting this in R Statistical software, although I
> am pretty new to both R and GRASS.

Which GRASS version do you use?

> While reading forums on this model, I came across this post. I see that Mr.
> Maning Sambale and Mr. Stefan Sylla have developed an AddOn for GRASS:
> r.forestfrag (http://grasswiki.osgeo.org/wiki/GRASS_AddOns#r.forestfrag).
> But while trying to attempt this model on GRASS, it keeps showing an error
> (AddOn not found).

If you are using GRASS 6: the script has been written for GRASS 7. You
can install both 6 and 7 in parallel on the same computer. However, I
see a link to a GRASS 6 version of the script in the Wiki page.

If you are using GRASS 7: The reason is that it has not been
implemented as Python script but as Shell script. To my knowledge
g.extention (the extension manager) is not capable to deal with shell
scripts.

> If you could please share the script? I was trying to implement the script
> mentioned in this post, but without much luck. In R, with the help of this
> script, replicating Pf is easier, however, for Pff, implementing the code in
> cardinal direction only, is more complex.

You can simply download the script from here:
http://trac.osgeo.org/grass/browser/grass-addons/grass7/raster/r.forestfrag/r.forestfrag?format=txt
and put it into the search path in your GRASS session. Be sure to set
the "executable" flag.

Please let us know how it goes.

Markus

_______________________________________________
grass-user mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-user

_______________________________________________
grass-user mailing list
[hidden email]
http://lists.osgeo.org/mailman/listinfo/grass-user

[anusheema]

CONTENTS DELETED

The author has deleted this message.