/* *$Id: patch.c,v 2.3.4.1 2006/10/06 04:00:08 hamish Exp $ ************************************************************ * MODULE: r.le.patch/patch.c * * Version 5.0 Nov. 1, 2001 * * * * AUTHOR: W.L. Baker, University of Wyoming * * BAKERWL@UWYO.EDU * * * * PURPOSE: To analyze attributes of patches in a landscape * * patch.c computes and saves the patch measures * * * * COPYRIGHT: (C) 2001 by W.L. Baker * * * * This program is free software under the GNU General * * Public License(>=v2). Read the file COPYING that comes * * with GRASS for details * * * ************************************************************/ #include #include #include "patch.h" extern struct CHOICE *choice; extern int total_patches, ntype, n_scale, n_unit, *recl_count; extern float *shape_PA, *shape_CPA, *shape_RCC, **recl_tb, *size_cl; FILE *a5,*a6,*a1_4,*a7,*a8,*c1_4,*c5,*c6,*c7,*c8,*c9c,*c9e,*c10c, *c10e,*s3,*s4,*s5,*s6,*s7_8,*s1_2,*h3,*h4,*h5,*h6,*h1_2,*n1_4, *p4,*p5,*p6,*p1_3, *outfile; int i; /************************************/ /* RUN THE DEFAULT PATCH MEASURES */ /************************************/ void df_patch (PATCH *patch_list) { PATCH *tmp = patch_list; int *type_dens, type_coh; char path[30]; if(!total_patches) return; /* Write the scale and unit in each file */ /***********************************************************/ if (choice->att[1] || choice->att[2] || choice->att[3] || choice->att[4]) { a1_4 = fopen0("r.le.out/a1-4.out", "a"); fprintf(a1_4, "%5d %5d", n_scale, n_unit); } if (choice->att[5]) { a5 = fopen0("r.le.out/a5.out", "a"); fprintf(a5, "%5d %5d", n_scale, n_unit); } if (choice->att[6]) { a6 = fopen0("r.le.out/a6.out", "a"); fprintf(a6, "%5d %5d", n_scale, n_unit); } if (choice->att[7]) { a7 = fopen0("r.le.out/a7.out", "a"); fprintf(a7, "%5d %5d", n_scale, n_unit); } if (choice->att[8]) { a8 = fopen0("r.le.out/a8.out", "a"); fprintf(a8, "%5d %5d", n_scale, n_unit); } /***********************************************************/ if (choice->size[1] || choice->size[2]) { s1_2 = fopen0("r.le.out/s1-2.out", "a"); fprintf(s1_2, "%5d %5d", n_scale, n_unit); } if (choice->size[3]) { s3 = fopen0("r.le.out/s3.out", "a"); fprintf(s3, "%5d %5d", n_scale, n_unit); } if (choice->size[4]) { s4 = fopen0("r.le.out/s4.out", "a"); fprintf(s4, "%5d %5d", n_scale, n_unit); } if (choice->size[5]) { s5 = fopen0("r.le.out/s5.out", "a"); fprintf(s5, "%5d %5d", n_scale, n_unit); } if (choice->size[6]) { s6 = fopen0("r.le.out/s6.out", "a"); fprintf(s6, "%5d %5d\n", n_scale, n_unit); } if (choice->size[7] || choice->size[8]) { s7_8 = fopen0("r.le.out/s7-8.out", "a"); fprintf(s7_8, "%5d %5d", n_scale, n_unit); } /***********************************************************/ if (choice->core[1] || choice->core[2] || choice->core[3] || choice->core[4]) { c1_4 = fopen0("r.le.out/c1-4.out", "a"); fprintf(c1_4, "%5d %5d", n_scale, n_unit); } if (choice->core[5]) { c5 = fopen0("r.le.out/c5.out", "a"); fprintf(c5, "%5d %5d", n_scale, n_unit); } if (choice->core[6]) { c6 = fopen0("r.le.out/c6.out", "a"); fprintf(c6, "%5d %5d", n_scale, n_unit); } if (choice->core[7]) { c7 = fopen0("r.le.out/c7.out", "a"); fprintf(c7, "%5d %5d", n_scale, n_unit); } if (choice->core[8]) { c8 = fopen0("r.le.out/c8.out", "a"); fprintf(c8, "%5d %5d", n_scale, n_unit); } if (choice->core[9]) { c9c = fopen0("r.le.out/c9c.out", "a"); fprintf(c9c, "%5d %5d", n_scale, n_unit); c9e = fopen0("r.le.out/c9e.out", "a"); fprintf(c9e, "%5d %5d", n_scale, n_unit); } if (choice->core[10]) { c10c = fopen0("r.le.out/c10c.out", "a"); fprintf(c10c, "%5d %5d\n", n_scale, n_unit); c10e = fopen0("r.le.out/c10e.out", "a"); fprintf(c10e, "%5d %5d\n", n_scale, n_unit); } /************************************************************/ if (choice->shape[1] || choice->shape[2]) { h1_2 = fopen0("r.le.out/h1-2.out", "a"); fprintf(h1_2, "%5d %5d", n_scale, n_unit); } if (choice->shape[3]) { h3 = fopen0("r.le.out/h3.out", "a"); fprintf(h3, "%5d %5d", n_scale, n_unit); } if (choice->shape[4]) { h4 = fopen0("r.le.out/h4.out", "a"); fprintf(h4, "%5d %5d", n_scale, n_unit); } if (choice->shape[5]) { h5 = fopen0("r.le.out/h5.out", "a"); fprintf(h5, "%5d %5d", n_scale, n_unit); } if (choice->shape[6]) { h6 = fopen0("r.le.out/h6.out", "a"); fprintf(h6, "%5d %5d\n", n_scale, n_unit); } /************************************************************/ if (choice->boundary[1] || choice->boundary[2] || choice->boundary[3] || choice->boundary[4]) { n1_4 = fopen0("r.le.out/n1-4.out", "a"); fprintf(n1_4, "%5d %5d", n_scale, n_unit); } /************************************************************/ if (choice->perim[1] || choice->perim[2] || choice->perim[3]) { p1_3 = fopen0("r.le.out/p1-3.out", "a"); fprintf(p1_3, "%5d %5d", n_scale, n_unit); } if (choice->perim[4]) { p4 = fopen0("r.le.out/p4.out", "a"); fprintf(p4, "%5d %5d", n_scale, n_unit); } if (choice->perim[5]) { p5 = fopen0("r.le.out/p5.out", "a"); fprintf(p5, "%5d %5d", n_scale, n_unit); } if (choice->perim[6]) { p6 = fopen0("r.le.out/p6.out", "a"); fprintf(p6, "%5d %5d", n_scale, n_unit); } type_dens = (int *)G_calloc(25, sizeof(int)); for (i = 0; i < 25; i++) type_dens[i] = 0; /* for each patch on the patch list */ while(tmp) { if ((type_coh = recl_coh(tmp->att)) >= 0) type_dens[type_coh] ++; if (choice->att[0]) df_att(tmp, type_coh, type_dens); if (choice->core[0]) df_core(tmp, type_coh, type_dens); if (choice->size[0]) df_size(tmp, type_coh, type_dens); if (choice->shape[0]) df_shape(tmp, type_coh, type_dens); if (choice->perim[0]) df_perim(tmp, type_coh, type_dens); if (choice->boundary[0]) df_boundary(tmp); if (strcmp(choice->out,"") && choice->wrum != 'm') { sprintf(path, "r.le.out/%s", choice->out); outfile = fopen0(path, "a"); fprintf(outfile, "%3d %3d %6d %7.1f %4.0f %4.0f %8.0f \ %8.0f %8.0f %8.0f %6.3f %6.3f %6.3f %8d %4.3f\n", n_scale, n_unit, tmp->num, tmp->att, tmp->c_row, tmp->c_col, tmp->area, tmp->core, tmp->edge, tmp->perim, tmp->perim/tmp->area, 0.282*tmp->perim/sqrt(tmp->area), 2.0*sqrt(tmp->area/PI)/tmp->long_axis, tmp->twist, tmp->omega); fclose(outfile); } tmp = tmp->next; } if (choice->att[1] || choice->att[2] || choice->att[3] || choice->att[4]) fclose(a1_4); if (choice->att[5]) fclose(a5); if (choice->att[6]) fclose(a6); if (choice->att[7]) fclose(a7); if (choice->att[8]) fclose(a8); if (choice->core[1] || choice->core[2] || choice->core[3] || choice->core[4]) fclose(c1_4); if (choice->core[5]) fclose(c5); if (choice->core[6]) fclose(c6); if (choice->core[7]) fclose(c7); if (choice->core[8]) fclose(c8); if (choice->core[9]) { fclose(c9c); fclose(c9e); } if (choice->core[10]) { fclose(c10c); fclose(c10e); } if (choice->size[1] || choice->size[2]) fclose(s1_2); if (choice->size[3]) fclose(s3); if (choice->size[4]) fclose(s4); if (choice->size[5]) fclose(s5); if (choice->size[6]) fclose(s6); if (choice->size[7] || choice->size[8]) fclose(s7_8); if (choice->shape[1] || choice->shape[2]) fclose(h1_2); if (choice->shape[3]) fclose(h3); if (choice->shape[4]) fclose(h4); if (choice->shape[5]) fclose(h5); if (choice->shape[6]) fclose(h6); if (choice->boundary[1] || choice->boundary[2] || choice->boundary[3] || choice->boundary[4]) fclose(n1_4); if (choice->perim[1] || choice->perim[2] || choice->perim[3]) fclose(p1_3); if (choice->perim[4]) fclose(p4); if (choice->perim[5]) fclose(p5); if (choice->perim[6]) fclose(p6); G_free (type_dens); total_patches = 0; return; } /************************************/ /* COMPUTE THE ATTRIBUTE MEASURES */ /************************************/ void df_att (PATCH *tmp, int type_coh, int *type_dens) { static double sumx=0.0, sumx2=0.0, w_att=0.0, w_att2=0.0, total=0.0; static double *area, total2=0.0; /* variables: IN: tmp = the next patch on the patch list type_coh = identification no. for the group for this patch *type_dens = array of no. of patches by group INTERNAL: sumx = sum of all the patch attributes sumx2 = sum of all the patch attributes squared w_att = sum of (patch attributes x areas) w_att2 = sum of (patch attributes x areas squared) total = sum of all the patch areas *area = array of total areas by gp */ if (tmp->num == 1) area = (double *)G_calloc(25, sizeof(double)); sumx += tmp->att; sumx2 += tmp->att * tmp->att; w_att += tmp->area * tmp->att; w_att2 += tmp->area * tmp->att * tmp->att; total += tmp->area; total2 += tmp->area * tmp->area; area[type_coh] += tmp->area; if (!tmp->next) { save_att(w_att, w_att2, total, total2, sumx, sumx2, type_dens, area); w_att = w_att2 = total = sumx = sumx2 = 0.0; G_free (area); } return; } /************************************/ /* SAVE THE ATTRIBUTE MEASURES */ /************************************/ void save_att (double w_att, double w_att2, double t_size, double t_size2, double sum, double sum2, int *density, double *area) { register int i; double wm, wstdv, m, stdv; /* variables: IN: from df_att w_att = sum of (patch attributes x areas) w_att2 = sum of (patch attributes x areas squared) t_size = sum of all the patch areas t_size2 = sum of all the patch areas squared sum = sum of all the patch attributes sum2 = sum of all the patch attributes squared *area = array of total areas by group *density = array of no. of patches by group INTERNAL: wm = mean pixel attribute (a1) wstdv = st. dev. pixel attribute (a2) m = mean patch attribute (a3) stdv = st. dev. patch attribute (a4) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ wm = w_att/t_size; wstdv = w_att2/t_size - wm*wm; if (wstdv > 0) wstdv = sqrt(wstdv); else wstdv = 0; m = sum/total_patches; stdv = sum2/total_patches - m*m; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0; /* write a1=mn. pix. att., a2=s.d. pix. att., a3=mn. patch att., a4=s.d. patch att. */ if (choice->att[1] || choice->att[2] || choice->att[3] || choice->att[4]) { fprintf(a1_4, " %15.3f %15.3f %15.3f %15.3f\n", wm, wstdv, m, stdv); } /* write a5 = cover by gp */ if (choice->att[5]) { for(i = 0; i < ntype; i++) fprintf(a5, " %11.3f", area[i]/t_size); fprintf(a5, "\n"); } /* write a6 = density by gp */ if (choice->att[6]) { for(i = 0; i < ntype; i++) fprintf(a6, " %10d", *(density + i)); fprintf(a6, "\n"); } /* write a7 = total patches */ if (choice->att[7]) fprintf(a7, " %11d\n", total_patches); /* write a8 = eff. mesh no. */ if (choice->att[8]) { if (t_size2 > 0.0) fprintf(a8, " %11.3f\n", (t_size*t_size)/t_size2); else fprintf(a8, " %11.3f\n", t_size2); } return; } /************************************/ /* COMPUTE THE CORE MEASURES */ /************************************/ void df_core (PATCH *tmp, int type_coh, int *type_dens) { static int **density1c=NULL, **density1e=NULL, *densityc=NULL, *densitye=NULL, first=1; static double mcore=0.0, medge=0.0, *mcore1=NULL, *medge1=NULL, sumc2=0.0, sume2=0.0, *sum22c, *sum22e; register int i; int core_coh, edge_coh; /* variables: IN: tmp = the next patch on the patch list type_coh = identif. no. for the group for this patch type_dens[] = array of no. of patches by group INTERNAL: mcore = sum of patch core sizes medge = sum of patch edge sizes sumc2 = sum of patch cores squared sume2 = sum of patch edges squared mcore1[] = array of patch cores by group medge1[] = array of patch edges by group sum22c[] = array of patch cores squared by group sum22e[] = array of patch edges squared by group densityc[] = array of no. of patches by core size class densitye[] = array of no. of patches by edge size class density1c[] = array of no. of patches by core size class by group density1e[] = array of no. of patches by edge size class by group */ if (first){ densityc = (int *)G_calloc(25, sizeof(int)); densitye = (int *)G_calloc(25, sizeof(int)); sum22c = (double *)G_calloc(25, sizeof(double)); sum22e = (double *)G_calloc(25, sizeof(double)); mcore1 = (double *)G_calloc(25, sizeof(double)); medge1 = (double *)G_calloc(25, sizeof(double)); } /* if output is by size class (c9 & c10) determine which size class the current patch is in */ if (choice->core[9] || choice->core[10]) { core_coh = index_coh(tmp->core, size_cl); densityc[core_coh] ++; edge_coh = index_coh(tmp->edge, size_cl); densitye[edge_coh] ++; } mcore += tmp->core; medge += tmp->edge; sumc2 += tmp->core * tmp->core; sume2 += tmp->edge * tmp->edge; if (type_coh >= 0) { mcore1[type_coh] += tmp->core; medge1[type_coh] += tmp->edge; sum22c[type_coh] += tmp->core * tmp->core; sum22e[type_coh] += tmp->edge * tmp->edge; } /* if c10 */ if (choice->core2) { if (first) { density1c = (int **)G_calloc(25, sizeof(int *)); for(i = 0; i < 25; i++) density1c[i] = (int *)G_calloc(25, sizeof(int)); density1e = (int **)G_calloc(25, sizeof(int *)); for(i = 0; i < 25; i++) density1e[i] = (int *)G_calloc(25, sizeof(int)); } if (type_coh >= 0) { if (core_coh >= 0) density1c[type_coh][core_coh] ++; if (edge_coh >= 0) density1e[type_coh][edge_coh] ++; } } if (first) first = 0; if (!tmp->next) { save_core(sumc2, sume2, mcore, medge, mcore1, medge1, sum22c, sum22e, densityc, densitye, type_dens, density1c, density1e); mcore = medge = sumc2 = sume2 = 0; G_free (densityc); G_free (densitye); G_free (sum22c); G_free (sum22e); G_free (mcore1); G_free (medge1); if (choice->core2) { for(i = 0; i < 25; i++) { G_free (density1c[i]); G_free (density1e[i]); } G_free (density1c); G_free (density1e); } first = 1; } return; } /************************************/ /* SAVE THE CORE MEASURES */ /************************************/ void save_core (double sumc2, double sume2, double mcore, double medge, double *mcore1, double *medge1, double *sum22c, double *sum22e, int *densityc, int *densitye, int *type_dens, int **density1c, int **density1e) { register int i, j; double tmpc, tmpe, stdvc, stdve; /* variables: IN: sumc2 = sum of patch cores squared sume2 = sum of patch edges squared mcore = sum of patch cores medge = sum of patch edges mcore1[] = array of sums of patch cores by gp medge1[] = array of sums of patch edges by gp sum22c[] = array of patch cores squared by gp sum22e[] = array of patch edges squared by gp densityc[] = array of no. of patches by core size class densitye[] = array of no. of patches by edge size class type_dens[] = array of no. of patches by gp ? INTERNAL: tmpc = mean patch core (c1 & c5) stdvc = st. dev. patch core (c2 & c6) tmpe = mean patch edge (c3 & c7) stdve = st. dev. patch edge (c4 & c8) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ /* calc. & write c1=mn. core, c2=s.d. core, c3=mn. edge, c4=s.d. edge */ if (choice->core[1] || choice->core[2] || choice->core[3] || choice->core[4]) { tmpc = mcore/total_patches; stdvc = sumc2/total_patches - tmpc*tmpc; if (stdvc > 0) stdvc = sqrt(stdvc); else stdvc = 0; tmpe = medge/total_patches; stdve = sume2/total_patches - tmpe*tmpe; if (stdve > 0) stdve = sqrt(stdve); else stdve = 0; fprintf(c1_4, " %15.3f %15.3f %15.3f %15.3f\n", tmpc, stdvc, tmpe, stdve); } /* calc. & write c5=mn. core by gp */ if (choice->core[5]) { for(i = 0; i < ntype; i++){ if (tmpc = type_dens[i]) tmpc = mcore1[i] / tmpc; fprintf(c5, " %11.3f", tmpc); } fprintf(c5, "\n"); } /* calc. & write c6=s.d. core by gp */ if (choice->core[6]) { for(i = 0; i < ntype; i++){ stdvc = 0; if (type_dens[i]){ tmpc = mcore1[i] / type_dens[i]; stdvc = sum22c[i]/type_dens[i] - tmpc*tmpc; if (stdvc > 0) stdvc = sqrt(stdvc); else stdvc = 0; } fprintf(c6, " %11.3f", stdvc); } fprintf(c6, "\n"); } /* calc. & write c7=mn. edge by gp */ if (choice->core[7]) { for(i = 0; i < ntype; i++){ if (tmpe = type_dens[i]) tmpe = medge1[i] / tmpe; fprintf(c7, " %11.3f", tmpe); } fprintf(c7, "\n"); } /* calc. & write c8=s.d. edge by gp */ if (choice->core[8]) { for(i = 0; i < ntype; i++){ stdve = 0; if (type_dens[i]){ tmpe = medge1[i] / type_dens[i]; stdve = sum22e[i]/type_dens[i] - tmpe*tmpe; if(stdve > 0) stdve = sqrt(stdve); else stdve = 0; } fprintf(c8, " %11.3f", stdve); } fprintf(c8, "\n"); } /* write c9=no. by size class */ if (choice->core[9]) { for(i = 0; i < size_cl[0] - 1; i++) fprintf(c9c, " %11d", *(densityc+i)); fprintf(c9c, "\n"); for(i = 0; i < size_cl[0] - 1; i++) fprintf(c9e, " %11d", *(densitye+i)); fprintf(c9e, "\n"); } /* write c10=no. by size class by gp */ if (choice->core2){ for(i = 0; i < ntype; i++){ fprintf(c10c, " Gp[%2d]",i + 1); for(j = 0; j < *size_cl - 1; j++) fprintf(c10c, " %11d", density1c[i][j]); fprintf(c10c, "\n"); } for(i = 0; i < ntype; i++){ fprintf(c10e, " Gp[%2d]",i + 1); for(j = 0; j < *size_cl - 1; j++) fprintf(c10e, " %11d", density1e[i][j]); fprintf(c10e, "\n"); } } return; } /************************************/ /* COMPUTE THE SIZE MEASURES */ /************************************/ void df_size (PATCH *tmp, int type_coh, int *type_dens) { static int **density1=NULL, *density=NULL, first=1; static double msize=0.0, *msize1=NULL, sum2=0.0, *sum22; register int i; int size_coh; /* variables: IN: tmp = the next patch on the patch list type_coh = identif. no. for the gp for this patch type_dens[] = array of no. of patches by gp INTERNAL: msize = sum of patch sizes msize1[] = array of patch sizes by gp sum2 = sum of patch sizes squared sum22[] = array of patch sizes squared by gp density[] = array of no. of patches by size class */ if (first){ density = (int *)G_calloc(25, sizeof(int)); sum22 = (double *)G_calloc(25, sizeof(double)); msize1 = (double *)G_calloc(25, sizeof(double)); } /* if output is by size class (s5 & s6) determine which size class the current patch is in */ if (choice->size[5] || choice->size[6]) { size_coh = index_coh(tmp->area, size_cl); density[size_coh] ++; } msize += tmp->area; sum2 += tmp->area * tmp->area; if (type_coh >= 0) { msize1[type_coh] += tmp->area; sum22[type_coh] += tmp->area * tmp->area; } if (choice->size2) { if(first) { density1 = (int **)G_calloc(25, sizeof(int *)); for(i = 0; i < 25; i++) density1[i] = (int *)G_calloc(25, sizeof(int)); } if (type_coh >= 0 && size_coh >= 0) density1[type_coh][size_coh] ++; } if (first) first = 0; if (!tmp->next) { save_size(sum2, msize, msize1, sum22, density, type_dens, density1); msize = sum2 = 0.0; G_free (density); G_free (msize1); G_free (sum22); if (density1) { for(i = 0; i < 25; i++) G_free (density1[i]); G_free (density1); } first = 1; } return; } /************************************/ /* SAVE THE SIZE MEASURES */ /************************************/ void save_size (double sum2, double msize, double *msize1, double *sum22, int *density, int *type_dens, int **density1) { register int i, j; double tmp, stdv; /* variables: IN: sum2 = sum of patch sizes squared msize = sum of patch sizes msize1[] = array of sums of patch sizes by gp sum22[] = array of patch sizes squared by gp density[] = array of no. of patches by size class type_dens[] = array of no. of patches by gp INTERNAL: tmp = mean patch size (s1 & s3) stdv = st. dev. patch size (s2 & s4) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ /* calc. & write s1=mn. size, s2=s.d. size */ if (choice->size[1] || choice->size[2]) { tmp = msize/total_patches; stdv = sum2/total_patches - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; fprintf(s1_2, " %15.3f %15.3f\n", tmp, stdv); } /* calc. & write s3=mn. size by gp */ if (choice->size[3]) { for(i = 0; i < ntype; i++) { if (tmp = type_dens[i]) tmp = msize1[i] / tmp; fprintf(s3, " %11.3f", tmp); } fprintf(s3, "\n"); } /* calc. & write s4=s.d. size by gp */ if (choice->size[4]) { for(i = 0; i < ntype; i++) { stdv = 0.0; if (type_dens[i]){ tmp = msize1[i] / type_dens[i]; stdv = sum22[i]/type_dens[i] - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; } fprintf(s4, " %11.3f", stdv); } fprintf(s4, "\n"); } /* write s5=no. by size class */ if (choice->size[5]) { for(i = 0; i < size_cl[0] - 1; i++) fprintf(s5, " %11d", *(density+i)); fprintf(s5, "\n"); } /* write s6=no. by size class by gp */ if (choice->size2) { /* if(!(fp = fopen("r.le.out/s6.out", "a"))) G_fatal_error("Can't write file s6.out; do you have write permission?\n"); fprintf(fp, "%5d %5d\n", n_scale, n_unit); */ for(i = 0; i < ntype; i++) { fprintf(s6, " Gp[%2d]",i + 1); for(j = 0; j < *size_cl - 1; j++) fprintf(s6, " %11d", density1[i][j]); fprintf(s6, "\n"); } } /* calculate & write s7 (eff. mesh size) & s8 (deg. landsc. division) */ if (choice->size[7] || choice->size[8]) fprintf(s7_8, " %15.3f %15.3f\n", (1.0 / msize) * sum2, 1.0 - (sum2 / (msize * msize))); return; } /************************************/ /* COMPUTE THE SHAPE MEASURES */ /************************************/ void df_shape (PATCH *tmp, int type_coh, int *type_dens) { static int *den1, *den2, *den3, **density1=NULL, **density2=NULL, **density3, new=1; static double mshape=0.0, *mshape1=NULL, mshape_p=0.0, *mshape1_p, *sqr11, *sqr21, *sqr31, mshape_r=0.0, *mshape1_r, sq1=0.0, sq2=0.0, sq3=0.0; register int i; int shape_coh1=0, shape_coh2=0, shape_coh3=0; double shp1, shp2, shp3; /* variables: IN: tmp = the next patch on the patch list type_coh = identif. no. for the gp for this patch type_dens[] = array of no. of patches by gp INTERNAL: shp1 = CPA shape index shp2 = PA shape index shp3 = RCC shape index mshape = sum of CPA shape indices for all patches mshape_p = sum of PA shape indices for all patches mshape_r = sum of RCC shape indices for all patches mshape1[] = array of CPA indices by gp mshape1_p[] = array of PA indices by gp mshape1_r[] = array of RCC indices by gp sq1 = sum of CPA indices squared sq2 = sum of PA indices squared sq3 = sum of RCC indices squared sqr11[] = array of CPA indices squared by gp sqr21[] = array of PA indices squared by gp sqr31[] = array of RCC indices squared by gp den1[] = array of CPA indices by index class den2[] = array of PA indices by index class den3[] = array of RCC indices by index class density1[] = array of CPA indices by index class by gp density2[] = array of PA indices by index class by gp density3[] = array of RCC indices by index class by gp */ shp1 = 0.282 * tmp->perim / sqrt(tmp->area); /* CPA method m2 */ shp2 = tmp->perim / tmp->area; /* PA method m1 */ shp3 = 2.0 * sqrt(tmp->area / PI) / tmp->long_axis; /* RCC method m3 */ if (new) { mshape1 = (double *)G_calloc(25, sizeof(double)); mshape1_p = (double *)G_calloc(25, sizeof(double)); mshape1_r = (double *)G_calloc(25, sizeof(double)); sqr11 = (double *)G_calloc(25, sizeof(double)); sqr21 = (double *)G_calloc(25, sizeof(double)); sqr31 = (double *)G_calloc(25, sizeof(double)); den1 = (int *)G_calloc(25, sizeof(int)); den2 = (int *)G_calloc(25, sizeof(int)); den3 = (int *)G_calloc(25, sizeof(int)); } mshape += shp1; mshape_p += shp2; mshape_r += shp3; sq1 += shp1 * shp1; sq2 += shp2 * shp2; sq3 += shp3 * shp3; if (type_coh >= 0) { mshape1[type_coh] += shp1; mshape1_p[type_coh] += shp2; mshape1_r[type_coh] += shp3; sqr11[type_coh] += shp1*shp1; sqr21[type_coh] += shp2*shp2; sqr31[type_coh] += shp3*shp3; } /* if sh2=h5 or h6 */ if (choice->shape[5] || choice->shape[6]) { /* if sh1=m1 */ if (choice->Mx[1]) { if (0 <= (shape_coh2 = index_coh(shp2, shape_PA))) den2[shape_coh2] ++; } /* if sh1=m2 */ if (choice->Mx[2]) { if (0 <= (shape_coh1 = index_coh(shp1, shape_CPA))) den1[shape_coh1] ++; } /* if sh1=m3 */ if (choice->Mx[3]) { if (0 <= (shape_coh3 = index_coh(shp3, shape_RCC))) den3[shape_coh3] ++; } } /* if sh2=h6 */ if (choice->shape2) { if (new) { density1 = (int **)G_calloc(25, sizeof(int *)); density2 = (int **)G_calloc(25, sizeof(int *)); density3 = (int **)G_calloc(25, sizeof(int *)); for(i = 0; i < 25; i++) if (!(density1[i] = (int *)G_calloc (25, sizeof(int))) || !(density2[i] = (int *)G_calloc (25, sizeof(int))) || !(density3[i] = (int *)G_calloc (25, sizeof(int))) ) G_fatal_error("Failure to allocate memory for sh2, exit.\n"); } if (type_coh >= 0) { if (shape_coh1 >= 0) density1[type_coh][shape_coh1] ++; if (shape_coh2 >= 0) density2[type_coh][shape_coh2] ++; if (shape_coh3 >= 0) density3[type_coh][shape_coh3] ++; } } if (new) new = 0; if (!tmp->next) { /** if all the patches are done **/ save_shape(sq1, sq2, sq3, sqr11, sqr21, sqr31, mshape, mshape_p, mshape_r, mshape1, mshape1_p, mshape1_r, type_dens, den1, den2, den3, density1, density2, density3); mshape = sq1 = sq2 =sq3 = mshape_p = mshape_r = 0; G_free (mshape1); G_free (mshape1_p); G_free (mshape1_r); G_free (sqr11); G_free (sqr21); G_free (sqr31); G_free (den1); G_free (den2); G_free (den3); if (density1) { for(i = 0; i < 25; i++) { G_free (density1[i]); G_free (density2[i]); G_free (density3[i]); } G_free (density1); G_free (density2); G_free (density3); } new = 1; } return; } /************************************/ /* SAVE THE SHAPE MEASURES */ /************************************/ void save_shape (double sq1, double sq2, double sq3, double *sqr11, double *sqr21, double *sqr31, double mshape, double mshape_p, double mshape_r, double *mshape1, double *mshape1_p, double *mshape1_r, int *type_dens, int *den1, int *den2, int *den3, int **density1, int **density2, int **density3) { register int i, j; double tmp, stdv; /* variables: IN: sq1 = sum of CPA indices squared sq2 = sum of PA indices squared sq3 = sum of RCC indices squared sqr11[] = array of CPA indices squared by gp sqr21[] = array of PA indices squared by gp sqr31[] = array of RCC indices squared by gp mshape = sum of CPA shape indices for all patches mshape_p = sum of PA shape indices for all patches mshape_r = sum of RCC shape indices for all patches mshape1[] = array of CPA indices by gp mshape1_p[] = array of PA indices by gp mshape1_r[] = array of RCC indices by gp type_dens[] = array of no. of patches by gp ? den1[] = array of CPA indices by index class den2[] = array of PA indices by index class den3[] = array of RCC indices by index class density1[] = array of CPA indices by index class by gp density2[] = array of PA indices by index class by gp density3[] = array of RCC indices by index class by gp INTERNAL: tmp = mean shape index (h1) stdv = st. dev. shape index (h2) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ /*CALC. & WRITE h1-h5 FOR CPA INDEX (m2) */ /* calc. & write h1=mn. shape, h2=s.d. shape for CPA index (m2) */ if ((choice->shape[1] || choice->shape[1]) && choice->Mx[2]) { tmp = mshape / total_patches; stdv = sq1/total_patches - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; fprintf(h1_2, " %15.3f %15.3f\n", tmp, stdv); } /* calc. & write h3=mn. shape by gp for CPA index (m2) */ if (choice->shape[3] && choice->Mx[2]) { for(i = 0; i < ntype; i++) { if (tmp = type_dens[i]) tmp = mshape1[i] / tmp; fprintf(h3, " %10.3f", tmp); } fprintf(h3, "\n"); } /* calc. & write h4=s.d. shape by gp for CPA index (m2) */ if (choice->shape[4] && choice->Mx[2]) { for(i = 0; i < ntype; i++) { stdv = 0.0; if (type_dens[i] > 1) { tmp = mshape1[i] / type_dens[i]; stdv = sqr11[i]/type_dens[i] - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; } fprintf(h4, " %10.3f", stdv); } fprintf(h4, "\n"); } /* write h5=no. by shape index class for CPA index (m2) */ if (choice->shape[5] && choice->Mx[2]) { for(j = 0; j < *shape_CPA - 1; j++) fprintf(h5, " %10d", den1[j]); fprintf(h5, "\n"); } /*CALC. & WRITE h1-h5 FOR PA INDEX (m1) */ /* calc. & write h1=mn. shape, h2=s.d. shape for PA index (m1) */ if ((choice->shape[1] || choice->shape[2]) && choice->Mx[1]) { tmp = mshape_p / total_patches; stdv = sq2/total_patches - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; fprintf(h1_2, " %15.3f %15.3f\n", tmp, stdv); } /* calc. & write h3=mn. shape by gp for PA index (m1) */ if (choice->shape[3] && choice->Mx[1]) { for(i = 0; i < ntype; i++) { if (tmp = type_dens[i]) tmp = mshape1_p[i] / tmp; fprintf(h3, " %10.3f", tmp); } fprintf(h3, "\n"); } /* calc. & write h4=s.d. shape by gp for PA index (m1) */ if (choice->shape[4] && choice->Mx[1]) { for(i = 0; i < ntype; i++){ stdv = 0.0; if (type_dens[i] > 1){ tmp = mshape1_p[i] / type_dens[i]; stdv = sqr21[i]/type_dens[i] - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0; } fprintf(h4, " %10.3f", stdv); } fprintf(h4, "\n"); } /* write h5=no. by shape index class for PA index (m1) */ if (choice->shape[5] && choice->Mx[1]) { for(j = 0; j < *shape_PA - 1; j++) fprintf(h5, " %10d", den2[j]); fprintf(h5, "\n"); } /*CALC. & WRITE h1-h5 FOR RCC INDEX (m3) */ /* calc. & write h1=mn. shape, h2=s.d. shape for RCC index (m3) */ if ((choice->shape[1] || choice->shape[2]) && choice->Mx[3]) { tmp = mshape_r / total_patches; stdv = sq3/total_patches - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; fprintf(h1_2, " %15.3f %15.3f\n", tmp, stdv); } /* calc. & write h3=mn. shape by gp for RCC index (m3) */ if (choice->shape[3] && choice->Mx[3]) { for(i = 0; i < ntype; i++) { if (tmp = type_dens[i]) tmp = mshape1_r[i] / tmp; fprintf(h3, " %10.3f", tmp); } fprintf(h3, "\n"); } /* calc. & write h4=s.d. shape by gp for RCC index (m3) */ if (choice->shape[4] && choice->Mx[3]) { for(i=0; i 1) { tmp = mshape1_r[i] / type_dens[i]; stdv = sqr31[i]/type_dens[i] - tmp*tmp; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; } fprintf(h4, " %10.3f", stdv); } fprintf(h4, "\n"); } /* write h5=no. by shape index class for RCC index (m3) */ if (choice->shape[5] && choice->Mx[3]) { for(j = 0; j < *shape_RCC - 1; j++) fprintf(h5, " %10d", den3[j]); fprintf(h5, "\n"); } /* CALC. & WRITE h6 = NO. IN EA. SHAPE INDEX CLASS BY GP */ if (choice->shape[6]) { if (density1) { if (choice->Mx[1]) { for(i = 0; i < ntype; i++) { fprintf(h6, " Gp[%2d]",i + 1); for(j = 0; j < *shape_PA - 1; j++) fprintf(h6, " %10d", density2[i][j]); fprintf(h6, "\n"); } } if (choice->Mx[2]) { for(i = 0; i < ntype; i++) { fprintf(h6, " Gp[%2d]",i + 1); for(j = 0; j < *shape_CPA - 1; j++) fprintf(h6, " %10d", density1[i][j]); fprintf(h6, "\n"); } } if (choice->Mx[3]) { for(i = 0; i < ntype; i++) { fprintf(h6, " Gp[%2d]",i + 1); for(j = 0; j < *shape_RCC - 1; j++) fprintf(h6, " %10d", density3[i][j]); fprintf(h6, "\n"); } } } } return; } /****************************************************/ /* COMPUTE AND SAVE BOUNDARY COMPLEXITY MEASURES */ /****************************************************/ void df_boundary (PATCH *tmp) { static double sumomega=0.0, sumomega2=0.0; static int sumtwist=0, sumtwist2=0; register int i; double meantwist=0.0, stdvtwist=0.0, meanomega=0.0, stdvomega=0.0; /* variables: IN: tmp = the next patch in the patch list INTERNAL: meantwist = mean twist number (n1) stdvtwist = st. dev. twist number (n2) meanomega = mean omega index (n3) stdvomega = st. dev. omega index (n4) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ sumtwist += tmp->twist; sumtwist2 += tmp->twist * tmp->twist; sumomega += tmp->omega; sumomega2 += tmp->omega * tmp->omega; if (!tmp->next) { meantwist = (double)sumtwist/total_patches; stdvtwist = (double)sumtwist2/total_patches - meantwist*meantwist; if (stdvtwist > 0.0) stdvtwist = sqrt(stdvtwist); else stdvtwist = 0.0; meanomega = (double)sumomega/total_patches; stdvomega = (double)sumomega2/total_patches - meanomega*meanomega; if (stdvomega > 0.0) stdvomega = sqrt(stdvomega); else stdvomega = 0.0; /* write n1=mn. twist, n2=s.d. twist, n3=mn. omega, n4=s.d. omega */ if (choice->boundary[1] || choice->boundary[2] || choice->boundary[3] || choice->boundary[4]) fprintf(n1_4, " %15.3f %15.3f %15.3f %15.3f\n", meantwist, stdvtwist, meanomega, stdvomega); meantwist = stdvtwist = meanomega = stdvomega = 0.0; } return; } /*************************************/ /* COMPUTE & SAVE PERIMETER MEASURES */ /*************************************/ void df_perim (PATCH *tmp, int type_coh, int *type_dens) { static double perim=0.0, *perim1, sum2=0.0, *sum21, first=1.0; register int i; double mean, stdv; if (first){ perim1 = (double *)G_calloc(25, sizeof(double)); sum21 = (double *)G_calloc(25, sizeof(double)); first = 0.0; } /* variables: IN: tmp = the next patch in the patch list type_coh = identif. no. of the gp for this patch type_dens[] = array of no. of patches by gp INTERNAL: perim = sum of perimeters (p1) sum2 = sum of perimeters squared perim1[] = array of sum of perims by gp (pp4) sum21[] = array of sum of perims squared by gp mean = mean perim. (p2) stdv = st. dev. perim. (p3) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ perim += tmp->perim; sum2 += tmp->perim*tmp->perim; if (type_coh >= 0) { perim1[type_coh] += tmp->perim; sum21[type_coh] += tmp->perim*tmp->perim; } if (!tmp->next) { /** save the perimeter measures **/ mean = perim/total_patches; stdv = sum2/total_patches - mean*mean; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; /* write p1=sum per., p2=mn. per.,p3=s.d. per. */ if (choice->perim[1] || choice->perim[2] || choice->perim[3]) fprintf(p1_3, " %15.3f %15.3f %15.3f\n", perim, mean, stdv); /* write p4=sum per. by gp */ if (choice->perim[4]) { for(i = 0; i < ntype; i++) fprintf(p4, " %11.3f", perim1[i]); fputs("\n", p4); } /* write p5=mn. per. by gp */ if (choice->perim[5]) { for(i = 0; i < ntype; i++) { if (type_dens[i]) mean = perim1[i]/type_dens[i]; else mean = 0.0; fprintf(p5, " %11.3f", mean); } fputs("\n", p5); } /* calc. & write p6=s.d. per. by gp */ if (choice->perim[6]) { for(i = 0; i < ntype; i++) { stdv = 0; if (type_dens[i]) { mean = perim1[i]/type_dens[i]; stdv = sum21[i]/type_dens[i] - mean*mean; if (stdv > 0) stdv = sqrt(stdv); else stdv = 0.0; } fprintf(p6, " %11.3f", stdv); } fputs("\n", p6); } G_free (perim1); G_free (sum21); first = 1; perim = sum2 = 0; } return; } /************************************/ /* RUN THE MOVING WINDOW PATCH */ /* MEASURES */ /************************************/ void mv_patch (PATCH *patch_list, double **value, int index) { PATCH *tmp = patch_list; /* Variables: IN: patch_list = the list of patches value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: *tmp = pointer to the list of patches */ if (!total_patches) return; while(tmp) { if (choice->att[0]) m_att(tmp, value, index); if (choice->size[0]) m_size(tmp, value, index); if (choice->core[0]) m_core(tmp, value, index); if (choice->shape[0] && choice->Mx[1]) m_shape(tmp, 1, value, index); if (choice->shape[0] && choice->Mx[2]) m_shape(tmp, 2, value, index); if (choice->shape[0] && choice->Mx[3]) m_shape(tmp, 3, value, index); if (choice->boundary) m_boundary(tmp, value, index); if (choice->perim[0]) m_perim(tmp, value, index); tmp = tmp->next; } total_patches = 0; return; } /************************************/ /* FOR DF_PATCH PROGRAM (NOT MOVING */ /* WINDOW) DETERMINE WHICH GROUP */ /* ATT BELONGS TO */ /************************************/ int recl_coh (double att) { register int i; extern int ntype; extern float **recl_tb; for(i = 0; i < ntype; i++) { if (in_group(att, recl_tb[i], i)) { return i; } } return -999; } /************************************/ /* DETERMINE WHETHER ATT BELONGS TO */ /* THE CHOSEN rTH GP IN THE RECLASS */ /* TABLE */ /************************************/ int in_group (double att, float *group, int r) { register int i; /*Variables EXTERNAL recl_count = an array containg the number of elements in each row of the reclass table IN att = the attribute value group = an array containing the rth row of the reclass table r = the chosen row of the reclass table INTERNAL i = index */ /* For each element i in the rth row of the reclass table. When this program is called by the moving window programs r is always 0 */ for(i = 1; i < recl_count[r]; i++) { /* if the i element of the rth row of the reclass table is -999, this indicates "thru" in the reclass table, so check to see if the attribute lies in the range expressed in the "thru" statement */ if (-999 == *(group + i)) { if (*(group + i) && *(group + i - 1) <= att && att <= *(group + i + 1)) return 1; else i++; } /* if the rth row of the reclass table does not contain a "thru" statement, then just check whether the i element of the rth row of the reclass table is equal to the attribute */ else if ((double)*(group + i) == att) { return 1; } } return 0; } /************************************/ /* DETERMINE WHICH INDEX CLASS ATT */ /* BELONGS TO */ /************************************/ int index_coh (double att, float *group) { register int i; /* Variables: IN att = the attribute to be checked group = an array of index class limits (e.g., size classes) */ for(i = (int)*group - 1; i >= 1; i--) { if ((double)*(group + i) <= att) return i-1; } return -999; } /************************************/ /* MOVING WINDOW ATTRIBUTE MEASURES */ /************************************/ void m_att (PATCH *tmp, double **value, int index) { static double sum1=0.0, sum12=0.0, sum2=0.0, sum22=0.0, sum32=0.0, total1=0.0, total2=0.0, area=0.0, area2=0.0; static int density=0; double mean, stdv; /* choice->att 1 = mean pixel attrib. (a1) 2 = st. dev. pixel attrib. (a2) 3 = mean patch attrib. (a3) 4 = st. dev. patch attrib. (a4) 5 = cover in gp 1 (a5) 6 = density in gp 1 (a6) 7 = total density (a7) 8 = eff. mesh no. (a8) Variables: IN: tmp = the list of patches value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: sum1 = sum of patch area x patch attrib. sum2 = sum of patch attributes sum12 = sum of patch area x patch attrib. squared sum22 = sum of patch attributes squared sum32 = sum of patch areas in gp 1 total1 = sum of patch areas for a1 and a2 total2 = sum of patch areas for a5 area = sum of patch areas for a8 area2 = sum of patch areas squared density = no. of patches in gp 1 value = output value for selected att measure GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ if (choice->att[1] || choice->att[2]) { sum1 += tmp->area * tmp->att; total1 += tmp->area; if (choice->att[2]) sum12 += tmp->area * tmp->att * tmp->att; } if (choice->att[3] || choice->att[4]) { sum2 += tmp->att; if (choice->att[4]) sum22 += tmp->att * tmp->att; } if (choice->att[5]) { total2 += tmp->area; if (in_group(tmp->att, recl_tb[0], 0)) { sum32 += tmp->area; } } if (choice->att[6]) { if (in_group(tmp->att, recl_tb[0], 0)) density++; } if (choice->att[8]) { area += tmp->area; area2 += tmp->area * tmp->area; } if (!tmp->next) { /* calc. a1 = mn. pixel attrib. */ if (choice->att[1] && total1) { value[index][0] = sum1/total1; } /* calc. a2 = s.d. pixel attrib. */ if (choice->att[2] && total1) { mean = sum1/total1; stdv = sum12/total1 - mean*mean; if (stdv > 0) value[index][1] = sqrt(stdv); } /* calc. a3 = mn. patch attrib. */ if (choice->att[3] && total_patches) { value[index][2] = sum2/total_patches; } /* calc. a4 = s.d. patch attrib. */ if (choice->att[4] && total_patches) { mean = sum2/total_patches; stdv = sum22/total_patches - mean*mean; if (stdv > 0) value[index][3] = sqrt(stdv); } /* calc. a5 = cover in gp 1 */ if (choice->att[5] && total2) { value[index][4] = sum32/total2; } /* calc. a6 = density in gp 1 */ if (choice->att[6]) { value[index][5] = density; } /* calc. a7 = total density */ if (choice->att[7]) value[index][6] = total_patches; /* calc. a8 = eff. mesh no. */ if (choice->att[8]) { value[index][36] = (area * area)/area2; } total1 = total2 = sum1 = sum12 = 0.0; sum2 = sum22 = sum32 = area = area2 = 0.0; density = 0; } return; } /************************************/ /* MOVING WINDOW SIZE MEASURES */ /************************************/ void m_size (PATCH *tmp, double **value, int index) { static double sum1=0.0, sum12=0.0, sum2=0.0, sum22=0.0; static int density1=0, density2=0, density3=0; double mean, stdv; /* choice->size == 1 - mean patch size (s1) 2 - st. dev. patch size (s2) 3 - mean patch size by gp 1 (s3) 4 - st. dev. patch size by gp 1 (s4) 5 - no. by size class 1 (s5) 6 - no. by size class 1 by gp 1 (s6) Variables: IN: tmp = the list of patches value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: sum1 = sum of patch sizes sum12 = sum of patch sizes squared sum2 = sum of patch sizes in gp 1 sum22 = sum of patch sizes in gp 1 squared density1 = no. of patches in gp 1 density2 = no. of patches in size class 1 density3 = no. of patches in gp 1 and size class 1 */ if (choice->size[1] || choice->size[2] || choice->size[7] || choice->size[8]) { sum1 += tmp->area; if (choice->size[2] || choice->size[7] || choice->size[8]) sum12 += tmp->area * tmp->area; } if (choice->size[3] || choice->size[4]) { if (in_group(tmp->att, recl_tb[0], 0)) { density1 ++; sum2 += tmp->area; if (choice->size[4]) sum22 += tmp->area * tmp->area; } } if (choice->size[5] && tmp->area < size_cl[2]) { density2 ++; } if (choice->size[6] && tmp->area < size_cl[2] && in_group(tmp->att, recl_tb[0], 0)) { density3 ++; } if (!tmp->next){ /* calc. s1 = mn. patch size */ if (choice->size[1] && total_patches) { value[index][7] = sum1/total_patches; } /* calc. s2 = s.d. patch size */ if (choice->size[2] && total_patches) { mean = sum1/total_patches; stdv = sum12/total_patches - mean*mean; if (stdv > 0) value[index][8] = sqrt(stdv); } /* calc. s3 = mn. patch size by gp 1 */ if (choice->size[3] && density1) { value[index][9] = sum2/density1; } /* calc. s4 = s.d. patch size by gp 1 */ if (choice->size[4] && density1 > 1) { mean = sum2/density1; stdv = sum22/density1 - mean*mean; if (stdv > 0) value[index][10] = sqrt(stdv); } /* calc. s5 = no. by size class 1 */ if (choice->size[5]) { value[index][11] = density2; } /* calc. s6 = no. by size class 1 by gp 1*/ if (choice->size[6]) { value[index][12] = density3; } /* calc. s7 = eff. mesh size */ if (choice->size[7]) { value[index][37] = (1.0 / sum1) * sum12; } /* calc. s8 = deg. landsc. division */ if (choice->size[8]) { value[index][38] = 1.0 - (sum12 / (sum1 * sum1)); } density1 = density2 = density3 = 0; sum1 = sum12 = sum2 = sum22 = 0.0; } return; } /************************************/ /* MOVING WINDOW CORE MEASURES */ /************************************/ void m_core (PATCH *tmp, double **value, int index) { static double sum1c=0.0, sum1e=0.0, sum2c=0.0, sum2e=0.0, sum12c=0.0, sum12e=0.0, sum22c=0.0, sum22e=0.0; static int density1c=0, density1e=0, density2c=0, density2e=0, density3c=0, density3e=0; double meanc, stdvc, meane, stdve; /* choice->core == 1 - mean core size (c1) 2 - st. dev. core size (c2) 3 - mean edge size (c3) 4 - st. dev. edge size (c4) 5 - mean core size by gp (c5) 6 - st. dev. core size by gp (c6) 7 - mean edge size by gp (c7) 8 - st. dev. edge size by gp (c8) 9 - no. by size class (c9) 10 - no. by size class by gp (c10) variables: IN: tmp = the list of patches value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: sum1c = sum of core sizes sum12c = sum of core sizes squared sum1e = sum of edge sizes sum12e = sum of edge sizes squared sum2c = sum of core sizes in gp 1 sum22c = sum of core sizes in gp 1 squared sum2e = sum of edge sizes in gp 1 sum22e = sum of edge sizes in gp 1 squared density1c = no. of cores in gp 1 density1e = no. of edges in gp 1 density2c = no. of cores in size class 1 density2e = no. of edges in size class 1 density3c = no. of cores in size class 1 and gp 1 density3e = no. of edges in size class 1 and gp 1 meanc = mean core size stdvc = standard deviation of mean core size meane = mean edge size stdve = standard deviation of mean edge size */ if (choice->core[1] || choice->core[2]) { sum1c += tmp->core; if (choice->core[2]) sum12c += tmp->core * tmp->core; } if (choice->core[3] || choice->core[4]) { sum1e += (int)tmp->edge; if (choice->core[4]) sum12e += tmp->edge * tmp->edge; } if (choice->core[5] || choice->core[6] || choice->core[7] || choice->core[8]) if (in_group(tmp->att, recl_tb[0], 0)) { if (choice->core[5] || choice->core[6]) { density1c ++; sum2c += tmp->core; if (choice->core[6]) sum22c += tmp->core * tmp->core; } if (choice->core[7] || choice->core[8]) { density1e ++; sum2e += tmp->edge; if (choice->core[8]) sum22e += tmp->edge * tmp->edge; } } if (choice->core[9]) { if (tmp->core < size_cl[2]) density2c ++; if (tmp->edge < size_cl[2]) density2e ++; } if (choice->core[10]) { if (tmp->core < size_cl[2] && in_group(tmp->att, recl_tb[0], 0)) density3c ++; if (tmp->edge < size_cl[2] && in_group(tmp->att, recl_tb[0], 0)) density3e ++; } if (!tmp->next){ /* calc. c1 = mn. core size */ if (choice->core[1] && total_patches) { value[index][13] = sum1c/total_patches; } /* calc. c2 = s.d. core size */ if (choice->core[2] && total_patches) { meanc = sum1c/total_patches; stdvc = sum12c/total_patches - meanc*meanc; if (stdvc > 0) value[index][14] = sqrt(stdvc); } /* calc. c3 = mn. edge size */ if (choice->core[3] && total_patches) { value[index][15] = sum1e/total_patches; } /* calc. c4 = s.d. edge size */ if (choice->core[4] && total_patches) { meane = sum1e/total_patches; stdve = sum12e/total_patches - meane*meane; if (stdve > 0) value[index][16] = sqrt(stdve); } /* calc. c5 = mn. core size by gp 1 */ if (choice->core[5] && density1c) { value[index][17] = sum2c/density1c; } /* calc. c6 = s.d. core size by gp 1 */ if (choice->core[6] && density1c > 1) { meanc = sum2c/density1c; stdvc = sum22c/density1c - meanc*meanc; if (stdvc > 0) value[index][18] = sqrt(stdvc); } /* calc. c7 = mn. edge size by gp 1 */ if (choice->core[7] && density1e) { value[index][19] = sum2e/density1e; } /* calc. c8 = s.d. edge size by gp 1 */ if (choice->core[8] && density1e > 1) { meane = sum2e/density1e; stdve = sum22e/density1e - meane*meane; if (stdve > 0) value[index][20] = sqrt(stdve); } /* calc. c9 = no. by size class 1 */ if (choice->core[9]) { value[index][21] = density2c; } /* calc. c10 = no. by size class 1 by gp 1 */ if (choice->core[10]) { value[index][22] = density3c; } density1c = density1e = 0; density2c = density2e = 0; density3c = density3e = 0; sum1c = sum1e = sum2c = sum2e = 0.0; sum12c = sum12e = sum22c = sum22e = 0.0; } return; } /************************************/ /* MOVING WINDOW SHAPE MEASURES */ /************************************/ void m_shape (PATCH *tmp, int way, double **value, int index) { static double sum1=0, sum12=0, sum2=0, sum22=0; static int density1=0, density2=0, density3=0; double mean, shp, stdv; /* choice->shape 1 = mean patch shape (h1) 2 = st.dev. patch shape (h2) 3 = mean patch shape by gp 1 (h3) 4 = st.dev. shape by gp 1 (h4) 5 = number by shape class 1 (h5) 6 = number by shape class 1 by gp 1 (h6) Variables: IN: tmp = the list of patches way = shape index choice (see shp) value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: shp = P/A (way=1) CPA (way=2) RCC (way=3) sum1 = sum of shape indices sum12 = sum of shape indices squared sum2 = sum of shape indices in gp 1 sum22 = sum of shape indices in gp 1 squared density1 = no. of patches in gp 1 density2 = no. of patches in shape class 1 density3 = no. of patches in shape class 1 and gp 1 */ if (way == 1) { if (tmp->area) shp = tmp->perim/tmp->area; else shp = 0.0; } else if (way == 2) { if (tmp->area) shp = (0.282 * tmp->perim)/sqrt(tmp->area); else shp = 0.0; } else { if (tmp->long_axis) shp = 2 * sqrt(tmp->area/PI)/tmp->long_axis; else shp = 0.0; } if (choice->shape[1] || choice->shape[2]) { sum1 += shp; if (choice->shape[2]) sum12 += shp * shp; } if (choice->shape[3] || choice->shape[4]) { if (in_group(tmp->att, recl_tb[0], 0)) { density1 ++; sum2 += shp; if (choice->shape[4]) sum22 += shp * shp; } } if (choice->shape[5] || choice->shape[6]) if ((way == 1 && (shp < shape_PA[2] && shp >= shape_PA[1])) || (way == 2 && (shp < shape_CPA[2] && shp >= shape_CPA[1])) || (way == 3 && (shp < shape_RCC[2] && shp >= shape_RCC[1]))) { if (choice->shape[5]) density2 ++; else if (in_group(tmp->att, recl_tb[0], 0)) density3 ++; } if (!tmp->next){ /* calc. h1=mn. patch shape */ if (choice->shape[1] && total_patches) value[index][23] = sum1/total_patches; /* calc. h2=s.d. patch shape */ if (choice->shape[2] && total_patches > 1) { mean = sum1/total_patches; stdv = sum12/total_patches - mean*mean; if (stdv > 0) value[index][24] = sqrt(stdv); } /* calc. h3=mn. patch shape in gp 1 */ if (choice->shape[3] && density1) value[index][25] = sum2/density1; /* calc. h4=s.d. patch shape in gp 1 */ if (choice->shape[4] && density1 > 1) { mean = sum2/density1; stdv = sum22/density1 - mean*mean; if (stdv > 0) value[index][26] = sqrt(stdv); } /* calc. h5=no. in shape index class 1 */ if (choice->shape[5]) value[index][27] = density2; /* calc. h6=no. in shape index class 1 & gp 1 */ if (choice->shape[6]) value[index][28] = density3; density1 = density2 = density3 = 0; sum1 = sum12 = sum2 = sum22 = 0.0; } return; } /******************************************/ /* COMPUTE AND SAVE BOUNDARY COMPLEXITY */ /* MEASURES FOR MOVING WINDOW ANALYSIS */ /******************************************/ void m_boundary (PATCH *tmp, double **value, int index) { static double sumomega=0.0, sumomega2=0.0; static int sumtwist=0, sumtwist2=0; register int i; double meantwist=0.0, stdvtwist=0.0, meanomega=0.0, stdvomega=0.0; /* Variables: IN: tmp = the next patch in the patch list value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: meantwist = mean twist number (n1) stdvtwist = st. dev. twist number (n2) meanomega = mean omega index (n3) stdvomega = st. dev. omega index (n4) GLOBAL: total_patches = tot. no. patches in sampling area (fm. trace.c) */ sumtwist += tmp->twist; sumtwist2 += tmp->twist * tmp->twist; sumomega += tmp->omega; sumomega2 += tmp->omega * tmp->omega; if (!tmp->next) { if (choice->boundary[1] || choice->boundary[2]) { meantwist = (double)sumtwist/total_patches; stdvtwist = (double)sumtwist2/total_patches - meantwist*meantwist; if (stdvtwist > 0.0) stdvtwist = sqrt(stdvtwist); else stdvtwist = 0.0; } if (choice->boundary[3] || choice->boundary[4]) { meanomega = (double)sumomega/total_patches; stdvomega = (double)sumomega2/total_patches - meanomega*meanomega; if (stdvomega > 0.0) stdvomega = sqrt(stdvomega); else stdvomega = 0.0; } /* write n1=mean twist number */ if (choice->boundary[1]) { value[index][29] = meantwist; } /* write n2=st. dev. twist number */ if (choice->boundary[2]) { value[index][39] = stdvtwist; } /* write n3=mean omega index */ if (choice->boundary[3]) { value[index][40] = meanomega; } /* write n4=st. dev. omega index */ if (choice->boundary[4]) { value[index][41] = stdvomega; } sumtwist = sumtwist2 = meantwist = stdvtwist = 0.0; sumomega = sumomega2 = meanomega = stdvomega = 0.0; } return; } /************************************/ /* MOVING WINDOW PERIMETER MEASURES */ /************************************/ void m_perim (PATCH *tmp, double **value, int index) { static double sum1=0.0, sum12=0.0, sum2=0.0, sum22=0.0; static int density=0; double mean, stdv; /* choice->perim == 1 - sum of perimeters (p1) 2 - mean perimeter (p2) 3 - st. dev. perimeters (p3) 4 - sum of perimeters in gp 1 (p4) 5 - mean perimeter in gp 1 (p5) 6 - st. dev. perimeters in gp 1 (p6) Variables: IN: tmp = the list of patches value = buffer containing a full row of the results of the moving window if a moving window map, otherwise 0 index = number of the column in the moving window INTERNAL: sum1 = sum of patch perimeters sum12 = sum of patch perimeters squared sum2 = sum of patch perimeters in gp 1 sum22 = sum of patch perimeters in gp 1 squared density = no. of patches in gp 1 */ if (choice->perim[1] || choice->perim[2] || choice->perim[3]) { sum1 += tmp->perim; if (choice->perim[3]) sum12 += tmp->perim * tmp->perim; } if (choice->perim[4] || choice->perim[5] || choice->perim[6]) { if (in_group(tmp->att, recl_tb[0], 0)){ sum2 += tmp->perim; if (choice->perim[5] || choice->perim[6]){ density ++; if (choice->perim[6]) sum22 += tmp->perim * tmp->perim; } } } if (!tmp->next){ /* calc. p1 = sum of perims. */ if (choice->perim[1]) value[index][30] = sum1; /* calc. p2 = mn. perim. */ if (choice->perim[2] && total_patches) value[index][31] = sum1/total_patches; /* calc. p3 = s.d. perim. */ if (choice->perim[3] && total_patches > 1){ mean = sum1/total_patches; stdv = sum12/total_patches - mean*mean; if (stdv > 0) value[index][32] = sqrt(stdv); } /* calc. p4 = sum of perims. in gp 1 */ if (choice->perim[4]) value[index][33] = sum2; /* calc. p5 = mn. perim. in gp 1 */ if (choice->perim[5] && density) value[index][34] = sum2/density; /* calc. p6 = s.d. perims. in gp 1 */ if (choice->perim[6] && density > 1){ mean = sum2/density; stdv = sum22/density - mean*mean; if (stdv > 0) value[index][35] = sqrt(stdv); } density = 0; sum1 = sum12 = sum2 = sum22 = 0.0; } return; } double eu_d (double x1, double y1, double x2, double y2) { return( sqrt( (x1-x2)*(x1-x2) + (y1-y2)*(y1-y2)) ); }