[GRASS-SVN] r51216 - in grass-addons/grass6/raster: . r.cva
svn_grass at osgeo.org
svn_grass at osgeo.org
Sat Mar 31 03:58:03 EDT 2012
Author: amuriy
Date: 2012-03-31 00:58:03 -0700 (Sat, 31 Mar 2012)
New Revision: 51216
Added:
grass-addons/grass6/raster/r.cva/
grass-addons/grass6/raster/r.cva/Makefile
grass-addons/grass6/raster/r.cva/README
grass-addons/grass6/raster/r.cva/TODO
grass-addons/grass6/raster/r.cva/authors
grass-addons/grass6/raster/r.cva/bugs
grass-addons/grass6/raster/r.cva/config.h
grass-addons/grass6/raster/r.cva/count_visible_cells.c
grass-addons/grass6/raster/r.cva/count_visible_cells.h
grass-addons/grass6/raster/r.cva/cumulate_visible_cells.c
grass-addons/grass6/raster/r.cva/cumulate_visible_cells.h
grass-addons/grass6/raster/r.cva/delete3.c
grass-addons/grass6/raster/r.cva/delete_lists.c
grass-addons/grass6/raster/r.cva/delete_lists.h
grass-addons/grass6/raster/r.cva/delete_one_list.c
grass-addons/grass6/raster/r.cva/description
grass-addons/grass6/raster/r.cva/description.html
grass-addons/grass6/raster/r.cva/global_vars.h
grass-addons/grass6/raster/r.cva/init_segment_file.c
grass-addons/grass6/raster/r.cva/init_segment_file.h
grass-addons/grass6/raster/r.cva/line_of_sight.c
grass-addons/grass6/raster/r.cva/line_of_sight.h
grass-addons/grass6/raster/r.cva/main.c
grass-addons/grass6/raster/r.cva/make_list3.c
grass-addons/grass6/raster/r.cva/make_point3.c
grass-addons/grass6/raster/r.cva/point.h
grass-addons/grass6/raster/r.cva/point_sample.c
grass-addons/grass6/raster/r.cva/point_sample.h
grass-addons/grass6/raster/r.cva/pts_elim3.c
grass-addons/grass6/raster/r.cva/pts_elim3.h
grass-addons/grass6/raster/r.cva/r.cva.tmp.html
grass-addons/grass6/raster/r.cva/radians.h
grass-addons/grass6/raster/r.cva/random_int.c
grass-addons/grass6/raster/r.cva/random_sample.c
grass-addons/grass6/raster/r.cva/random_sample.h
grass-addons/grass6/raster/r.cva/scan_all_cells.c
grass-addons/grass6/raster/r.cva/scan_all_cells.h
grass-addons/grass6/raster/r.cva/segment3.c
grass-addons/grass6/raster/r.cva/segment3.h
grass-addons/grass6/raster/r.cva/segment_file.c
grass-addons/grass6/raster/r.cva/segment_file.h
grass-addons/grass6/raster/r.cva/systematic_sample.c
grass-addons/grass6/raster/r.cva/systematic_sample.h
grass-addons/grass6/raster/r.cva/version
grass-addons/grass6/raster/r.cva/zufall.c
grass-addons/grass6/raster/r.cva/zufall.h
grass-addons/grass6/raster/r.cva/zufalli.c
Modified:
grass-addons/grass6/raster/Makefile
Log:
<r.cva> added to grass6/raster; grass6/raster/Makefile edited
Modified: grass-addons/grass6/raster/Makefile
===================================================================
--- grass-addons/grass6/raster/Makefile 2012-03-30 23:00:09 UTC (rev 51215)
+++ grass-addons/grass6/raster/Makefile 2012-03-31 07:58:03 UTC (rev 51216)
@@ -1,31 +1,53 @@
MODULE_TOPDIR = ..
-# r.convert
-# r.roughness
SUBDIRS = \
LandDyn \
mcda \
+ r.area \
+ r.basin \
r.boxcount \
r.boxcount.sh \
r.burn.frict \
r.clim \
+ r.clump2 \
+ r.colors.out_sld \
+ r.colors.out_vtk \
r.colors.quantiles \
+ r.colors.tools \
+ r.convergence \
+ r.convert \
+ r.csr \
+ r.cva \
r.denoise \
+ r.diversity \
+ r.flip \
+ r.fuzzy \
+ r.fuzzy.logic \
+ r.fuzzy.system \
r.game_of_life \
+ r.hazard.flood \
+ r.in.swisstopo \
r.inund.fluv \
r.in.xyz.auto \
r.ipso \
r.isoregions \
+ r.maxent.lambdas \
r.out.gmap \
r.out.gmt \
r.out.gmt2 \
- r.out.kml \
+ r.out.kap_template \
+ r.out.kml \
+ r.out.maxent_swd \
r.out.netcdf \
r.pack \
r.pi \
r.prominence \
r.rast4d \
+ r.roughness \
+ r.seg \
r.soils.texture \
+ r.stack \
+ r.stream.angle \
r.stream.basins \
r.stream.del \
r.stream.distance \
@@ -37,11 +59,14 @@
r.surf.volcano \
r.terracost \
r.terraflow.tools \
- r.univar2.zonal \
+ r.threshold \
+ r.univar.zonal \
+ r.unpack \
r.viewshed \
r.wf \
- r.xtent
+ r.xtent
+
include $(MODULE_TOPDIR)/include/Make/Dir.make
default: htmldir
Added: grass-addons/grass6/raster/r.cva/Makefile
===================================================================
--- grass-addons/grass6/raster/r.cva/Makefile (rev 0)
+++ grass-addons/grass6/raster/r.cva/Makefile 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,15 @@
+MODULE_TOPDIR = ../..
+
+PGM = r.cva
+
+LIBES = $(GISLIB) $(DATETIMELIB) $(SEGMENTLIB) $(VECTLIB) $(DBMILIB) $(RASTERLIB)
+
+DEPENDENCIES = $(RASTERDEP) $(VECTDEP) $(DBMIDEP) $(GISDEP)
+
+EXTRA_CFLAGS = $(VECT_CFLAGS)
+EXTRA_INC = $(VECT_INC)
+
+
+include $(MODULE_TOPDIR)/include/Make/Module.make
+
+default: cmd
Added: grass-addons/grass6/raster/r.cva/README
===================================================================
--- grass-addons/grass6/raster/r.cva/README (rev 0)
+++ grass-addons/grass6/raster/r.cva/README 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,2 @@
+This is a GRASS 6 extension that contains additional functionality
+for cumulative viewshed analysis.
Added: grass-addons/grass6/raster/r.cva/TODO
===================================================================
--- grass-addons/grass6/raster/r.cva/TODO (rev 0)
+++ grass-addons/grass6/raster/r.cva/TODO 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,121 @@
+HAS THE PATCH FOR r.los BEEN APPLIED IN GRASS 6.1 CVS (DOES IT DO ANYTHING IN r.los, anyway?)
+
+
+BUGS:
+ dist=sqrt(del_x * del_x + del_y * del_y)*window.ns_res;
+ in pts_elim3.c seems to indicate, that this module always
+ assumes N/S res to be = E/W res for distance calculations!
+
+
+CHECK:
+ Check that OFFSETB now gives the desired results.
+ (pts_elim3.c)
+
+ Check that target_elev is really just handled as an alias
+ for OFFSETB (main.c)
+
+ Check correctness of Earth Curvature Correction
+
+ point_sample.c: L149: it seems that db_init_string (&sql) is called iteratively
+ in this loop w/o releasing any memory resulting in a mem leak
+
+PARAMS:
+ the default operation should NOT be reflex viewsheds!
+ Some of the parameters are redundant (obs height, tgt height)
+ or should at least be more clearly named (e.g.: "global OFFSET A")
+
+ make sure the word "sites" disappeards from all titles and parameter names
+
+ unify writing of all parameter names: capitalise first word, NO period at end
+
+ - parameters must be declared GRASS 6.0.0 style, so that the right
+ selector buttons will be provided!
+
+
+DOCS:
+ HTML man page for r.prominence
+
+ state, that missing values are not allowed in site parameters!
+
+ A separate section to clarify:
+ clarify NULL handling: discarded in systematic samples
+ not valid in random samples
+ discarded in site sampling
+ discarded in total sampling
+ achieved sample size discounts NULL cells
+ what does the "convert -1 to NULL" option do?
+ state that NULL cells do not count for
+ the reflexive viewshed!
+
+
+
+
+* Store visibility count of sites in an attribute in the input vector map,
+(count case only)
+if attribute with name VCOUNT exists (or create such an attribute).
+
+Add a "z" factor for height scaling.
+
+Convenience: user supplies a map of obstructing objects as vector objects
+with a field HEIGHT that gives relative height of object. This is
+rasterised into tmp1 map, the height is added to dem and the result
+output as tmp2 and fed into r.cva instead of the original dem.
+
+Rename "all" option to "total", as it is less ambiguous. Change default option
+to "sites".
+
+Implement fuzzy viewsheds that can attribute for unknown variations in terrain
+heights, due to error or vegetation. Add a random offset within a user-defined
+range to observer point,
+target point and all blocking points. Do this n times and the final result
+is the average visibility over all trials.
+Offset size can be determined from error specifications of DEM and
+average height of vegetation (?).
+
+WE MUST BE ABLE TO MASK OUT CERTAIN AREAS FROM VEGETATION OFFSET, WHICH
+ARE CLEARED AREAS (?).
+
+Alternative, systematic fuzzy viewsheds:
+The user supplies a set of maps with different configurations of spatial
+obstacles. All of these configurations are taken into account and the final
+result is the average visibility over all configurations (e.g. 4 of 5 times=80%).
+
+* Implement probabilistic viewsheds: user defines a distance D1, at which objects
+of interest are no longer clearly identifiable from the obersvation Point P and
+distance D2, at which they are no longer observable at all.
+The program computes a corresponding interpolation curve from "1" (P) to "0.5"
+(P+D1) and on to "0" (P+D2). Instead of "0"/"1" viewsheds,
+we can then have a continuous probability. For the cumulative case, we would
+multiply probabilities instead of adding. We can also have functions that
+relate to viewing angles. Secondly, use a visibility map with values in the
+range of 0..1 that encode probability of objects being spotted on a certain
+type of terrain. With these functions, we are moving away from terrain
+visibility and towards hypothetical object visibility.
+This needs a float output map. Two alternatives for the cumulative case:
+directly work on an FCELL map (SLOW, but can link in cache!) or create another
+segment file of FCELL type. Option one should work, as operations are fairly
+localised.
+
+The terrain adjacent to a viewing location also has a major impact on visibility
+target locations! This can be catered for by adding the height of obstructing
+objects to the terrain.
+
+Another field has to be added to the linked list of points storing
+its probability in point.h [DONE: probability]
+
+This probability has to be updated according to the function linking it to
+distance and viewing angle:
+Add a parameter probability to function definition in make_point3.c and
+make_point3.h.
+Add probability calculation function to make_list3.c, at the end, but before
+call to make_point().
+
+It has also got to be updated using floats in ONE map that encodes visibility
+for e.g. different types of terrain (e.g. green objects are harder to spot
+on grassland than on rock surface).
+
+For the two (*) above: files cumulate_visible_cells.c (-f option) and
+count_visible_cells.c as a starting point.
+
+Is the mask option really needed or can the same effect be achieved by setting
+the GRASS global mask?
Added: grass-addons/grass6/raster/r.cva/authors
===================================================================
--- grass-addons/grass6/raster/r.cva/authors (rev 0)
+++ grass-addons/grass6/raster/r.cva/authors 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,3 @@
+#comment
+Mark Lake (University College London)<br>
+Benjamin Ducke (University of Kiel, Germany)
Added: grass-addons/grass6/raster/r.cva/bugs
===================================================================
--- grass-addons/grass6/raster/r.cva/bugs (rev 0)
+++ grass-addons/grass6/raster/r.cva/bugs 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,6 @@
+# Please describe all known bugs and shortcomings of the current implementation
+# of your extension.
+<strong>r.cva:</strong><br>
+Earth curvature correction has not been sufficiently tested.<br>
+The 'offsetb' option (=target height) has not been sufficiently tested.<br>
+Please feel free to test and report!<br>
Added: grass-addons/grass6/raster/r.cva/config.h
===================================================================
--- grass-addons/grass6/raster/r.cva/config.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/config.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,207 @@
+/* include/config.h. Generated automatically by configure. */
+
+/*
+ * config.h.in
+ */
+
+#ifndef _config_h
+#define _config_h
+
+#define GDEBUG 1
+
+/* define if standard C headers exist: headers are
+ * stdlib.h, stdarg.h, string.h, and float.h */
+
+/* define if curses.h exists */
+/* #undef HAVE_CURSES_H */
+
+/* define if keypad in lib[n]curses */
+/* #undef HAVE_KEYPAD */
+
+/* define if limits.h exists */
+#define HAVE_LIMITS_H 1
+
+/* define if termio.h exists */
+#define HAVE_TERMIO_H 1
+
+/* define if termios.h exists */
+#define HAVE_TERMIOS_H 1
+
+/* define if unistd.h exists */
+#define HAVE_UNISTD_H 1
+
+/* define if values.h exists */
+#define HAVE_VALUES_H 1
+
+/* define if zlib.h exists */
+#define HAVE_ZLIB_H 1
+
+/* define if sys/ioctl.h exists */
+#define HAVE_SYS_IOCTL_H 1
+
+/* define if sys/mtio.h exists */
+#define HAVE_SYS_MTIO_H 1
+
+/* define if sys/resource.h exists */
+#define HAVE_SYS_RESOURCE_H 1
+
+/* define if sys/time.h exists */
+#define HAVE_SYS_TIME_H 1
+
+/* define if time.h and sys/time.h can be included together */
+#define TIME_WITH_SYS_TIME 1
+
+/* define if sys/timeb.h exists */
+#define HAVE_SYS_TIMEB_H 1
+
+/* define if sys/types.h exists */
+#define HAVE_SYS_TYPES_H 1
+
+/* define if sys/utsname.h exists */
+#define HAVE_SYS_UTSNAME_H 1
+
+/* define if g2c.h exists */
+#define HAVE_G2C_H 1
+
+/* define if f2c.h exists */
+/* #undef HAVE_F2C_H */
+
+/* define gid_t type */
+/* #undef gid_t */
+
+/* define off_t type */
+/* #undef off_t */
+
+/* define uid_t type */
+/* #undef uid_t */
+
+/* Define the return type of signal handlers */
+#define RETSIGTYPE void
+
+/* define curses.h WINDOW structure component */
+#define CURSES_MAXY @CURSES_MAXY@
+
+/* define if ftime() exists */
+#define HAVE_FTIME 1
+
+/* define if gethostname() exists */
+#define HAVE_GETHOSTNAME 1
+
+/* define if gettimeofday() exists */
+#define HAVE_GETTIMEOFDAY 1
+
+/* define if lseek() exists */
+#define HAVE_LSEEK 1
+
+/* define if time() exists */
+#define HAVE_TIME 1
+
+/* define if uname() exists */
+#define HAVE_UNAME 1
+
+/* define if seteuid() exists */
+#define HAVE_SETEUID 1
+
+/* define if setpriority() exists */
+#define HAVE_SETPRIORITY 1
+
+/* define if setreuid() exists */
+#define HAVE_SETREUID 1
+
+/* define if setruid() exists */
+/* #undef HAVE_SETRUID */
+
+/* define if setpgrp() takes no argument */
+#define SETPGRP_VOID 1
+
+/* define if drand48() exists */
+#define HAVE_DRAND48 1
+
+/* define if postgres is to be used */
+/* #undef HAVE_POSTGRES */
+
+/* define if OGR is to be used */
+#define HAVE_OGR 1
+
+/* define if postgres client header exists */
+/* #undef HAVE_LIBPQ_FE_H */
+
+/* define if PQcmdTuples in lpq */
+/* #undef HAVE_PQCMDTUPLES */
+
+/* define if ODBC exists */
+/* #undef HAVE_SQL_H */
+
+/* define if PNG support exists */
+/* #undef HAVE_PNG */
+
+/* define if jpeglib.h exists */
+/* #undef HAVE_JPEGLIB_H */
+
+/* define if fftw.h exists */
+/* #undef HAVE_FFTW_H */
+
+/* define if dfftw.h exists */
+/* #undef HAVE_DFFTW_H */
+
+/* define if BLAS exists */
+/* #undef HAVE_LIBBLAS */
+
+/* define if LAPACK exists */
+/* #undef HAVE_LIBLAPACK */
+
+/* define if dbm.h exists */
+/* #undef HAVE_DBM_H */
+
+/* define if readline exists */
+/* #undef HAVE_READLINE_READLINE_H */
+
+/* define if ft2build.h exists */
+/* #undef HAVE_FT2BUILD_H */
+
+/* Whether or not we are using G_socks for display communications */
+#define USE_G_SOCKS 1
+
+/* define if X is disabled or unavailable */
+/* #undef X_DISPLAY_MISSING */
+
+/* define if libintl.h exists */
+#define HAVE_LIBINTL_H 1
+
+/* define if iconv.h exists */
+#define HAVE_ICONV_H 1
+
+/* define if NLS requested */
+/* #undef USE_NLS */
+
+/* define if <GL/GLwMDrawA.h> exists */
+/* #undef HAVE_GL_GLWMDRAWA_H */
+
+/* define if <X11/GLw/GLwMDrawA.h> exists */
+/* #undef HAVE_X11_GLW_GLWMDRAWA_H */
+
+/* define if TclTk exists */
+/* #undef HAVE_TCLTK */
+
+/* define if putenv() exists */
+#define HAVE_PUTENV 1
+
+/* define if setenv() exists */
+#define HAVE_SETENV 1
+
+/* define if glXCreatePbuffer exists */
+/* #undef HAVE_PBUFFERS */
+
+/* define if glXCreateGLXPixmap exists */
+/* #undef HAVE_PIXMAPS */
+
+/*
+ * configuration information solely dependent on the above
+ * nothing below this point should need changing
+ */
+
+#if defined(HAVE_VALUES_H) && !defined(HAVE_LIMITS_H)
+#define INT_MIN -MAXINT
+#endif
+
+#endif /* _config_h */
Added: grass-addons/grass6/raster/r.cva/count_visible_cells.c
===================================================================
--- grass-addons/grass6/raster/r.cva/count_visible_cells.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/count_visible_cells.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,108 @@
+/****************************************************************
+ * count_visible_cells.c
+ *
+ * Counts the number of cells visible from the current
+ * viewpoint, i.e. counts the number of cells on the
+ * 16 segment lists. Writes the current viewpoint number
+ * to the temporary segment file after a cell has been
+ * counted in order prevent cells which occur on more than
+ * one list from being counted more than once. Writes the
+ * number of visible cells to the output segment file
+ * and updates the minimum and maximum number of visible
+ * cells.
+ *
+ ****************************************************************/
+
+/* Written by Mark Lake on 1/3/96 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Count the number of cells visible from the current
+ * viewpoint.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) r.los does not perform this function.
+ *
+ * Called by:
+ *
+ * 1) scan_all_cells
+ * 2) random_sample
+ * 3) systematic_sample
+ * 4) point_sample
+ *
+ * Calls:
+ *
+ * 1) None.
+ */
+
+#include <grass/segment.h>
+#include <grass/gis.h>
+
+#include "config.h"
+#include "point.h"
+#include "global_vars.h"
+
+#define NEXT_SEARCH_PT SEARCH_PT->next
+
+void count_visible_cells (CELL cell_no, int row_viewpt, int col_viewpt,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p, struct point *heads[])
+
+/* struct point *heads[] is same as "struct point **heads" */
+
+{
+ int segment_no;
+ char *value_read, *value_to_write;
+ CELL visible_cells;
+ CELL data_read, data_to_write;
+ struct point *SEARCH_PT;
+ int mask;
+
+ visible_cells = 0;
+ value_read = (char*) &data_read;
+ value_to_write = (char*) &data_to_write;
+ mask = 1;
+
+
+ /* Count cells visible from current viewpoint */
+
+ for (segment_no = 1; segment_no <= 16; segment_no++)
+ {
+ data_to_write = cell_no;
+ SEARCH_PT = heads[segment_no-1];
+ while(SEARCH_PT != NULL)
+ {
+ segment_get (seg_out_1_p, value_read,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+ if (data_read == 0) /* if not already counted or deleted */
+ {
+ visible_cells++;
+
+
+ /* Write cell_no to temporary segment file to show
+ that cell has already been counted */
+
+ segment_put (seg_out_1_p, value_to_write,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+ }
+ SEARCH_PT = NEXT_SEARCH_PT;
+ }
+ }
+
+ /* Include viewpoint itself if it is not a masked destination */
+
+ if (patt_flag)
+ {
+ segment_get (seg_patt_p, value_read, row_viewpt, col_viewpt);
+ mask = (int) data_read;
+ if (mask)
+ visible_cells += 1; /* include viewpoint */
+ }
+
+ value_to_write = (char*) &visible_cells;
+ segment_put (seg_out_2_p, value_to_write, row_viewpt, col_viewpt);
+}
+
Added: grass-addons/grass6/raster/r.cva/count_visible_cells.h
===================================================================
--- grass-addons/grass6/raster/r.cva/count_visible_cells.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/count_visible_cells.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,10 @@
+#ifndef __COUNT_VISIBLE_CELLS_H__
+#define __COUNT_VISIBLE_CELLS_H__
+
+#include <grass/segment.h>
+
+void count_visible_cells (CELL cell_no, int row_viewpt, int col_viewpt,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p, struct point *heads[]);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/cumulate_visible_cells.c
===================================================================
--- grass-addons/grass6/raster/r.cva/cumulate_visible_cells.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/cumulate_visible_cells.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,125 @@
+/****************************************************************
+ * cumulate_visible_cells.c
+ *
+ * Marks the cells visible from the current viewpoint
+ * by adding 1 to the values of visible cells on the
+ * output segment file. Writes the current viewpoint
+ * number to the temporary segment file to prevent
+ * cells which occur on more than one segment list from
+ * being marked more than once for a given viewpoint.
+ * The end effect is to produce a map showing from how
+ * viewpoints a given cell can be seen. Updates the
+ * minimum and maximum number of visible cells.
+ *
+ ****************************************************************/
+
+/* Written by Mark Lake on 1/3/96 to:-
+ *
+ * 1) Cumulate the number of visible cells.
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * Modifications to r.los:-
+ *
+ * 1) r.los does not perform this function.
+ *
+ * Called by:-
+ *
+ * Called by:
+ *
+ * 1) scan_all_cells
+ * 2) random_sample
+ * 3) systematic_sample
+ * 4) point_sample
+ *
+ * Calls:
+ *
+ * 1) None.
+ */
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "global_vars.h"
+
+#include "point.h"
+
+#define NEXT_SEARCH_PT SEARCH_PT->next
+
+void cumulate_visible_cells (CELL cell_no,
+ int row_viewpt, int col_viewpt,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p,
+ struct point *heads[])
+
+/* struct point *heads[] is same as "struct point **heads" */
+
+{
+ int segment_no;
+ int mask;
+ char *value_read, *value_to_write;
+ CELL visible_cells;
+ CELL data_read, data_to_write;
+ struct point *SEARCH_PT;
+
+ /* Count cells visible from current viewpoint */
+ visible_cells = 0;
+ value_read = (char*) &data_read;
+ value_to_write = (char*) &data_to_write;
+ mask = 1;
+
+ for (segment_no = 1; segment_no <= 16; segment_no++) {
+ SEARCH_PT = heads[segment_no-1];
+
+#ifdef DEBUG
+ printf ("\nCumulate_cells SEARCH_PT = %p", SEARCH_PT);
+ fflush (stdout);
+#endif
+
+ while(SEARCH_PT != NULL) {
+ segment_get (seg_out_1_p, value_read,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+ /* if not already counted or deleted */
+ if (data_read < cell_no ) {
+ visible_cells++;
+
+ /* Add one to cell on output segment file */
+ segment_get (seg_out_2_p, value_read,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+ if (data_read == background) /* NULL value */
+ data_read = 0;
+
+ data_to_write = data_read + 1;
+ segment_put (seg_out_2_p, value_to_write,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+
+ /* Write cell_no to temporary segment file to show
+ that cell has already been counted */
+ data_to_write = cell_no;
+ segment_put (seg_out_1_p, value_to_write,
+ row_viewpt - SEARCH_PT->y,
+ col_viewpt + SEARCH_PT->x);
+ }
+ SEARCH_PT = NEXT_SEARCH_PT;
+ }
+ } /* END (iterate through segments */
+
+ /* Include viewpoint itself if it is not a masked destination */
+ if (patt_flag) {
+ segment_get (seg_patt_p, value_read, row_viewpt, col_viewpt);
+ mask = (int) data_read;
+ }
+
+ if (mask) {
+ segment_get (seg_out_2_p, value_read, row_viewpt, col_viewpt);
+ if (data_read == background)
+ data_read = 0;
+
+ data_to_write = data_read + 1;
+ segment_put (seg_out_2_p, value_to_write, row_viewpt, col_viewpt);
+ }
+}
Added: grass-addons/grass6/raster/r.cva/cumulate_visible_cells.h
===================================================================
--- grass-addons/grass6/raster/r.cva/cumulate_visible_cells.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/cumulate_visible_cells.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,12 @@
+#ifndef __CUMULATE_VISIBLE_CELLS_H__
+#define __CUMULATE_VISIBLE_CELLS_H__
+
+#include <grass/segment.h>
+
+void cumulate_visible_cells (CELL cell_no,
+ int row_viewpt, int col_viewpt,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p,
+ struct point *heads[]);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/delete3.c
===================================================================
--- grass-addons/grass6/raster/r.cva/delete3.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/delete3.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,106 @@
+/****************************************************************/
+/* */
+/* delete3.c */
+/* */
+/* This function detaches a point data structure from */
+/* the linked list and frees memory allocated for it. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by MWL on 14/12/95 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Detach a point structure from the linked list.
+ *
+ * Version 3.0, differs from v.1.0 as follows:-
+ *
+ * 1) Recieves cell_no and writes to segment file at deleted point.
+ * This ensures that make_point3.c can differentiate between a
+ * point that has been deleted during the current los analysis,
+ * and one which was deleted or counted during los
+ * analysis from a preceeding viewpoint.
+ *
+ * Implicated files:-
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) Writes cell_no to seg_out_p instead of value 1.
+ *
+ * Called by:-
+ *
+ * 1) hidden_point_elimination()
+ *
+ * Calls:-
+ *
+ * 1) None. */
+
+#include <stdlib.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+
+#define PT_TO_DELETE_X PT_TO_DELETE->x
+#define PT_TO_DELETE_Y PT_TO_DELETE->y
+#define PT_NEXT_TO_DELETED PT_TO_DELETE->next
+#define PT_PREVIOUS_TO_DELETED PT_TO_DELETE->previous
+#define NEXT_PT_BACK_PTR PT_TO_DELETE->next->previous
+#define PREVIOUS_PT_NEXT_PTR PT_TO_DELETE->previous->next
+
+struct point *delete(struct point *PT_TO_DELETE, struct point *head,
+ SEGMENT *seg_out_p, int row_viewpt, int col_viewpt,
+ CELL cell_no)
+{
+ CELL data;
+ char *value;
+
+ data = cell_no;
+ value = (char *) &data;
+ segment_put(seg_out_p,value,
+ row_viewpt-PT_TO_DELETE_Y,PT_TO_DELETE_X+col_viewpt);
+
+
+ /* If first and last point. This fixes a bug in r.los.
+ See pts_elim3.c for details */
+
+ if((PT_TO_DELETE==head) && (PT_NEXT_TO_DELETED == NULL))
+ {
+ head = PT_NEXT_TO_DELETED;
+ free(PT_TO_DELETE);
+ return(head);
+ }
+
+ if(PT_TO_DELETE==head) /* first one ? */
+ {
+ NEXT_PT_BACK_PTR = NULL;
+ head = PT_NEXT_TO_DELETED;
+ free(PT_TO_DELETE);
+ return(head);
+ }
+
+
+ /* If last point. This fixes a bug in r.los.
+ See pts_elim3.c for details */
+
+ if (PT_NEXT_TO_DELETED == NULL)
+ {
+ PREVIOUS_PT_NEXT_PTR = PT_NEXT_TO_DELETED;
+ free(PT_TO_DELETE);
+ return (head);
+ }
+
+ /* otherwise */
+
+ NEXT_PT_BACK_PTR = PT_PREVIOUS_TO_DELETED;
+ PREVIOUS_PT_NEXT_PTR = PT_NEXT_TO_DELETED;
+ free(PT_TO_DELETE);
+
+ return(head);
+ }
+
+/************* END OF FUNCTION "DELETE" *************************/
Added: grass-addons/grass6/raster/r.cva/delete_lists.c
===================================================================
--- grass-addons/grass6/raster/r.cva/delete_lists.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/delete_lists.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,59 @@
+/****************************************************************
+ * delete_lists.c
+ *
+ * Delete the 16 segment lists used to record the cells
+ * that are visible from each viewpoint.
+ *
+ ****************************************************************/
+
+/* Modified from r.cumlos
+ * by MWL on 1/3/96 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Delete the 16 segment lists used to record the cells
+ * that are visible from each viewpoint.
+ *
+ * Modifications to r.cumlos version 3.0:-
+ *
+ * 1) Taken out of main() and placed in function delete_lists().
+ *
+ * Modifications to r.los:-
+ *
+ * 1) r.los does not perform this function.
+ *
+ * Called by:-
+ *
+ * 1) scan_all_cells()
+ * 2) random_sample()
+ * 3) systematic_sample()
+ * 4) from_sites()
+ *
+ * Calls:-
+ *
+ * 1) delete_one_list()
+ */
+
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+
+void delete_lists (struct point *heads[])
+
+/* struct point *heads[] is same as "struct point **heads" */
+
+{
+ int segment_no;
+ void delete_one_list ();
+
+ /* delete the list of visible cells */
+
+ for (segment_no = 1; segment_no <= 16; segment_no++)
+ {
+ /* printf ("\nDeleting list for segment = %d ", segment_no);
+ fflush (stdout);*/
+ delete_one_list (heads[segment_no-1]);
+ }
+}
+
Added: grass-addons/grass6/raster/r.cva/delete_lists.h
===================================================================
--- grass-addons/grass6/raster/r.cva/delete_lists.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/delete_lists.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,8 @@
+#ifndef __DELETE_LISTS_H__
+#define __DELETE_LISTS_H__
+
+#include "point.h"
+
+void delete_lists (struct point *heads[]);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/delete_one_list.c
===================================================================
--- grass-addons/grass6/raster/r.cva/delete_one_list.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/delete_one_list.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,71 @@
+/****************************************************************/
+/* */
+/* delete_list.c */
+/* */
+/* This function frees the memory used in a segment list. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by MWL on 14/12/95 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Free the memory used in a single segment list.
+ *
+ * Implicated files:-
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) r.los does not perform this function.
+ *
+ * Called by:-
+ *
+ * 1) delete_lists()
+ *
+ * Calls:-
+ *
+ * 1) None. */
+
+#include <stdlib.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+
+
+#define PT_NEXT_TO_DELETED PT_TO_DELETE->next
+#define PT_PREVIOUS_TO_DELETED PT_TO_DELETE->previous
+#define NEXT_PT_BACK_PTR PT_TO_DELETE->next->previous
+#define PREVIOUS_PT_NEXT_PTR PT_TO_DELETE->previous->next
+
+void delete_one_list (struct point *head)
+{
+ struct point *PT_TO_DELETE;
+
+/* PT_TO_DELETE = head;
+ while (PT_TO_DELETE != NULL)
+ {
+ printf ("\nPOINT = %d ", PT_TO_DELETE); fflush (stdout);
+ PT_TO_DELETE = PT_NEXT_TO_DELETED;
+ }*/
+
+ PT_TO_DELETE = head;
+ while (PT_TO_DELETE != NULL)
+ {
+/*printf ("\nPT_TO_DELETE = %d ", PT_TO_DELETE); fflush (stdout);*/
+ if (PT_NEXT_TO_DELETED != NULL)
+ NEXT_PT_BACK_PTR = NULL;
+/*printf ("\nPT_NEXT = %d ", PT_NEXT_TO_DELETED); fflush (stdout);*/
+ head = PT_NEXT_TO_DELETED;
+ free(PT_TO_DELETE);
+/*printf ("\nDELETED");*/
+ PT_TO_DELETE = head;
+ }
+}
+
+/************* END OF FUNCTION "DELETE_LIST" *************************/
+
Added: grass-addons/grass6/raster/r.cva/description
===================================================================
--- grass-addons/grass6/raster/r.cva/description (rev 0)
+++ grass-addons/grass6/raster/r.cva/description 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,6 @@
+# This file should contain a very brief description of what your extension does.
+# All lines will be output as you write them in here, including linebreaks.
+# commented lines and parts of lines after a comment mark ('#') will be skipped.
+This is a GRASS 6.1.0 (and above) extension to perform cumulative viewshed analysis.
+
+If you have GIS Manager installed (which is true per default for a GRASS 6.1.0 installation), you should be able to locate the new modules in a submenu called "Viewshed Analysis" under the "Xtns" menu in the main menu bar.
Added: grass-addons/grass6/raster/r.cva/description.html
===================================================================
--- grass-addons/grass6/raster/r.cva/description.html (rev 0)
+++ grass-addons/grass6/raster/r.cva/description.html 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,428 @@
+<TITLE>r.cva</TITLE>
+
+<body bgcolor="white">
+
+<H2>NAME</H2>
+
+<EM><b>r.cva</b></EM> - Cumulative viewshed analysis program.
+<BR>
+<EM>(GRASS Raster Program)</EM>
+
+<H2>SYNOPSIS</H2>
+<B>r.cva</B>
+<BR>
+<B>r.cva help</B>
+<BR>
+<B>r.cva</B>
+[<B>-afmnoqrs</B>]
+<B>input=</B><EM>name</EM>
+[<B>target_mask=</B><EM>name</EM>]
+[<B>viewpt_mask=</B><EM>name</EM>]
+[<B>sites=</B><EM>name</EM>]
+<B>output=</B><EM>name</EM>
+[<B>obs_elev=</B><EM>value</EM>]
+[<B>target_elev=</B><EM>value</EM>]
+[<B>max_dist=</B><EM>value</EM>]
+[<B>seed=</B><EM>value</EM>]
+[<B>sample=</B><EM>value</EM>]
+<B>type=</B><EM>value</EM>
+
+<H2>QUICKSTART</H2>
+To perform the most common cumulative viewshed analysis (for each cell of
+the input elevation raster map: count the number of observer sites from which it
+can be seen) type in the following:
+
+<p>
+
+<B>r.cva</B> <B>-f</B> <B>input=</B><EM>elevation</EM> <B>sites=</B><EM>observers</EM> <B>output=</B><EM>result</EM> <B>type=</B><EM>sites</EM>
+<P>
+'elevation' can be a raster map containing elevation values in CELL, FCELL or DCELL format. 'observers' must a map containing vector points.
+'result' will be an output map in CELL format containing the cumulative visibility counts for each cell.
+
+<A NAME = "description"><H2>DESCRIPTION</H2></A>
+
+<EM>r.cva</EM> generates a cumulative viewshed map. When run without
+the <EM>-f</EM> flag each viewpoint is marked with the number of cells in its
+viewshed (i.e. that can be seen from it). In this mode the only cells
+that contain genuine data in the output map are the viewpoints.
+Alternatively, if the <EM>-f</EM> flag is set then each cell in the
+<EM>output</EM> map is marked with the number of viewpoints from which
+it can be seen. In either case, viewpoints may be chosen in one of
+four ways. The first (`all') treats every cell as a viewpoint. This
+is so computationally demanding that it is unlikely to be feasible
+unless either the <EM>input</EM> map or <EM>max_dist</EM> is very
+small. The second (`systematic') uses a systematic sample of cells
+arranged on a rectangular grid. The third (`random') uses a random
+sample. By default the random sampling algorithm does not allow
+replacement (i.e. a cell can be treated as a viewpoint only once), but
+this behaviour can be altered by setting the <EM>-r</EM> flag. The
+fourth option (`sites') uses the points in the specified <EM>vector</EM>
+map. The options `systematic' and `random' allow the user to specify
+any sampling frequency between 0% and 100%.
+
+<p>
+
+The cells used in viewshed analysis can be constrained by one or both
+of two masks. Both masks should be binary raster maps in which the
+areas of interest are coded `1' and all other cells `0'.
+<EM>target_mask</EM> determines which cells are targets: a cell which
+is not a target will be returned as not visible from a given viewpoint
+even if it is actually visible according to the line-of-sight analysis
+(i.e. this mask performs the same function as <EM>patt_map</EM> in
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>). <EM>viewpt_mask</EM>
+constrains which cells may be treated as viewpoints. For example,
+suppose that the raster map `no_sea' codes sea (elevation = 0) as `0'
+and land (elevation > 0) as `1'. In this case, specifying `no_sea' as
+the <EM>target_mask</EM> would cause <EM>r.cva</EM> to count cells on
+land that are visible from any viewpoint. Alternatively, specifying
+`no_sea' as the <EM>viewpt_mask</EM> would count any cells visible
+from viewpoints on land. Finally, specifying `no_sea' as both the
+<EM>target_mask</EM> and <EM>viewpt_mask</EM> would count cells on
+land that are visible from viewpoints on land.
+
+<p>
+
+The geometry of the potential line of sight between each pair of cells
+may be altered by means of two offsets, <EM>obs_elev</EM> and
+<EM>target_elev</EM>. The height <EM>obs_elev</EM> is added to the
+viewpoint and the height <EM>target_elev</EM> to each target cell.
+Note that by swapping these values it is possible to produce maps
+showing the number of cells from which a viewpoint is visible, or if
+the <EM>-f</EM> flag is set, the number of viewpoints which it is
+possible to see from each map cell.
+
+<p>
+
+A particular set of viewpoint locations in a random sample can be
+repeated by specifying the same value for <EM>seed</EM>.
+Alternatively, different values can be used to produce independent
+samples for statistical testing.
+
+<p>
+
+The program can be run either non-interactively or interactively. To
+run <EM>r.cva</EM> non-interactively the user must minimally specify
+the <B>input</B> and <B>output</B> maps and the <B>type</B> of
+analysis. If the chosen <B>type</B> is `sites' then the
+<B>sites</B> list must also be specified. Any remaining parameters
+whose values have not been specified on the command line will be set
+to their default values (see <A HREF = "#parameters">below</A>). Non-interactive usage format is:
+
+<DL>
+<DD>
+<B>r.cva</B>
+[<B>-afmnoqrs</B>]
+<B>input=</B><EM>name</EM>
+[<B>target_mask=</B><EM>name</EM>]
+[<B>viewpt_mask=</B><EM>name</EM>]
+[<B>sites=</B><EM>name</EM>]
+<B>output=</B><EM>name</EM>
+[<B>obs_elev=</B><EM>value</EM>]
+[<B>target_elev=</B><EM>value</EM>]
+[<B>max_dist=</B><EM>value</EM>]
+[<B>seed=</B><EM>value</EM>]
+[<B>sample=</B><EM>value</EM>]
+<B>type=</B><EM>value</EM>
+</DL>
+
+<H2>OPTIONS</H2>
+
+<H3>Flags:</H3>
+<DL>
+
+<DT><B>-a</B>
+<DD> Treat the specified sample size as an absolute number of cells
+rather than a percentage. Note, this flag is only for use with random
+sampling.
+
+<DT><B>-f</B>
+<DD> Mark each cell in the output map with the number of viewpoints
+from which it can be seen.
+
+<DT><B>-h</B>
+<DD> Adjust height of oberserver points so that they will never fall
+below surface level. This is only useful in conjunction with the
+<b>spot</b> parameter or attribute (see <A HREF ="#atributes">below</A>).
+
+<DT><B>-i</B>
+<DD> Do not interpret control attributes in input vector map as
+described <A HREF ="#atributes">below</A>).
+
+<DT><B>-m</B> <DD> Differentiate between genuine zero data values and
+null data in the output map. When the <EM>-f</EM> flag has not been
+set this ensures that all non-viewpoint cells are coded `-1'. When
+the <EM>-f</EM> flag has been set then all cells coded `0' in the
+<EM>target_mask</EM> (if there is one) are coded `-1'. See <A HREF =
+"#notes">below</A> for guidance when to set this flag.
+
+<DT><B>-n</B> <DD> In conjunction with <B>-m</B>: write NULL values instead of '-1' to output map.
+
+<DT><B>-o</B>
+<DD> Overwrite the output raster map if it already exists.
+
+<DT><B>-q</B>
+<DD> Be quiet about progress, confirming only the output file and
+sampling frequency. This is useful when <EM>r.cva</EM> is invoked from a
+script with standard output redirected to a disk file.
+
+<DT><B>-r</B>
+<DD> Allow replacement during random sampling (i.e. allow a given cell
+to be picked as a viewpoint more than once). Note that no cumulation
+occurs for replacement viewpoints; this ensures that the resulting
+cumulative viewshed will be as though the viewpoint was selected once
+only (the final output map will, of course, be different because
+replacement effectively results in the use of fewer viewpoints for a
+given sampling frequency).
+
+<DT><B>-s</B> <DD>Be silent about progress, in which case it is
+strongly recommended that you use <EM><A
+HREF="../../html/r.info.html">r.info</A></EM> to confirm what sampling frequency
+was actually achieved.
+
+<A NAME = "parameters"><H3>Parameters:</H3></A>
+
+<DT><B>input=</B><EM>name</EM>
+<DD>Name of a raster map layer containing elevation
+data, used as program input.
+
+<DT><B>target_mask=</B><EM>name</EM> <DD>Name of a binary raster map
+layer in which target cells within the areas of interest are assigned
+the category value `1' and all other cells are assigned the category
+value `0'. A cell is treated as a target when it is viewed from a
+viewpoint (note, therefore, that viewpoints are targets when viewed
+from other viewpoints).
+
+<DT><B>viewpt_mask=</B><EM>name</EM>
+<DD>Name of a binary raster map layer in which viewpoints within the areas
+of interest are assigned the category value `1' and all other cells
+are assigned the category value `0'.
+
+<DT><B>sites=</B><EM>name</EM>
+<DD>The name of a vector map containing the sites (vector points) that are to be treated as
+viewpoints. Sites falling outside the current region are ignored.
+
+<DT><B>obs_elev=</B><EM>value</EM>
+<DD>Height of the observer (in metres) above the elevation of the
+viewpoint. Default: 1.75 m
+
+<DT><B>target_elev=</B><EM>value</EM>
+<DD>Height of the object of interest (in metres) above the elevation
+of the target cell. Default: 0.0 m
+
+<DT><B>max_dist=</B><EM>value</EM>
+<DD>Maximum distance (in metres) from the viewing point within which
+the line-of-sight analysis will be performed. Options: 0-99999
+(stated in map units). Default: 100
+
+<DT><B>seed=</B><EM>value</EM>
+<DD>The seed for the random number generator used during random
+sampling. This option allows exact reproduction of a previous result.
+Options: 0-32767. Default: 1
+
+<DT><B>sample=</B><EM>value</EM>
+<DD>The sampling frequency as a percentage of the number of map cells
+in the current region. Options: 0.0-100.0. Default: 10.0
+
+<DT><B>type=</B><EM>value</EM> <DD>The type of sampling regime used to
+select viewpoints for analysis. See <A HREF = "#description">above</A>
+for details. Options: `all', `systematic', `random', `sites'.
+
+<DT><B>spot, offseta, offsetb, azimuth1, azimuth2, vert1, vert2, radius1, radius2</B>
+<DD> See <A HREF ="#atributes">below</A>) for a description of what they do.
+
+<DT><B>curvc=</B><EM>value</EM>
+<DD>Earth curvature correction will be applied if the distance from the observer
+point exceeds the given threshold. Set this to 0.0 to disable curvature correction
+(default). While correction slows down calculations, it is recommend for
+viewing distances of more than about 3 km. The correction uses a simplified
+formula which is accurate to 1 in 10,000 for distances up to 40 km:<p>
+<i>d(height)=dist^2/(2*Re)</i><p>
+With <i>Re</i> being the radius of the earth (ca. 6,356,766 m).
+</DL>
+
+
+
+
+<A NAME = "attributes"><H3>Vector map attributes:</H3>
+
+When using 'sites' type operation mode, <EM>r.cva</EM> checks for the
+presence of a number of attributes in the input vector map. The names
+of these attributes are all uppercase and match those used by the
+viewshed module of ESRI's ArcGIS Spatial Extension (trademarks of ESRI Corp.).
+They control <EM>r.cva</EM> operation in an analogous way.
+Each attribute, if present, must store a double value. The following
+attribute names and value ranges are recognised and interpreted:
+
+<dl>
+<dt>SPOT</dt>
+ <dd>This determines absolute observer height. If this attribute is given,
+ observer height will no longer be determined from the DEM.
+ <EM>r.cva</EM> can automatically raise this value to that of the DEM
+ point below the oberserver's position, if the value in SPOT would
+ cause the observer to get below ground level (use the <b>-h</b> flag).</dd>
+<dt>OFFSETA</dt>
+ <dd>Use this to add an offset to each oberver point's height. This is
+ just an alias for the <b>obs_elev</b> option and exists for
+ compatibility reasons.</dd>
+<dt>OFFSETB</dt>
+ <dd>The value stored in 'OFFSETB' will be added to each target point
+ at the moment it is considered for visibility. All blocking points
+ on the line-of-sight between observer and target will be left untouched.
+ This is different from the <b>target_elev</b> parameter. It can be
+ used to simulate visibility of targets with a certain height.</dd>
+<dt>AZIMUTH1, AZIMUTH2</dt>
+ <dd>These attributes must be given as values between 0.0 and 360.0 and
+ can be used to limit the horizontal scan range of the viewshed.
+ Limits are measured in degrees with 0.0 in the North and values
+ increasing clock-wise. Thus, if you wanted to limit visibility scans
+ to eastward directions only, you would set AZIMUTH1=1 and AZIMUTH2=180.
+ </dd>
+<dt>VERT1, VERT2</dt>
+ <dd>In analogy to 'AZIMUTH1' and 'AZIMUTH2', these attributes can be
+ used to limit the vertical range of the scan. The range is -90 to 0 for
+ 'VERT1' and 0 to 90 for 'VERT2'.
+ At -90 degrees, the observer looks straight down to the ground and at
+ 90 degress up into the sky.
+ </dd>
+<dt>RADIUS1, RADIUS2</dt>
+ <dd>Finally, two attributes exist to limit the distance of the scan.
+ 'RADIUS1' determines the minimum distance from the observer and 'RADIUS2'
+ the maximum distance. The latter is thus equivalent to using the
+ <b>max_dist</b> parameter and exists for compatibility reasons.
+ </dd>
+</dl>
+
+Interpretation of these attributes can be turned off by using the <b>-i</b> flag.
+<p>
+Note that you can use the same control values in operation modes
+other than 'sites', if you specify the equivalent parameters
+(see <A HREF = "#parameters">above</A>). The latter will also work
+as global settings in 'sites' mode. Attributes stored in the vector
+map table will override these global settings unless the <b>-i</b>
+flag is specified. The meaning of <b>OFFSETA, OFFSETB</b> and <b>RADIUS1</b> is
+the same as <b>obs_elev, target_elev</b> and <b>max_dist</b>. The latter still
+exist for compatibility with earlier versions of <em>r.cva</em>.
+
+
+<A NAME = "notes"><H2>NOTES</H2></A>
+
+<EM>r.cva</EM> works in the current geographic region, but ignores the
+current mask. Masking is achieved by specifying a target_mask and/or a
+viewpt_mask. For accurate results, the program must be run with the
+resolution of the geographic region set equal to the resolution of the
+data (see <EM><A HREF="../../html/g.region.html">g.region</A></EM>).
+<P>
+This program assumes all height and distance measurements in meters.
+Legacy units must be converted to meters by the user before running <EM>r.cva</EM>.
+<P>
+It is also assumed, that N-S and E-W resolution of the current region are the
+same. If they differ, distances calculated for limiting the viewshed will
+be imprecise and results affected accordingly.
+<P>
+The number of cells visible from a given viewpoint includes the
+viewpoint itself only if the viewpoint is a target cell.
+Consequently, if neither the viewpoint nor any other cell within
+<EM>max_dist</EM> of it is a target cell then the size of that
+viewpoint's viewshed will be zero. Set the <EM>-m</EM> flag to locate
+such viewpoints (they will show as cells coded `0' as opposed to
+`-1'). Similarly, set the <EM>-m</EM> flag when the <EM>-f</EM> flag has
+been set in order to differentiate cells which are genuinely not
+visible from any viewpoint as opposed to those that are not target
+cells.
+<P>
+<EM>r.cva</EM> always informs the user what sampling frequency was
+attempted and subsequently achieved. This information is printed to
+standard output and is also placed in the output map's history file
+(which can be viewed using <EM><A HREF="../../html/r.info.html">r.info</A></EM>).
+The sample attempted may differ slightly from that requested owing to
+the need to sample a whole number of cells. Any such discrepancy will
+be greatest in the case of small maps of coarse resolution. The
+sample achieved may be further reduced when a <EM>viewpt_mask</EM> has been
+specified. This will not occur with the option `random', because its
+algorithm repeatedly allocates viewpoints until the number that fall
+within mask cells coded `1' produces the frequency attempted. Such a
+reduction can, however, occur with the options `all', `systematic' and
+`sites', because their algorithms allocate viewpoints without
+reference to the mask and then simply ignore those that fall in mask
+cells coded `0'.
+<P>
+
+<H2>SEE ALSO</H2>
+
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>
+<br>
+<EM><A HREF="../../html/r.info.html">r.info</A></EM>
+<br>
+<EM><A HREF="../../html/g.region.html">g.region</A></EM>
+
+<P>
+
+Cumulative viewshed analysis was defined by D. Wheatley, 1995,
+`Cumulative viewshed analysis: a GIS-based method for investigating
+intervisibility, and its archaeological application', in G. R. Lock
+and Z. Stancic (eds.) <EM>Archaeology and Geographic Information
+Systems: A European Perspective</EM>, London: Taylor and Francis,
+pp. 171-186.
+<P>
+Lake et al. (M. W. Lake, P. E. Woodman and S. J. Mithen, 1998,
+`Tailoring GIS software for archaeological applications: an example
+concerning viewshed analysis', <EM>Journal of Archaeological
+Science</EM> 25: 27-38) discuss the first version of <EM>r.cva</EM>
+and provide some timing estimates. They also discuss the importance
+of eliminating the edge effect that occurs when calculating the number
+of cells visible from a viewpoint that is closer than
+<EM>max_dist</EM> to the edge of the map.
+<P>
+Fisher et al. (P. Fisher, C. Farrelly, A. Maddocks and C. Ruggles,
+1997, `Spatial Analysis of Visible Areas from the Bronze Age Cairns of
+Mull', <EM>Journal of Archaeological Science</EM> 24: 581--592)
+provide a useful review of some of the pitfalls that may trap unwary
+users of cumulative viewshed analysis. Most of these are easily
+avoided using <EM>r.cva</EM>.
+
+<H2>AUTHOR</H2>
+
+Mark Lake, Institute of Archaeology, University College London (the author).
+<P>
+With contributions by:
+Benjamin Ducke, Institute of Archaeology, University of Bamberg, Germany.
+
+<H2>ACKNOWLEDGEMENTS</H2>
+
+This version of <EM>r.cva</EM> was completed during the author's
+tenure of a Leverhulme Trust Special Research Fellowship at the
+Institute of Archaeology, University College London, U.K.
+<P>
+The first version of <EM>r.cva</EM> was written by the author as part of the
+MAGICAL Project directed by Dr. Steven Mithen (Dept. of Archaeology,
+University of Reading, U.K.). The MAGICAL Project was made possible
+by a Natural Environment Research Council award (NERC GR3/9540) to
+Dr. Mithen.
+<P>
+<EM>r.cva</EM> draws heavily on the code for <EM><A
+HREF="../../html/r.los.html">r.los</A></EM>, which was written by Kewan
+Q. Khawaja, Intelligent Engineering Systems Laboratory, M.I.T.
+<P>
+Adaptation to GRASS 6 vector and raster formats, additional functionality,
+documentation updates by Benjamin Ducke (2005).
+
+<H2>BUGS</H2>
+
+Missing values in site attributes are currently not handled.
+<p>
+
+Earth curvature correction has not been sufficiently tested. The 'offsetb' option (=target_elev)
+has not been sufficiently tested. Please feel free to test and report!
+<p>
+
+When <EM>max_dist</EM> is set to 100, the N-S and E-W resolutions are
+50, and a <EM>target_mask</EM> has been specified, the cells
+immediately NE, SE, SW and NW of the viewpoint may not be masked when
+then should be. Please report any others that you find to the author.
+<P>
+
+Note that <EM>r.cva</EM> fixes two bugs in the current version of
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>: that <EM>obs_elev</EM> is
+truncated to its integer component and that the target mask
+(<EM>patt_map</EM> in <EM><A HREF="../../html/r.los.html">r.los</A></EM>) does
+not apply to the viewpoint or its eight immediate neighbours.
Added: grass-addons/grass6/raster/r.cva/global_vars.h
===================================================================
--- grass-addons/grass6/raster/r.cva/global_vars.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/global_vars.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,46 @@
+/***********************************************************************/
+/*
+ global_vars.h
+
+ CONTAINS
+
+ Global variables */
+
+/***********************************************************************/
+
+#ifdef MAIN
+#define GLOBAL
+#else
+#define GLOBAL extern
+#endif
+
+GLOBAL double Re;
+
+GLOBAL double obs_elev, target_elev;
+GLOBAL double max_dist;
+GLOBAL double sample;
+GLOBAL int background, seed, from_cells;
+GLOBAL int patt_flag, patt_flag_v;
+GLOBAL int show_mask;
+GLOBAL int fract;
+GLOBAL struct Cell_head window; /* Region info used for sites option */
+
+/* distance above which to apply earth curvature correction */
+GLOBAL double DIST_CC;
+
+/* these are set if the vector points support the ArcGIS attributes */
+GLOBAL double SPOT;
+GLOBAL int SPOT_SET;
+GLOBAL double OFFSETB;
+GLOBAL int OFFSETB_SET;
+GLOBAL int ADJUST;
+GLOBAL double AZIMUTH1;
+GLOBAL int AZIMUTH1_SET;
+GLOBAL double AZIMUTH2;
+GLOBAL int AZIMUTH2_SET;
+GLOBAL double VERT1;
+GLOBAL int VERT1_SET;
+GLOBAL double VERT2;
+GLOBAL int VERT2_SET;
+GLOBAL double RADIUS1;
+
Added: grass-addons/grass6/raster/r.cva/init_segment_file.c
===================================================================
--- grass-addons/grass6/raster/r.cva/init_segment_file.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/init_segment_file.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,74 @@
+/***********************************************************************/
+/*
+ init_segment_file.c
+
+ CONTAINS
+
+ 1) init_segment_file */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ init_segment_file
+
+ Called from:
+
+ 1) scan_all_cells scan_all_cells.c
+ 2) systematic_sample systematic_sample.c
+ 3) random_sample random_sample.c
+ 4) point_sample point_sample.c
+
+ Calls:
+
+ None */
+
+/***********************************************************************/
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "global_vars.h"
+
+void init_segment_file (int nrows, int ncols, SEGMENT *seg_to_init_p,
+ SEGMENT *seg_patt_p)
+{
+ int row, col;
+ char *value;
+ CELL data;
+
+ value = (char*) &data;
+
+ if (from_cells && patt_flag && (background == -1))
+ {
+ CELL data_read;
+ char *value_to_get;
+
+ value_to_get = (char*) &data_read;
+ for (row = 0; row < nrows; row ++)
+ {
+ for (col = 0; col < ncols; col ++)
+ {
+ segment_get (seg_patt_p, value_to_get,
+ row, col);
+ if (! data_read)
+ data = background;
+ else
+ data = 0;
+ segment_put (seg_to_init_p, value,
+ row, col);
+ }
+ }
+ }
+ else
+ {
+ data = background;
+ for (row = 0; row < nrows; row ++)
+ {
+ for (col = 0; col < ncols; col ++)
+ {
+ segment_put (seg_to_init_p, value, row, col);
+ }
+ }
+ }
+}
Added: grass-addons/grass6/raster/r.cva/init_segment_file.h
===================================================================
--- grass-addons/grass6/raster/r.cva/init_segment_file.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/init_segment_file.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,10 @@
+#ifndef __INIT_SEGMENT_FILE_H__
+#define __INIT_SEGMENT_FILE_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void init_segment_file (int nrows, int ncols, SEGMENT *seg_to_init_p,
+ SEGMENT *seg_patt_p);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/line_of_sight.c
===================================================================
--- grass-addons/grass6/raster/r.cva/line_of_sight.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/line_of_sight.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,104 @@
+/***********************************************************************/
+/*
+ line_of_site.c
+
+ Updated by Mark Lake on 17/8/01 for GRASS 5.x
+
+ CONTAINS
+
+ 1) line_of_site */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ line_of_site
+
+ Called from:
+
+ 1) scan_all_cells scan_all_cells.c
+ 2) random_sample random_sample.c
+ 3) systematic_sample systematic_sample.c
+ 4) point_sample point_sample.c
+
+ Calls:
+
+ 1) segment() segment3.c */
+
+/***********************************************************************/
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+#include "segment3.h"
+
+void line_of_sight (CELL cell_no, int row_viewpt, int col_viewpt,
+ int nrows, int ncols, double viewpt_elev,
+ struct point *heads[],
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p, SEGMENT *seg_patt_p,
+ int patt_flag, RASTER_MAP_TYPE data_type)
+
+ /* struct point *heads[] is same as "struct point **heads" */
+
+{
+ int a,b;
+ int segment_no, flip, xmax, ymax, sign_on_y, sign_on_x;
+ double slope_1, slope_2;
+
+ /* Do los analysis from current viewpoint for 16 segments */
+ for(segment_no=1; segment_no<=16; segment_no++)
+ {
+ sign_on_y= 1- (segment_no-1)/8 * 2;
+ if(segment_no>4 && segment_no<13)
+ sign_on_x= -1;
+ else sign_on_x=1;
+
+
+ /* Calc slopes for bounding rays of this segment */
+
+ if(segment_no==1 || segment_no==4 || segment_no==5 ||
+ segment_no==8 || segment_no==9 || segment_no==12 ||
+ segment_no==13 || segment_no==16)
+ {
+ slope_1= 0.0;
+ slope_2= 0.5;
+ }
+ else
+ {
+ slope_1= 0.5;
+ slope_2= 1.0;
+ }
+ if(segment_no==1 || segment_no==2 || segment_no==7 ||
+ segment_no==8 || segment_no==9 || segment_no==10 ||
+ segment_no==15 || segment_no==16)
+ flip = 0;
+ else flip = 1;
+
+
+ /* Calculate max and min 'x' and 'y' */
+
+ a= ((ncols-1)*(sign_on_x+1)/2 - sign_on_x * col_viewpt);
+ b= (1-sign_on_y)/2*(nrows-1) + sign_on_y*row_viewpt;
+ if(flip==0)
+ {
+ xmax=a;
+ ymax=b;
+ }
+ else
+ {
+ xmax=b;
+ ymax=a;
+ }
+
+ /* Perform los analysis for this segment */
+ /* this function is in file 'segment3.c' */
+ heads[segment_no-1] = segment(segment_no,xmax,ymax,
+ slope_1,slope_2,flip,sign_on_y,sign_on_x,
+ viewpt_elev,
+ seg_in_p,seg_out_1_p, seg_patt_p,
+ row_viewpt,col_viewpt,patt_flag,
+ cell_no, data_type);
+ }
+}
Added: grass-addons/grass6/raster/r.cva/line_of_sight.h
===================================================================
--- grass-addons/grass6/raster/r.cva/line_of_sight.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/line_of_sight.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,13 @@
+#ifndef __LINE_OF_SIGHT_H__
+#define __LINE_OF_SIGHT_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void line_of_sight (CELL cell_no, int row_viewpt, int col_viewpt,
+ int nrows, int ncols, double viewpt_elev,
+ struct point *heads[],
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p, SEGMENT *seg_patt_p,
+ int patt_flag, RASTER_MAP_TYPE data_type);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/main.c
===================================================================
--- grass-addons/grass6/raster/r.cva/main.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/main.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,687 @@
+/***********************************************************************/
+/*
+ r.cva in $GIS/src.contrib/IoA
+
+ *****
+
+ Mark Lake 30/6/99 (pdated 17/8/01 for GRASS 5.x)
+
+ University College London
+ Institute of Archaeology
+ 31-34 Gordon Square
+ London. WC1H 0PY
+
+ Email: mark.lake at ucl.ac.uk
+
+ *****
+
+ Adaptations for new GRASS 5 site and raster formats
+ Adaptations for new GRASS 6 vector format
+
+ by Benjamin Ducke 2004
+ benducke at compuserve.de
+
+ Summary of changes:
+ - updated to GRASS 5 sites API
+ - use GRASS routines to convert Site coordinates to raster coords
+ - handle CELL, FCELL or DCELL as elevation input map format
+ - handle NULL values in input map by copying them through to result
+ - renamed flag '-n' to '-a'
+ - new flag '-n': convert -1 to NULL in output map
+ - max viewing distance now set to 2000 m by default
+ - handle all progress display through 'G_percent()'
+ - all output now goes to 'stdout', so it can be piped to a file
+ - made source code more readable
+ - updated manual pages
+ - update tcltkgrass menu file
+
+ - updated to GRASS 6
+ - use new vector points format instead of sites
+ - read and interpret special attributes in the vector map's DB table
+ the same way ArcInfo(tm) Spatial Analyst Extension(tm I guess) does
+ - same attributes can also be set as global options
+ - earth curvature correction
+ - fixed nasty FP handling bugs that would cause the previous GRASS 5 version to
+ crash when compiled with newer versions of GCC
+ - lots of fixes for NULL data handling in DEM
+ - updated HTML manual page to reflect changes
+
+ *****
+
+ ACKNOWLEDGEMENTS
+
+ The first version of r.cva was written by the author as part of the
+ MAGICAL Project directed by Dr. Steven Mithen (University of Reading,
+ U.K.). The MAGICAL Project was made possible by a Natural Environment
+ Research Council award (NERC GR3/9540) to Dr. Mithen.
+
+ This version of r.cva was significantly enhanced during the author's
+ tenure of a Leverhulme Trust Special Research Fellowship at the
+ Institute of Archaeology, U.C.L., U.K.
+
+ r.cva draws heavily on the code for r.los, which was written by
+ Kewan Q. Khawaja, Intelligent Engineering Systems Laboratory, M.I.T.
+
+ *****
+
+ PURPOSE
+
+ 1) Perform cumulative viewshed (CVA) analysis.
+
+ OPTIONS
+
+ 1) Record the number of cells visible from:
+
+ all every map cell;
+ systematic a systematic sample of map cells;
+ random a random sample of map cells;
+ sites a list of sites.
+
+ 2) Record how often each cell in the map is visible from:
+
+ all every map cell;
+ systematic a systematic sample of map cells;
+ random a random sample of map cells;
+ sites a list of sites.
+
+ 3) Random sampling may occur with or without replacement.
+
+ 4) All types of sampling may be subject to binary masks such that:
+
+ i) Not all possible destination cells (cells to be `looked at')
+ are of interest;
+
+ ii) Not all possible viewpoints are of interest.
+
+ 5) It is possible to specify the height of the observer.
+
+ 6) It is possible to specify the height of the target object
+ above ground.
+
+ NOTE
+
+ r.cva includes fixes for two bugs in the current version of r.los.
+ These are that the observer elevation is cast to an integer, and that
+ the patt_map does not effect the viewpoint and its immediately
+ neighbouring cells */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/* main
+
+ Called from:
+
+ 1) None
+
+ Calls functions in:
+
+ 1) scan_all_cells.c
+ 2) random_sample.c
+ 3) systematic_sample.c
+ 4) point_sample.c
+ 5) segment_file.c */
+
+/***********************************************************************/
+
+#include <stdlib.h>
+#include <string.h>
+#include "config.h"
+#include <grass/segment.h>
+#include <grass/gis.h>
+#define MAIN
+#include "global_vars.h"
+#include "point.h"
+
+#include "segment_file.h"
+#include "scan_all_cells.h"
+#include "systematic_sample.h"
+#include "random_sample.h"
+#include "point_sample.h"
+
+#define kSEGFRACT 4 /* Controls no. of segments
+ DO NOT CHANGE */
+
+
+int main(int argc, char *argv[])
+{
+ int nrows, ncols;
+ int new,old,patt,patt_v,in_fd,out_fd_1, out_fd_2;
+ int patt_fd, patt_fd_v, rnd_fd;
+ int same_mask;
+ char old_mapset [50], patt_mapset [50], patt_mapset_v [50];
+ char current_mapset [50];
+ char in_name [128], out_name_1 [128], out_name_2 [128];
+ char patt_name [128], patt_name_v [128], rnd_name [128];
+ double attempted_sample, actual_sample;
+ struct Categories cats;
+ struct History history;
+ extern struct Cell_head window;
+ CELL *cell = NULL;
+ DCELL *dcell = NULL;
+ FCELL *fcell = NULL;
+ SEGMENT seg_in, seg_out_1, seg_out_2, seg_patt, seg_patt_v, seg_rnd;
+ SEGMENT *seg_patt_p_v = NULL;
+ struct point *segment();
+ struct GModule *module;
+ struct Option *input, *patt_map, *patt_map_v, *type;
+ struct Option *opt_obs_elev, *opt_target_elev;
+ struct Option *opt_max_dist, *output, *opt_seed, *opt_sample;
+ struct Option *sites;
+ struct Option *curvature_corr;
+ struct Option *spot, *offseta, *offsetb, *azimuth1, *azimuth2, *vert1, *vert2, *radius1, *radius2;
+ struct Flag *replace, *from, *overwrite, *quiet, *mask, *silent;
+ struct Flag *adjust, *absolute, *nulls, *ignore;
+ /* stores cell data mode (CELL_TYPE, FCELL_TYPE or DCELL_TYPE) */
+ int cell_mode;
+ int make_nulls = 0; /* this signals whether -1 in output map should be
+ converted to NULL */
+
+
+ /***********************************************************************/
+ /* Preliminary things */
+ /***********************************************************************/
+
+ /* Initialize the GIS calls and sort out error handling */
+
+ G_gisinit (argv[0]);
+ G_sleep_on_error (0);
+
+
+ /* Define the options */
+
+ module = G_define_module();
+ module->description = "Multiple/cumulative viewshed analysis";
+
+ input = G_define_standard_option(G_OPT_R_INPUT) ;
+ input->key = "input";
+ input->type = TYPE_STRING;
+ input->required = YES;
+ input->gisprompt = "old,fcell,dcell,cell,raster" ;
+ input->description= "Elevation raster map" ;
+
+ patt_map = G_define_standard_option(G_OPT_R_INPUT);
+ patt_map->key = "target_mask";
+ patt_map->type = TYPE_STRING;
+ patt_map->required = NO;
+ patt_map->description= "Binary raster map for target cells";
+ patt_map->gisprompt = "old,cell,raster" ;
+
+ patt_map_v = G_define_standard_option(G_OPT_R_INPUT);
+ patt_map_v->key = "viewpt_mask";
+ patt_map_v->type = TYPE_STRING;
+ patt_map_v->required = NO;
+ patt_map_v->description= "Binary raster map for viewpoints";
+ patt_map_v->gisprompt = "old,cell,raster" ;
+
+ sites = G_define_standard_option(G_OPT_V_MAP);
+ sites->key = "sites";
+ sites->type = TYPE_STRING;
+ sites->required = NO;
+ sites->description= "Vector map with sites points" ;
+
+ output = G_define_standard_option(G_OPT_R_OUTPUT);
+ output->key = "output";
+ output->type = TYPE_STRING;
+ output->required = YES;
+ output->gisprompt = "cell,raster" ;
+ output->description= "Raster map name for storing results";
+
+ opt_obs_elev = G_define_option() ;
+ opt_obs_elev->key = "obs_elev";
+ opt_obs_elev->type = TYPE_DOUBLE;
+ opt_obs_elev->required = NO;
+ opt_obs_elev->answer = "1.75";
+ opt_obs_elev->description= "Height of the viewing location";
+
+ opt_target_elev = G_define_option() ;
+ opt_target_elev->key = "target_elev";
+ opt_target_elev->type = TYPE_DOUBLE;
+ opt_target_elev->required = NO;
+ opt_target_elev->answer = "0.0";
+ opt_target_elev->description= "Height of the target location";
+
+ opt_max_dist = G_define_option() ;
+ opt_max_dist->key = "max_dist";
+ opt_max_dist->type = TYPE_DOUBLE;
+ opt_max_dist->required = NO;
+ opt_max_dist->answer = "2000";
+ opt_max_dist->options = "0.0-9999999.9" ;
+ opt_max_dist->description= "Max viewing distance in meters" ;
+
+ opt_seed = G_define_option() ;
+ opt_seed->key = "seed";
+ opt_seed->type = TYPE_INTEGER;
+ opt_seed->required = NO;
+ opt_seed->answer = "1";
+ opt_seed->options = "0-32767" ;
+ opt_seed->description= "Seed for the random number generator" ;
+
+ opt_sample = G_define_option() ;
+ opt_sample->key = "sample";
+ opt_sample->type = TYPE_DOUBLE;
+ opt_sample->required = NO;
+ opt_sample->answer = "10.0";
+ opt_sample->description= "Sampling frequency as a percentage of cells" ;
+
+ type = G_define_option();
+ type->key = "type";
+ type->type = TYPE_STRING;
+ type->required = YES;
+ type->answer = "sites";
+ type->options = "sites,systematic,random,all";
+ type->description= "Type of vieshed analysis";
+
+ spot = G_define_option();
+ spot->key = "spot";
+ spot->type = TYPE_DOUBLE;
+ spot->required = NO;
+ spot->description = "Absolute observer height";
+
+ offseta = G_define_option();
+ offseta->key = "offseta";
+ offseta->type = TYPE_DOUBLE;
+ offseta->required = NO;
+ offseta->description = "Observer offset (OFFSETA=obs_elev)";
+
+ offsetb = G_define_option();
+ offsetb->key = "offsetb";
+ offsetb->type = TYPE_DOUBLE;
+ offsetb->required = NO;
+ offsetb->description = "target offset (OFFSETB=target_elev)";
+
+ azimuth1 = G_define_option();
+ azimuth1->key = "azimuth1";
+ azimuth1->type = TYPE_DOUBLE;
+ azimuth1->required = NO;
+ azimuth1->options = "0.0-360.0";
+ azimuth1->description = "Minimum azimuth (AZIMUTH1)";
+
+ azimuth2 = G_define_option();
+ azimuth2->key = "azimuth2";
+ azimuth2->type = TYPE_DOUBLE;
+ azimuth2->required = NO;
+ azimuth2->options = "0.0-360.0";
+ azimuth2->description = "Maximum azimuth (AZIMUTH2)";
+
+ vert1 = G_define_option();
+ vert1->key = "vert1";
+ vert1->type = TYPE_DOUBLE;
+ vert1->required = NO;
+ vert1->options = "0.0-90.0";
+ vert1->description = "Minimum vertical angle (VERT1)";
+
+ vert2 = G_define_option();
+ vert2->key = "vert2";
+ vert2->type = TYPE_DOUBLE;
+ vert2->required = NO;
+ vert2->options = "-90.0-0.0";
+ vert2->description = "Maximum vertical angle (VERT2)";
+
+ radius1 = G_define_option();
+ radius1->key = "radius1";
+ radius1->type = TYPE_DOUBLE;
+ radius1->required = NO;
+ radius1->description = "Minimum distance from observer (RADIUS1)";
+
+ radius2 = G_define_option();
+ radius2->key = "radius2";
+ radius2->type = TYPE_DOUBLE;
+ radius2->required = NO;
+ radius2->description = "Maximum distance to observer (RADIUS2)";
+
+ curvature_corr = G_define_option ();
+ curvature_corr->key = "curvc";
+ curvature_corr->type = TYPE_DOUBLE;
+ curvature_corr->required = NO;
+ curvature_corr->answer = "0.0";
+ curvature_corr->description = "Earth curvature correction threshold (0.0 = off)";
+
+ from = G_define_flag ();
+ from->key = 'f';
+ from->description = "Calculate the `visibility from' rather than `viewsheds of' ";
+
+ mask = G_define_flag ();
+ mask->key = 'm';
+ mask->description = "Differentiate masked cells from data value zero ";
+
+ absolute = G_define_flag ();
+ absolute->key = 'a';
+ absolute->description = "Treat sample size as no. of cells ";
+
+ overwrite = G_define_flag ();
+ overwrite->key = 'o';
+ overwrite->description = "Overwrite output file if it exists ";
+
+ quiet = G_define_flag ();
+ quiet->key = 'q';
+ quiet->description = "Run quietly (i.e. do not report \% complete) ";
+
+ replace = G_define_flag ();
+ replace->key = 'r';
+ replace->description = "Allow replacement during random sampling ";
+
+ adjust = G_define_flag ();
+ adjust->key = 'h';
+ adjust->description = "Adjust spot heights below surface";
+
+ ignore = G_define_flag ();
+ ignore->key = 'i';
+ ignore->description = "Ignore site attributes";
+
+ nulls = G_define_flag ();
+ nulls->key = 'n';
+ nulls->description = "Convert -1 to NULL in output map ";
+
+ silent = G_define_flag ();
+ silent->key = 's';
+ silent->description = "Run silently ";
+
+ if (G_parser(argc, argv))
+ exit (-1);
+
+
+ /* Make numeric parameters globally available */
+ sscanf (opt_obs_elev->answer, "%lf", &obs_elev);
+ //sscanf (opt_target_elev->answer, "%lf", &target_elev);
+ if ( opt_target_elev->answer != NULL ) {
+ OFFSETB_SET = 1; /* OFFSETB is the same as old parameter target_elev */
+ OFFSETB = atof (opt_target_elev->answer);
+ }
+
+ sscanf (opt_max_dist->answer, "%lf", &max_dist);
+ sscanf (opt_seed->answer, "%d", &seed);
+ sscanf (opt_sample->answer, "%lf", &sample);
+ from_cells = from->answer;
+ if(patt_map->answer == NULL)
+ patt_flag = 0;
+ else
+ patt_flag = 1;
+ if(patt_map_v->answer == NULL)
+ patt_flag_v = 0;
+ else
+ patt_flag_v = 1;
+ if (mask->answer)
+ background = -1;
+ else
+ background = 0;
+
+ if (adjust->answer) {
+ ADJUST = 1;
+ } else {
+ ADJUST = 0;
+ }
+
+ if (nulls->answer)
+ make_nulls = 1;
+
+ fract = kSEGFRACT;
+
+ /* Check sample size is sensible */
+ if (!absolute->answer)
+ {
+ if ((sample < 0.0) || (sample > 100.0))
+ G_fatal_error ("Sample size must be in range 0 to 100 ");
+ }
+ else
+ if (sample < 0.0)
+ G_fatal_error ("Sample size must be a positive number ");
+
+ /* If both destination and viewpoint masks have been specified
+ check whether they are the same. If so set flag so that only
+ one is loaded into memory and the other is simply made a duplicate
+ set of pointers. */
+ if (patt_flag && patt_flag_v)
+ {
+ if (strcmp (patt_map->answer, patt_map_v->answer) == 0)
+ same_mask = 1;
+ else
+ same_mask = 0;
+ }
+ else same_mask = 0;
+
+ /* Must be quiet if silent */
+ if (silent->answer)
+ quiet->answer = 1;
+
+ /* Get current region */
+ G_get_window (&window);
+ nrows = G_window_rows();
+ ncols = G_window_cols();
+
+ /* determine mode of input map */
+ cell_mode = G_raster_map_type (input->answer,"");
+
+ /* Allocate buffer space for row-io to layer */
+ cell = G_allocate_c_raster_buf(); /* we need a cell buffer, anyhow */
+ if (cell_mode == FCELL_TYPE){ /* we might also need an fcell or dcell buffer for */
+ fcell = G_allocate_f_raster_buf(); /* elevation data */
+ }
+ if (cell_mode == DCELL_TYPE) {
+ dcell = G_allocate_d_raster_buf();
+ }
+
+ /* Open files and initialise memory for random access */
+ /* create segment files for all input raster maps and store */
+ /* raster maps in them */
+ /* input->answer is the elevation map */
+
+ /* store elevation map in a CELL, FCELL or DCELL segment file, */
+ /* according to its type */
+ if ( cell_mode == CELL_TYPE) {
+ Segment_infile (input->answer, old_mapset, "", &old,
+ in_name, &in_fd, &seg_in, cell, fract, cell_mode);
+ }
+ if ( cell_mode == FCELL_TYPE) {
+ Segment_infile (input->answer, old_mapset, "", &old,
+ in_name, &in_fd, &seg_in, fcell, fract, cell_mode);
+ }
+ if ( cell_mode == DCELL_TYPE) {
+ Segment_infile (input->answer, old_mapset, "", &old,
+ in_name, &in_fd, &seg_in, dcell, fract, cell_mode);
+ }
+
+ /* store mask maps in CELL segment files */
+ if (patt_flag)
+ Segment_infile (patt_map->answer, patt_mapset, "", &patt,
+ patt_name, &patt_fd, &seg_patt, cell, fract, CELL_TYPE);
+
+ if (patt_flag_v) {
+ if (same_mask)
+ seg_patt_p_v = &seg_patt;
+ else {
+ seg_patt_p_v = &seg_patt_v;
+ Segment_infile (patt_map_v->answer, patt_mapset_v, "", &patt_v,
+ patt_name_v, &patt_fd_v, seg_patt_p_v, cell, fract, CELL_TYPE);
+ }
+ }
+
+
+ Segment_named_outfile (output->answer, current_mapset, &new,
+ out_name_2, &out_fd_2, &seg_out_2,
+ overwrite->answer, quiet->answer, fract, CELL_TYPE);
+
+ Segment_tmpfile (out_name_1, &out_fd_1, &seg_out_1, fract, CELL_TYPE);
+
+ /* Record name of output map if in quiet mode but not silent mode */
+ if (!silent->answer)
+ {
+ if (quiet->answer)
+ {
+ fprintf (stdout, "\n\nOutput map is: '%s'", output->answer);
+ if (mask->answer)
+ fprintf (stdout, " (masked cells = -1)");
+ fflush (stdout);
+ }
+ }
+
+
+ /***********************************************************************/
+ /* Perform analysis of the type requested */
+ /***********************************************************************/
+
+ /* earth curvature correction threshold */
+ DIST_CC = atof (curvature_corr->answer);
+ Re = 6356766.0; /* radius of Earth in m */
+
+ /* set attributes to default values */
+ SPOT = 0.0;
+ SPOT_SET = 0;
+ OFFSETB = 0.0;
+ OFFSETB_SET = 0;
+ AZIMUTH1 = 0.0;
+ AZIMUTH1_SET = 0;
+ AZIMUTH2 = 360.0;
+ AZIMUTH2_SET = 0;
+ VERT1 = 90.0;
+ VERT1_SET = 0;
+ VERT2 = -90.0;
+ VERT2_SET = 0;
+ RADIUS1 = 0.0;
+
+ /* check for global attributes */
+ if (spot->answer != NULL) {
+ SPOT = atof (spot->answer);
+ SPOT_SET = 1;
+ }
+ if (offseta->answer != NULL) {
+ obs_elev = atof (offseta->answer);
+ }
+ if (offsetb->answer != NULL) {
+ OFFSETB_SET = 1;
+ OFFSETB = atof (offsetb->answer);
+ }
+ if (azimuth1->answer != NULL) {
+ AZIMUTH1 = atof (azimuth1->answer);
+ AZIMUTH1_SET = 1;
+ }
+ if (azimuth2->answer != NULL) {
+ AZIMUTH2 = atof (azimuth2->answer);
+ AZIMUTH2_SET = 1;
+ }
+ if (vert1->answer != NULL) {
+ VERT1 = atof (vert1->answer);
+ VERT1_SET = 1;
+ }
+ if (vert2->answer != NULL) {
+ VERT2 = atof (vert2->answer);
+ VERT2_SET = 1;
+ }
+ if (radius2->answer != NULL) {
+ max_dist = atof (radius2->answer);
+ }
+ if ((radius1->answer != NULL) && (atof (radius1->answer) < max_dist)) {
+ RADIUS1 = atof (radius1->answer);
+ }
+
+
+ /* use all raster cells in input map as observer locations (!) */
+ if (strcmp (type->answer, "all") == 0)
+ scan_all_cells (nrows, ncols,
+ &seg_in, &seg_out_1, &seg_out_2, &seg_patt,
+ seg_patt_p_v, &attempted_sample, &actual_sample,
+ quiet->answer, cell_mode);
+
+ /* use a gridded set of raster cells as observer locations */
+ if (strcmp (type->answer, "systematic") == 0)
+ systematic_sample (nrows, ncols,
+ &seg_in, &seg_out_1, &seg_out_2, &seg_patt,
+ seg_patt_p_v, &attempted_sample,
+ &actual_sample, quiet->answer,
+ absolute->answer, cell_mode);
+
+ /* use a random set of raster cells as observer locations */
+ if (strcmp (type->answer, "random") == 0) {
+ Segment_tmpfile (rnd_name, &rnd_fd, &seg_rnd, fract, CELL_TYPE);
+ random_sample (nrows, ncols,
+ &seg_in, &seg_out_1, &seg_out_2, &seg_patt,
+ seg_patt_p_v, &seg_rnd,
+ &attempted_sample, &actual_sample,
+ replace->answer, quiet->answer,
+ absolute->answer, cell_mode);
+ Close_segmented_tmpfile (rnd_name, rnd_fd, &seg_rnd);
+ }
+
+ /* do a point sample using site locations */
+ if (strcmp (type->answer, "sites") == 0)
+ point_sample (&window, nrows, ncols,
+ &seg_in, &seg_out_1, &seg_out_2, &seg_patt,
+ seg_patt_p_v, &attempted_sample, &actual_sample,
+ sites->answer, quiet->answer, cell_mode, ignore->answer);
+
+
+ /***********************************************************************/
+ /* Tidy up */
+ /***********************************************************************/
+
+ /* write result file */
+ Close_segmented_outfile (output->answer, new, out_name_2, out_fd_2,
+ &seg_out_2, cell, quiet->answer,
+ &seg_in, cell_mode, make_nulls);
+
+ /* Close remaining files */
+ Close_segmented_infile (old, in_name, in_fd, &seg_in);
+ if (patt_flag)
+ Close_segmented_infile (patt, patt_name, patt_fd, &seg_patt);
+ if ((patt_flag_v) && (! same_mask))
+ Close_segmented_infile (patt_v, patt_name_v, patt_fd_v,
+ seg_patt_p_v);
+ Close_segmented_tmpfile (out_name_1, out_fd_1, &seg_out_1);
+
+
+ /* Remind user about value used to indicate NULL data */
+ if ((! quiet->answer) && mask->answer)
+ fprintf (stdout, " (masked cells = -1)");
+
+
+ /* Provide user with info about actual sampling frequencies */
+ if (strcmp (type->answer, "sites") == 0)
+ sample = attempted_sample;
+ if (! quiet->answer)
+ {
+ fprintf (stdout, "\n");
+ fflush (stdout);
+ }
+ if (!silent->answer)
+ {
+ if (absolute->answer)
+ {
+ fprintf (stdout,
+ "\nSample: requested = %8.5lf attempted = %8.5lf actual = %8.5lf\n", sample, attempted_sample, actual_sample);
+ }
+ else
+ {
+ fprintf (stdout,
+ "\nSample: requested = %8.5lf%% attempted = %8.5lf%% actual = %8.5lf%%\n", sample, attempted_sample, actual_sample);
+ }
+ fflush (stdout);
+ }
+
+
+ /* Create support files */
+
+ if (! quiet->answer)
+ fprintf (stdout, "\nCreating support files "); fflush (stdout);
+
+
+ /* Create category file for output map */
+ G_read_cats (output->answer, current_mapset, &cats);
+ G_set_cats_fmt ("$1 cell$?s", 1.0, 0.0, 0.0, 0.0, &cats);
+ G_write_cats (output->answer, &cats);
+ G_free_cats (&cats);
+
+ /* Create history support file for output map */
+ G_short_history (output->answer, "raster", &history);
+ history.edlinecnt = 4;
+ sprintf (history.edhist [0], "Sample type = %s", type->answer);
+ sprintf (history.edhist[1], "Requested = %8.5lf%%", sample);
+ sprintf (history.edhist[2], "Attempted = %8.5lf%%", attempted_sample);
+ sprintf (history.edhist[3], "Actual = %8.5lf%%", actual_sample);
+ if (G_write_history (output->answer, &history) == -1)
+ G_warning ("Failed to write history file ");
+
+ if (! quiet->answer)
+ fprintf (stdout, "\nJob finished\n");
+
+ return (0);
+}
Added: grass-addons/grass6/raster/r.cva/make_list3.c
===================================================================
--- grass-addons/grass6/raster/r.cva/make_list3.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/make_list3.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,149 @@
+/****************************************************************/
+/* */
+/* make_list3.c */
+/* */
+/* This function adds a new point to the point list */
+/* for any segment of the map. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by Mark Lake on 14/12/95 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Add a new point to the list of cells that are
+ * potentially visible from the current viewpoint,
+ * i.e. which are within the maximum visible
+ * distance.
+ *
+ * Version 3.0, differs from v.1.0 as follows:-
+ *
+ * 1) Calls make_point.c version 3;
+ * 2) Passes cell_no to make_point.c.
+ *
+ * Implicated files:-
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) Recieves cell_no and passes it on to make_point.c.
+ *
+ * Called by:
+ *
+ * 1) segment
+ *
+ * Calls:
+ *
+ * 1) make_point
+ * 2) find_orientation
+ * 3) find_inclination
+ */
+
+#include <stdlib.h>
+
+#include <grass/segment.h>
+#include <grass/gis.h>
+
+#include "config.h"
+#include "point.h"
+
+#define NEXT_PT PRESENT_PT->next
+
+struct point *make_list(struct point *head, int y, int x,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_p,
+ double viewpt_elev, int quadrant,
+ int row_viewpt, int col_viewpt, CELL cell_no,
+ RASTER_MAP_TYPE data_type)
+{
+ double del_x,del_y,dist,sqrt(),orientation,inclination,
+ find_orientation(), find_inclination();
+ struct point *make_point();
+ static struct point *PRESENT_PT;
+ extern struct Cell_head window;
+ extern double AZIMUTH1;
+ extern double AZIMUTH2;
+ extern double VERT1;
+ extern double VERT2;
+ extern int AZIMUTH1_SET;
+ extern int AZIMUTH2_SET;
+ extern int VERT1_SET;
+ extern int VERT2_SET;
+ extern double max_dist;
+
+ double v_angle, h_angle;
+
+ del_x= abs(x) ;
+ del_y= abs(y) ;
+
+ dist = sqrt(del_x*del_x +del_y*del_y) * window.ns_res;
+
+ /* if distance from viewpt is greater than the max */
+ /* range specified, neglect that point */
+ if (dist > max_dist) {
+ return(head);
+ }
+
+ /* otherwise find orientation and inclination */
+ orientation = find_orientation(x,y,quadrant);
+ inclination = find_inclination(x,y,viewpt_elev,seg_in_p,
+ row_viewpt,col_viewpt, data_type);
+
+ /* check, if limits are imposed on orientation and inclination */
+ if ( AZIMUTH1_SET ) {
+ h_angle = orientation * (180.0 / 3.14159265359);
+ h_angle = 360 - h_angle;
+ h_angle = h_angle + 90;
+ if (h_angle > 360) {
+ h_angle = h_angle - 360;
+ }
+ if (h_angle < AZIMUTH1) {
+ return (head);
+ }
+ }
+ if ( AZIMUTH2_SET ) {
+ h_angle = orientation * (180.0 / 3.14159265359);
+ h_angle = 360 - h_angle;
+ h_angle = h_angle + 90;
+ if (h_angle > 360) {
+ h_angle = h_angle - 360;
+ }
+ if (h_angle > AZIMUTH2) {
+ return (head);
+ }
+ }
+ if ( VERT1_SET ) {
+ v_angle = inclination * (180.0 / 3.14159265359);
+ if (v_angle > VERT1) {
+ return (head);
+ }
+ }
+ if ( VERT2_SET ) {
+ v_angle = inclination * (180.0 / 3.14159265359);
+ if (v_angle < VERT2) {
+ return (head);
+ }
+ }
+
+
+#ifdef DEBUG
+ printf ("\nlist y = %d, x = %d", y, x);
+#endif
+
+ if(head == NULL)
+ { /* first point ? */
+ head= make_point(orientation,inclination,y,x,
+ row_viewpt,col_viewpt,seg_out_p,cell_no);
+ PRESENT_PT = head;
+ }
+ else
+ { /* add new point to tail of list */
+ NEXT_PT = make_point(orientation,inclination,y,x,
+ row_viewpt,col_viewpt,seg_out_p,cell_no);
+ PRESENT_PT = NEXT_PT ;
+ }
+ return(head);
+
+}
+
+/********** END OF FUNCTION "MAKE_LIST" *************************/
Added: grass-addons/grass6/raster/r.cva/make_point3.c
===================================================================
--- grass-addons/grass6/raster/r.cva/make_point3.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/make_point3.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,88 @@
+/****************************************************************/
+/* */
+/* make_point3.c */
+/* */
+/* This function allocates memory space for a new point, */
+/* initializes the various fields using the values of */
+/* the parameters passed and returns the address of this */
+/* new point so that it can be attached in the linked */
+/* list. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by Mark Lake on 14/12/95 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Create a new point on a segment list.
+ *
+ * Version 3.0, differs from v.1.0 as follows:-
+ *
+ * 1) Receives cell_no and uses it to determine whether or not
+ * to initialise segment file with value 0. Allows main3.c
+ * to count only those cells on a list which have the value 0.
+ * Since main3.c writes cell_no once a cell has been counted,
+ * this ensures that cells which occur on more than one list
+ * are counted only once. If the segment file already contains
+ * the current cell_no then the cell has already been deleted
+ * during the current los analysis and should not be
+ * reinitialised.
+ *
+ * Implicated files:-
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) Initialises segment file with value zero if not already
+ * deleted during the the current los analysis.
+ *
+ * Called by:
+ *
+ * 1) make_list
+ *
+ * Calls:
+ *
+ * 1) None.
+ */
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+#define NULL_HERE 0
+
+#define NEW_PT_X NEW_PT->x
+#define NEW_PT_Y NEW_PT->y
+#define NEW_PT_ORIENTATION NEW_PT->orientation
+#define NEW_PT_INCLINATION NEW_PT->inclination
+#define NEXT_NEW_PT NEW_PT->next
+
+struct point *make_point( double orientation, double inclination,
+ int y, int x, int row_viewpt, int col_viewpt,
+ SEGMENT *seg_out_p, CELL cell_no)
+{
+ struct point *NEW_PT;
+ CELL data_read, data_to_write;
+ char *value_read, *value_to_write;
+
+ NEW_PT = (struct point *) G_malloc(sizeof(struct point));
+ NEW_PT_ORIENTATION = orientation;
+ NEW_PT_INCLINATION = inclination;
+ NEW_PT_Y = y;
+ NEW_PT_X = x;
+ NEXT_NEW_PT = NULL_HERE;
+
+ value_read = (char*) &data_read;
+ value_to_write = (char*) &data_to_write;
+ data_to_write = 0;
+ segment_get (seg_out_p, value_read,
+ row_viewpt - y, x + col_viewpt);
+ if (data_read != cell_no)
+ segment_put (seg_out_p, value_to_write,
+ row_viewpt - y, x + col_viewpt);
+ return(NEW_PT);
+}
+
+/********* END OF FUNCTION "MAKE_POINT" *************************/
Added: grass-addons/grass6/raster/r.cva/point.h
===================================================================
--- grass-addons/grass6/raster/r.cva/point.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/point.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,38 @@
+#ifndef __POINT_H__
+#define __POINT_H__
+
+
+/****************************************************************/
+/* */
+/* point.h in ~/src/Glos */
+/* */
+/* This header file defines the data structure of a */
+/* point (structure containing various attributes of */
+/* a grid cell). */
+/* */
+/****************************************************************/
+
+ struct point {
+
+ double orientation;
+ /* horizontal angle(degrees) measured from +ve x-axis */
+
+ double inclination;
+ /* vertical angle(degrees) from the viewing point */
+
+ double probability;
+ /* for probabilistic viewsheds: this stores the probability
+ of the point being visible */
+
+ int x; /* x-coor measured from viewing point location */
+ int y; /* y-coor measured from viewing point location */
+
+ struct point *next; /* pointer to next point in list*/
+ struct point *previous; /* ptr to previous pt. in list */
+
+ };
+
+/****************************************************************/
+
+#endif
+
Added: grass-addons/grass6/raster/r.cva/point_sample.c
===================================================================
--- grass-addons/grass6/raster/r.cva/point_sample.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/point_sample.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,329 @@
+/***********************************************************************/
+/*
+ point_sample.c
+
+ Updated by Mark Lake on 17/8/01 for GRASS 5.x
+
+ CONTAINS
+
+ 1) point_sample */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ point_sample
+
+ Called from:
+
+ 1) main main.c
+
+ Calls:
+
+ 1) line_of_sight line_of_site.c
+ 2) count_visible_cells count_visible_cells.c
+ 3) cumulate_visible_cells cumulate_visible_cells.c
+ 4) delete_lists delete_lists.c */
+
+/***********************************************************************/
+
+#include <stdlib.h>
+#include <string.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+#include <grass/Vect.h>
+#include <grass/dbmi.h>
+
+#include "config.h"
+#include "point.h"
+#include "global_vars.h"
+
+#include "init_segment_file.h"
+#include "line_of_sight.h"
+#include "count_visible_cells.h"
+#include "cumulate_visible_cells.h"
+#include "delete_lists.h"
+
+
+void point_sample (struct Cell_head *window, int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ char *site_file, int terse, RASTER_MAP_TYPE data_type,
+ int ignore_att)
+{
+ int row_viewpt, col_viewpt;
+
+ /* 'value' stores a byte-sized (CELL) item from the mask segment file(s) */
+ void *value = NULL;
+
+ CELL cell_no;
+ double viewpt_elev = 0;
+ struct point *heads[16];
+ long int sites_in_region, sites_of_interest, cells_in_map;
+ int null_value;
+ char *site_mapset;
+ CELL mask = 0;
+ long num_sites = 0;
+ /* the following are used to store different raster map types */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+ struct Map_info in_vect_map;
+ struct line_pnts *vect_points;
+ struct line_cats *vect_cats;
+ int cur_type;
+ char errmsg [200];
+ double x,y,z;
+ int n_points = 1;
+
+ /* site attribute management */
+ int i;
+ struct field_info *field;
+ char buf[5000], *colname;
+ int dbcol, dbncols, ctype, sqltype, more;
+ dbString sql, str;
+ dbDriver *driver;
+ dbHandle handle;
+ dbCursor cursor;
+ dbTable *table;
+ dbColumn *column;
+ dbValue *dbvalue;
+
+
+ vect_points = Vect_new_line_struct ();
+ vect_cats = Vect_new_cats_struct ();
+
+ if ((site_mapset = G_find_vector2 (site_file, "")) == NULL) {
+ sprintf (errmsg, "Could not find input vector map %s\n", site_file);
+ G_fatal_error (errmsg);
+
+ }
+
+ Vect_set_open_level (1);
+ if (1 > Vect_open_old (&in_vect_map, site_file, site_mapset)) {
+ sprintf (errmsg, "Could not open input vector points.\n");
+ G_fatal_error (errmsg);
+ }
+
+
+ /* Initialise output segment file with appropriate data value */
+ init_segment_file (nrows, ncols, seg_out_2_p, seg_patt_p);
+
+ /* Initialise variables */
+ cells_in_map = nrows * ncols;
+ sites_in_region = 0;
+ sites_of_interest = 0;
+ cell_no = 0;
+
+ /* filter vector objects to current region and point types only */
+ Vect_set_constraint_region (&in_vect_map, window->north, window->south,
+ window->east, window->west, 0.0, 0.0);
+ Vect_set_constraint_type (&in_vect_map, GV_POINT);
+
+ /* calculate number of vector points in map and current region */
+ while ((cur_type = Vect_read_next_line (&in_vect_map, vect_points, NULL)) > 0) {
+ num_sites ++;
+ }
+
+ /* rewind vector points file: close and re-open */
+ Vect_close (&in_vect_map);
+ Vect_open_old(&in_vect_map, site_file, site_mapset);
+ Vect_set_constraint_type (&in_vect_map, GV_POINT);
+ Vect_set_constraint_region (&in_vect_map, window->north, window->south,
+ window->east, window->west, 0.0, 0.0);
+
+
+ /* Traverse site list */
+ while ((cur_type = Vect_read_next_line (&in_vect_map, vect_points, vect_cats)) > 0) {
+ mask = 0;
+ cell_no ++;
+
+ /* check for site attributes */
+ for (i = 0; i < vect_cats->n_cats; i++) {
+ field = Vect_get_field( &in_vect_map, vect_cats->field[i]);
+ if ((field != NULL) && (!ignore_att)) {
+ db_init_string (&sql);
+ driver = db_start_driver(field->driver);
+ db_init_handle (&handle);
+ db_set_handle (&handle, field->database, NULL);
+ if (db_open_database(driver, &handle) == DB_OK){
+ /* sql query */
+ sprintf (buf, "select * from %s where %s = %d", field->table,
+ field->key,
+ vect_cats->cat[i]);
+ db_set_string (&sql, buf);
+
+ db_open_select_cursor(driver, &sql, &cursor, DB_SEQUENTIAL);
+ table = db_get_cursor_table (&cursor);
+ db_fetch (&cursor, DB_NEXT, &more );
+ dbncols = db_get_table_number_of_columns (table);
+ for( dbcol = 0; dbcol < dbncols; dbcol++) {
+ column = db_get_table_column(table, dbcol);
+ sqltype = db_get_column_sqltype (column);
+ ctype = db_sqltype_to_Ctype(sqltype);
+ dbvalue = db_get_column_value(column);
+ db_convert_value_to_string( dbvalue, sqltype, &str);
+ colname = db_get_column_name (column);
+ /* check if we have some of the supported attributes */
+ if (ctype == DB_C_TYPE_DOUBLE) {
+ if (!strcmp (colname,"SPOT")){
+ SPOT = atof (db_get_string (&str));
+ SPOT_SET = 1;
+ }
+ if (!strcmp (colname,"OFFSETA")){
+ obs_elev = atof (db_get_string (&str));
+ }
+ if (!strcmp (colname,"OFFSETB")){
+ OFFSETB_SET =1;
+ OFFSETB = atof (db_get_string (&str));
+ }
+ if (!strcmp (colname,"AZIMUTH1")){
+ AZIMUTH1 = atof (db_get_string (&str));
+ AZIMUTH1_SET = 1;
+ }
+ if (!strcmp (colname,"AZIMUTH2")){
+ AZIMUTH2 = atof (db_get_string (&str));
+ AZIMUTH2_SET = 1;
+ }
+ if (!strcmp (colname,"VERT1")){
+ VERT1 = atof (db_get_string (&str));
+ VERT1_SET = 1;
+ }
+ if (!strcmp (colname,"VERT2")){
+ VERT2 = atof (db_get_string (&str));
+ VERT2_SET = 1;
+ }
+ if (!strcmp (colname,"RADIUS1")){
+ RADIUS1 = atof (db_get_string (&str));
+ }
+ if (!strcmp (colname,"RADIUS2")){
+ max_dist = atof (db_get_string (&str));
+ }
+ }
+ }
+ db_close_cursor(&cursor);
+ db_close_database(driver);
+ db_shutdown_driver(driver);
+ } else {
+ db_shutdown_driver(driver);
+ }
+ }
+ }
+
+ Vect_copy_pnts_to_xyz (vect_points, &x, &y, &z, &n_points);
+ /* If site falls within current region */
+ sites_in_region ++;
+
+ /* Get array coordinates of viewpoint */
+ /* row_viewpt = (window->north - site->north) / window->ns_res;
+ col_viewpt = (site->east - window->west) / window->ew_res;
+ */
+ row_viewpt = G_northing_to_row ( y, window );
+ col_viewpt = G_easting_to_col ( x, window );
+
+ /* Check whether viewpoint is in area of interest */
+ if (patt_flag_v) {
+ value = (char*) &mask;
+ segment_get(seg_patt_p_v, value, row_viewpt, col_viewpt);
+ } else
+ mask = 1;
+
+ /* Only include if viewpoint is in area of interest */
+ if (mask == 1) {
+ sites_of_interest ++;
+
+ /* We do not check for duplicates because it is possible
+ in cases of low resolution that two or more sites with
+ different map coordinates could have the same array
+ coordinates */
+
+ /* Read elevation of current view point */
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ null_value = 0;
+ segment_get (seg_in_p, value, row_viewpt, col_viewpt); /* height now in value */ /* and in read_viewpt_elev */
+
+ if ( data_type == CELL_TYPE ) {
+ viewpt_elev = c_value + obs_elev;
+ if (G_is_c_null_value (&c_value)) {
+ null_value = 1;
+ }
+ }
+ if ( data_type == FCELL_TYPE ) {
+ viewpt_elev = f_value + obs_elev;
+ if (G_is_f_null_value (&f_value)) {
+ null_value = 1;
+ }
+ }
+ if ( data_type == DCELL_TYPE ) {
+ viewpt_elev = d_value + obs_elev;
+ if (G_is_d_null_value (&d_value)) {
+ null_value = 1;
+ }
+ }
+
+ /* skip sites on NULL dem cells */
+ if (!null_value) {
+ if (SPOT_SET) {
+ if (ADJUST) {
+ if (SPOT >= viewpt_elev) {
+ viewpt_elev = SPOT + obs_elev;
+ }
+ } else {
+ viewpt_elev = SPOT;
+ }
+ }
+
+ /* Do line of sight analysis from current viewpoint */
+ /* this updates the information in 'heads', defined in 'point.h' */
+ line_of_sight (cell_no, row_viewpt, col_viewpt,
+ nrows, ncols, viewpt_elev,
+ heads, seg_in_p, seg_out_1_p, seg_patt_p,
+ patt_flag, data_type);
+
+ /* Calculate viewshed */
+ if (! from_cells) {
+ count_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+ } else
+ cumulate_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+
+ /* Free memory used in line-of-sight analysis */
+ delete_lists (heads);
+ }
+
+ /* Print progress message */
+ if (!terse) {
+ G_percent (sites_of_interest, num_sites, 1);
+ }
+ } /* END (if (mask == 1)) */
+ } /* END (loop through sites list) */
+ Vect_close (&in_vect_map);
+
+ /* Calculate_sample */
+ /* % of cells in region */
+ *attempted_sample = 100.0 / (double) cells_in_map * (double) sites_in_region;
+
+ /* % of cells in area of interest */
+ *actual_sample = 100.0 / (double) cells_in_map * (double) sites_of_interest;
+
+#ifdef DEBUG
+ fprintf (stdout, "\nSites in: region = %ld, area of interest = %ld",
+ sites_in_region,
+ sites_of_interest);
+ fflush (stdout);
+#endif
+}
Added: grass-addons/grass6/raster/r.cva/point_sample.h
===================================================================
--- grass-addons/grass6/raster/r.cva/point_sample.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/point_sample.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,15 @@
+#ifndef __POINT_SAMPLE_H__
+#define __POINT_SAMPLE_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void point_sample (struct Cell_head *window, int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ char *site_file, int terse, RASTER_MAP_TYPE data_type,
+ int ignore_att);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/pts_elim3.c
===================================================================
--- grass-addons/grass6/raster/r.cva/pts_elim3.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/pts_elim3.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,550 @@
+/****************************************************************/
+/* */
+/* pts_elim3.c */
+/* */
+/* This function prunes a linked list of all points */
+/* picked up from a map segment to leave only those */
+/* points in the list that are visible from the viewpt. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by Mark Lake on 14/12/95 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Eliminate invisible cells from a segment list.
+ *
+ * Version 3.0, differs from v.1.0 as follows:-
+ *
+ * 1) Calls delete.c version 3;
+ * 2) Receives cell_no and passes it to delete.c.
+ *
+ * Implicated files:-
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:-
+ *
+ * 1) Recieves cell_no and passes it on to delete.c.
+ *
+ * Called by:
+ *
+ * 1) segment
+ *
+ * Calls:
+ *
+ * 1) find_orientation
+ * 2) find_inclination
+ * 3) delete
+ */
+
+#include <stdlib.h>
+#include <math.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "pts_elim3.h"
+#include "point.h"
+#include "radians.h"
+#include "global_vars.h"
+
+#define SECOND_PT head->next
+
+#define NEXT_BLOCKING_PT BLOCKING_PT->next
+#define BLOCKING_PT_X BLOCKING_PT->x
+#define BLOCKING_PT_Y BLOCKING_PT->y
+#define BLOCKING_PT_INCLINATION BLOCKING_PT->inclination
+#define BLOCKING_PT_ORIENTATION BLOCKING_PT->orientation
+
+#define NEXT_CHECKED_PT CHECKED_PT->next
+#define CHECKED_PT_X CHECKED_PT->x
+#define CHECKED_PT_Y CHECKED_PT->y
+#define CHECKED_PT_INCLINATION CHECKED_PT->inclination
+#define CHECKED_PT_ORIENTATION CHECKED_PT->orientation
+
+
+/****************************************************************/
+/* */
+/* This function takes a linked list of points picked */
+/* up from a segment of the map and deletes all the points */
+/* from the list that are not visible from the viewing pt. */
+/* */
+/****************************************************************/
+
+struct point *hidden_point_elimination (struct point *head,
+ double viewpt_elev,
+ SEGMENT *seg_in_p,
+ SEGMENT *seg_out_p,
+ SEGMENT *seg_patt_p,
+ int quadrant,
+ int sign_on_y, int sign_on_x,
+ int row_viewpt, int col_viewpt,
+ int patt_flag, CELL cell_no,
+ RASTER_MAP_TYPE data_type)
+{
+ struct point *CHECKED_PT, *BLOCKING_PT, *delete();
+ double orientation_neighbor_1,orientation_neighbor_2,
+ find_orientation(),inclination_neighbor_1,
+ inclination_neighbor_2,interpolated_inclination,
+ find_inclination(),find_inclination2(),correct_neighbor_inclination,
+ correct_neighbor_orientation,fabs();
+ int correct_neighbor_x,correct_neighbor_y,neighbor_1_y,
+ neighbor_1_x,neighbor_2_x,neighbor_2_y, uu,vv;
+ CELL mask;
+ char *value;
+ int do_check = 1;
+ extern double RADIUS1;
+ extern int OFFSETB_SET;
+ double del_x, del_y,dist,atan(),sqrt();
+
+
+ uu = (sign_on_y + sign_on_x)/2;
+ vv = (sign_on_y - sign_on_x)/2;
+
+ value = (char *) &mask;
+
+
+#ifdef DEBUG
+ printf ("\nrow = %d, col = %d", row_viewpt, col_viewpt);
+ fflush (stdout);
+#endif
+
+ /* move blocking pt. from the 2nd pt till the end */
+ for(BLOCKING_PT = SECOND_PT;
+ BLOCKING_PT != NULL;
+ BLOCKING_PT = NEXT_BLOCKING_PT)
+ {
+
+#ifdef DEBUG
+ printf ("\nBlocking point y = %d, x = %d",
+ BLOCKING_PT_Y, BLOCKING_PT_X);
+ fflush (stdout);
+#endif
+
+ /* calc coors of the two immediate neighbors on either */
+ /* side of the blocking point */
+
+ if(BLOCKING_PT_X == 0 || BLOCKING_PT_Y == 0)
+ {
+ neighbor_1_x = BLOCKING_PT_X - vv;
+ neighbor_1_y = BLOCKING_PT_Y + uu;
+
+ neighbor_2_x = BLOCKING_PT_X + uu;
+ neighbor_2_y = BLOCKING_PT_Y + vv;
+ }
+ else
+ {
+ neighbor_1_x = BLOCKING_PT_X - uu;
+ neighbor_1_y = BLOCKING_PT_Y - vv;
+
+ neighbor_2_x = BLOCKING_PT_X + vv;
+ neighbor_2_y = BLOCKING_PT_Y - uu;
+ }
+
+ /* find orientation and inclination for both neighbors */
+ orientation_neighbor_1 =
+ find_orientation(neighbor_1_x,neighbor_1_y,quadrant);
+
+ orientation_neighbor_2 =
+ find_orientation(neighbor_2_x,neighbor_2_y,quadrant);
+
+ inclination_neighbor_1 =
+ find_inclination(neighbor_1_x,neighbor_1_y,viewpt_elev,
+ seg_in_p,row_viewpt,col_viewpt, data_type);
+
+ inclination_neighbor_2 =
+ find_inclination(neighbor_2_x,neighbor_2_y,viewpt_elev,
+ seg_in_p,row_viewpt,col_viewpt, data_type);
+
+ /* check all points behind the blocking point */
+ for(CHECKED_PT = head;
+ CHECKED_PT != BLOCKING_PT;
+ CHECKED_PT = NEXT_CHECKED_PT)
+ {
+
+#ifdef DEBUG
+ printf ("\nChecked point y = %d, x = %d",
+ CHECKED_PT_Y, CHECKED_PT_X);
+ fflush (stdout);
+#endif
+
+
+ /* if pattern layer specified, check to see if checked */
+ /* point is of interest. If not, delete it from list */
+ if(patt_flag == 1)
+ {
+ segment_get(seg_patt_p,value,
+ row_viewpt-CHECKED_PT_Y,col_viewpt+CHECKED_PT_X);
+
+ if(mask == 0)
+ {
+ head=delete(CHECKED_PT,head,seg_out_p,
+ row_viewpt,col_viewpt,cell_no);
+ goto next_iter;
+ }
+ }
+
+ /* if the OFFSETB parameter is set, we will raise the height of only the
+ target point under consideration, leaving all blocking points the
+ way they are. This improves visibility of current target point and
+ simulates an object of a known height */
+ if (OFFSETB_SET) {
+ CHECKED_PT_INCLINATION = find_inclination2(CHECKED_PT_X,CHECKED_PT_Y,viewpt_elev,
+ seg_in_p,row_viewpt,col_viewpt, data_type);
+ }
+
+ do_check = 1;
+
+ /* no need to check if this condition holds true */
+ if(BLOCKING_PT_INCLINATION <= CHECKED_PT_INCLINATION) {
+ do_check = 0;
+ }
+ if (BLOCKING_PT_INCLINATION == NULLPT) {
+ do_check = 0;
+ }
+
+ /* delete this point, if it lies closer than the minimum distance */
+ del_x = abs(CHECKED_PT_X);
+ del_y = abs(CHECKED_PT_Y);
+ dist=sqrt(del_x * del_x + del_y * del_y)*window.ns_res;
+ if ( dist < RADIUS1) {
+ head= delete(CHECKED_PT,head,seg_out_p,row_viewpt,col_viewpt,cell_no);
+ do_check = 0;
+ }
+
+ if (do_check == 1)
+
+ { /* otherwise, proceed to check */
+
+
+ if(CHECKED_PT_ORIENTATION == BLOCKING_PT_ORIENTATION)
+ {
+ head= delete(CHECKED_PT,head,seg_out_p,
+ row_viewpt,col_viewpt,cell_no);
+ /* if checked point directly behind, delete it */
+ }
+ else
+ { /* if checked point not directly behind, check */
+
+ /* find the coors of the actual neighbor that might be */
+ /* required for interpolation. */
+ if(CHECKED_PT_ORIENTATION > BLOCKING_PT_ORIENTATION)
+ {
+ correct_neighbor_x = neighbor_1_x;
+ correct_neighbor_y = neighbor_1_y;
+ correct_neighbor_inclination = inclination_neighbor_1;
+ correct_neighbor_orientation = orientation_neighbor_1;
+ }
+ else
+ {
+ correct_neighbor_x = neighbor_2_x;
+ correct_neighbor_y = neighbor_2_y;
+ correct_neighbor_inclination = inclination_neighbor_2;
+ correct_neighbor_orientation = orientation_neighbor_2;
+ }
+
+ if(fabs(BLOCKING_PT_ORIENTATION-CHECKED_PT_ORIENTATION) <
+ fabs(BLOCKING_PT_ORIENTATION-correct_neighbor_orientation))
+
+ { /* yes, the point neighboring the blocking point */
+ /* must be taken into consideration */
+
+ if(CHECKED_PT_Y == correct_neighbor_y &&
+ CHECKED_PT_X == correct_neighbor_x); /* same point */
+
+ else
+ { /* CHECK !! */
+
+
+ /* if the checked point's inclination is even lower */
+ /* than that of the blocking pt.'s neighbor, blocked */
+ if(CHECKED_PT_INCLINATION < correct_neighbor_inclination)
+ {
+ head= delete(CHECKED_PT,head,seg_out_p,
+ row_viewpt,col_viewpt,cell_no);
+ }
+
+ else
+ { /* INTERPOLATION */
+
+
+ interpolated_inclination= BLOCKING_PT_INCLINATION +
+ (CHECKED_PT_ORIENTATION - BLOCKING_PT_ORIENTATION)/
+ (correct_neighbor_orientation-BLOCKING_PT_ORIENTATION)*
+ (correct_neighbor_inclination - BLOCKING_PT_INCLINATION);
+
+ if(CHECKED_PT_INCLINATION < interpolated_inclination)
+ { /* interpolated point blocks */
+ head= delete(CHECKED_PT,head,seg_out_p,
+ row_viewpt,col_viewpt,cell_no);
+ }
+ }
+ }
+ }
+ }
+ }
+ next_iter:
+ ;
+ } /* end of loop over points to be checked */
+
+
+
+
+ /* if pattern layer specified, check if blocking point */
+ /* itself is an area of interest. If not, of no use */
+ if(patt_flag == 1)
+ {
+ segment_get(seg_patt_p, value,row_viewpt- BLOCKING_PT_Y,
+ col_viewpt+BLOCKING_PT_X);
+ if(mask == 0)
+ {
+
+ /* Commenting out the following fixes a bug in r.los.
+ In that program the 8 cells around the viewpoint
+ are marked as visible (when visible)
+ even if they fall outside the area of interest
+ specified by the patt_map. This occurs because
+ these cells are always the last blocking points
+ on the segment lists, and therefore don't get
+ deleted when they should. This fix allows them
+ to be deleted, but it required modifications
+ to delete, in delete3.c. MWL 25/6/99 */
+
+ /* if (NEXT_BLOCKING_PT != NULL) */
+
+ head = delete(BLOCKING_PT, head, seg_out_p,
+ row_viewpt, col_viewpt,cell_no);
+
+ }
+ }
+
+ } /* end of loop over blocking points */
+
+ return(head);
+
+}
+
+/*********** END OF FUNCTION "HIDDEN_POINT_ELIMINATION" *********/
+
+/* Called by:-
+ *
+ * 1) make_list()
+ * 2) hidden_point_elimination()
+ *
+ * Calls:-
+ *
+ * 1) None.
+ */
+
+
+/****************************************************************/
+/* */
+/* This function finds the orientation of a point if */
+/* provided with the number of the quadrant and the */
+/* coordinates of that point. */
+/* */
+/****************************************************************/
+
+double find_orientation(int x, int y, int quadrant)
+{
+ double del_x,del_y,atan(),angle;
+ int abs();
+
+ del_x = abs(x) ;
+ del_y = abs(y);
+
+ if(del_x == 0.0) angle = PIBYTWO;
+ else angle = atan(del_y/del_x) ;
+
+ switch(quadrant)
+ {
+ case 1 :
+ break;
+ case 2 :
+ angle = PI - angle;
+ break;
+ case 3 :
+ angle = PI + angle;
+ break;
+ case 4 :
+ angle = TWOPI - angle;
+ break;
+ default :
+ break;
+ }
+
+ return(angle);
+
+} /* END OF FUNCTION ANGLE */
+
+/************* END OF FUNCTION "FIND_ORIENTATION" ***************/
+
+/* Called by:-
+ *
+ * 1) make_list()
+ * 2) hidden_point_elimination()
+ *
+ * Calls:-
+ *
+ * 1) None.
+ */
+
+/****************************************************************/
+/* */
+/* This function calculates the vertical angle of a point */
+/* with respect to the viewing pt. */
+/* */
+/****************************************************************/
+
+double find_inclination(int x, int y, double viewpt_elev,
+ SEGMENT *seg_in_p,
+ int row_viewpt, int col_viewpt, RASTER_MAP_TYPE data_type)
+{
+ double del_x, del_y,dist,atan(),sqrt();
+ int abs();
+ double dest_elev = 0.0;
+ extern struct Cell_head window;
+ void *value = NULL;
+ /* these vars can store one data value from the elevation input map, */
+ /* which could be CELL, DCELL or FCELL type. */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+ del_x = abs(x) ;
+ del_y = abs(y) ;
+
+ dist=sqrt(del_x * del_x + del_y * del_y)*window.ns_res;
+
+ /* this takes care, that target elevation has the right format,
+ depending on type of input DEM */
+
+ if (data_type == CELL_TYPE) {
+ value = (CELL*) &c_value;
+ }
+ if (data_type == FCELL_TYPE) {
+ value = (FCELL*) &f_value;
+ }
+ if (data_type == DCELL_TYPE) {
+ value = (DCELL*) &d_value;
+ }
+
+ /* read value from elevation input map, convert to appropriate */
+ /* map type */
+ segment_get(seg_in_p,value,row_viewpt-y,x+col_viewpt);
+
+ if (data_type == CELL_TYPE) {
+ if ( G_is_c_null_value (&c_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = c_value;
+ }
+ }
+
+ if (data_type == FCELL_TYPE) {
+ if ( G_is_f_null_value (&f_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = f_value;
+ }
+ }
+
+ if (data_type == DCELL_TYPE) {
+ if ( G_is_d_null_value (&d_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = d_value;
+ }
+ }
+
+ /* CURVATURE CORRECTION */
+ /* decrease height of target point */
+ if (DIST_CC > 0.0) {
+ if (dist >= DIST_CC) {
+ dest_elev = dest_elev - ((dist*dist) / (2 * Re));
+ }
+ }
+
+
+ return(atan((dest_elev - viewpt_elev) / dist));
+}
+
+
+/* THIS IS CALLED IF ONLY THE TARGET POINT UNDER CURRENT CONSIDERATION IS TO
+ BE RAISED BY 'OFFSETB' AMOUNT */
+double find_inclination2(int x, int y, double viewpt_elev,
+ SEGMENT *seg_in_p,
+ int row_viewpt, int col_viewpt, RASTER_MAP_TYPE data_type)
+{
+ double del_x, del_y,dist,atan(),sqrt();
+ int abs();
+ double dest_elev = 0.0;
+ extern struct Cell_head window;
+ void *value = NULL;
+ /* these vars can store one data value from the elevation input map, */
+ /* which could be CELL, DCELL or FCELL type. */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+ del_x = abs(x) ;
+ del_y = abs(y) ;
+
+ dist=sqrt(del_x * del_x + del_y * del_y)*window.ns_res;
+
+ /* this takes care, that target elevation has the right format,
+ depending on type of input DEM */
+
+ if (data_type == CELL_TYPE) {
+ value = (CELL*) &c_value;
+ }
+ if (data_type == FCELL_TYPE) {
+ value = (FCELL*) &f_value;
+ }
+ if (data_type == DCELL_TYPE) {
+ value = (DCELL*) &d_value;
+ }
+
+ /* read value from elevation input map, convert to appropriate */
+ /* map type */
+ segment_get(seg_in_p,value,row_viewpt-y,x+col_viewpt);
+
+ if (data_type == CELL_TYPE) {
+ if ( G_is_c_null_value (&c_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = c_value + OFFSETB;
+ }
+ }
+
+ if (data_type == FCELL_TYPE) {
+ if ( G_is_f_null_value (&f_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = f_value + OFFSETB;
+ }
+ }
+
+ if (data_type == DCELL_TYPE) {
+ if ( G_is_d_null_value (&d_value) ) {
+ return (NULLPT);
+ } else {
+ dest_elev = d_value + OFFSETB;
+ }
+ }
+
+
+ /* CURVATURE CORRECTION */
+ /* decrease height of target point */
+ if (DIST_CC > 0.0) {
+ if (dist >= DIST_CC) {
+ dest_elev = dest_elev - ((dist*dist) / (2 * Re));
+ }
+ }
+
+ return(atan((dest_elev - viewpt_elev) / dist));
+}
+
+/************ END OF FUNCTION "FIND_INCLINATION"*****************/
Added: grass-addons/grass6/raster/r.cva/pts_elim3.h
===================================================================
--- grass-addons/grass6/raster/r.cva/pts_elim3.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/pts_elim3.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,17 @@
+#ifndef __PTS_ELIM3_H__
+#define __PTS_ELIM3_H__
+
+int NULLPT = 999999.999999;
+
+struct point *hidden_point_elimination (struct point *head,
+ double viewpt_elev,
+ SEGMENT *seg_in_p,
+ SEGMENT *seg_out_p,
+ SEGMENT *seg_patt_p,
+ int quadrant,
+ int sign_on_y, int sign_on_x,
+ int row_viewpt, int col_viewpt,
+ int patt_flag, CELL cell_no,
+ RASTER_MAP_TYPE data_type);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/r.cva.tmp.html
===================================================================
--- grass-addons/grass6/raster/r.cva/r.cva.tmp.html (rev 0)
+++ grass-addons/grass6/raster/r.cva/r.cva.tmp.html 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,445 @@
+<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
+<html>
+<head>
+<title>GRASS GIS Manual: r.cva</title>
+<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
+<link rel="stylesheet" href="grassdocs.css" type="text/css">
+</head>
+<body bgcolor="white">
+<img src="grass_logo.png" alt="GRASS logo"><hr align=center size=6 noshade>
+<h2>NAME</h2>
+<em><b>r.cva</b></em>
+
+<TITLE>r.cva</TITLE>
+
+<body bgcolor="white">
+
+<H2>NAME</H2>
+
+<EM><b>r.cva</b></EM> - Cumulative viewshed analysis program.
+<BR>
+<EM>(GRASS Raster Program)</EM>
+
+<H2>SYNOPSIS</H2>
+<B>r.cva</B>
+<BR>
+<B>r.cva help</B>
+<BR>
+<B>r.cva</B>
+[<B>-afmnoqrs</B>]
+<B>input=</B><EM>name</EM>
+[<B>target_mask=</B><EM>name</EM>]
+[<B>viewpt_mask=</B><EM>name</EM>]
+[<B>sites=</B><EM>name</EM>]
+<B>output=</B><EM>name</EM>
+[<B>obs_elev=</B><EM>value</EM>]
+[<B>target_elev=</B><EM>value</EM>]
+[<B>max_dist=</B><EM>value</EM>]
+[<B>seed=</B><EM>value</EM>]
+[<B>sample=</B><EM>value</EM>]
+<B>type=</B><EM>value</EM>
+
+<H2>QUICKSTART</H2>
+To perform the most common cumulative viewshed analysis (for each cell of
+the input elevation raster map: count the number of observer sites from which it
+can be seen) type in the following:
+
+<p>
+
+<B>r.cva</B> <B>-f</B> <B>input=</B><EM>elevation</EM> <B>sites=</B><EM>observers</EM> <B>output=</B><EM>result</EM> <B>type=</B><EM>sites</EM>
+<P>
+'elevation' can be a raster map containing elevation values in CELL, FCELL or DCELL format. 'observers' must a map containing vector points.
+'result' will be an output map in CELL format containing the cumulative visibility counts for each cell.
+
+<A NAME = "description"><H2>DESCRIPTION</H2></A>
+
+<EM>r.cva</EM> generates a cumulative viewshed map. When run without
+the <EM>-f</EM> flag each viewpoint is marked with the number of cells in its
+viewshed (i.e. that can be seen from it). In this mode the only cells
+that contain genuine data in the output map are the viewpoints.
+Alternatively, if the <EM>-f</EM> flag is set then each cell in the
+<EM>output</EM> map is marked with the number of viewpoints from which
+it can be seen. In either case, viewpoints may be chosen in one of
+four ways. The first (`all') treats every cell as a viewpoint. This
+is so computationally demanding that it is unlikely to be feasible
+unless either the <EM>input</EM> map or <EM>max_dist</EM> is very
+small. The second (`systematic') uses a systematic sample of cells
+arranged on a rectangular grid. The third (`random') uses a random
+sample. By default the random sampling algorithm does not allow
+replacement (i.e. a cell can be treated as a viewpoint only once), but
+this behaviour can be altered by setting the <EM>-r</EM> flag. The
+fourth option (`sites') uses the points in the specified <EM>vector</EM>
+map. The options `systematic' and `random' allow the user to specify
+any sampling frequency between 0% and 100%.
+
+<p>
+
+The cells used in viewshed analysis can be constrained by one or both
+of two masks. Both masks should be binary raster maps in which the
+areas of interest are coded `1' and all other cells `0'.
+<EM>target_mask</EM> determines which cells are targets: a cell which
+is not a target will be returned as not visible from a given viewpoint
+even if it is actually visible according to the line-of-sight analysis
+(i.e. this mask performs the same function as <EM>patt_map</EM> in
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>). <EM>viewpt_mask</EM>
+constrains which cells may be treated as viewpoints. For example,
+suppose that the raster map `no_sea' codes sea (elevation = 0) as `0'
+and land (elevation > 0) as `1'. In this case, specifying `no_sea' as
+the <EM>target_mask</EM> would cause <EM>r.cva</EM> to count cells on
+land that are visible from any viewpoint. Alternatively, specifying
+`no_sea' as the <EM>viewpt_mask</EM> would count any cells visible
+from viewpoints on land. Finally, specifying `no_sea' as both the
+<EM>target_mask</EM> and <EM>viewpt_mask</EM> would count cells on
+land that are visible from viewpoints on land.
+
+<p>
+
+The geometry of the potential line of sight between each pair of cells
+may be altered by means of two offsets, <EM>obs_elev</EM> and
+<EM>target_elev</EM>. The height <EM>obs_elev</EM> is added to the
+viewpoint and the height <EM>target_elev</EM> to each target cell.
+Note that by swapping these values it is possible to produce maps
+showing the number of cells from which a viewpoint is visible, or if
+the <EM>-f</EM> flag is set, the number of viewpoints which it is
+possible to see from each map cell.
+
+<p>
+
+A particular set of viewpoint locations in a random sample can be
+repeated by specifying the same value for <EM>seed</EM>.
+Alternatively, different values can be used to produce independent
+samples for statistical testing.
+
+<p>
+
+The program can be run either non-interactively or interactively. To
+run <EM>r.cva</EM> non-interactively the user must minimally specify
+the <B>input</B> and <B>output</B> maps and the <B>type</B> of
+analysis. If the chosen <B>type</B> is `sites' then the
+<B>sites</B> list must also be specified. Any remaining parameters
+whose values have not been specified on the command line will be set
+to their default values (see <A HREF = "#parameters">below</A>). Non-interactive usage format is:
+
+<DL>
+<DD>
+<B>r.cva</B>
+[<B>-afmnoqrs</B>]
+<B>input=</B><EM>name</EM>
+[<B>target_mask=</B><EM>name</EM>]
+[<B>viewpt_mask=</B><EM>name</EM>]
+[<B>sites=</B><EM>name</EM>]
+<B>output=</B><EM>name</EM>
+[<B>obs_elev=</B><EM>value</EM>]
+[<B>target_elev=</B><EM>value</EM>]
+[<B>max_dist=</B><EM>value</EM>]
+[<B>seed=</B><EM>value</EM>]
+[<B>sample=</B><EM>value</EM>]
+<B>type=</B><EM>value</EM>
+</DL>
+
+<H2>OPTIONS</H2>
+
+<H3>Flags:</H3>
+<DL>
+
+<DT><B>-a</B>
+<DD> Treat the specified sample size as an absolute number of cells
+rather than a percentage. Note, this flag is only for use with random
+sampling.
+
+<DT><B>-f</B>
+<DD> Mark each cell in the output map with the number of viewpoints
+from which it can be seen.
+
+<DT><B>-h</B>
+<DD> Adjust height of oberserver points so that they will never fall
+below surface level. This is only useful in conjunction with the
+<b>spot</b> parameter or attribute (see <A HREF ="#atributes">below</A>).
+
+<DT><B>-i</B>
+<DD> Do not interpret control attributes in input vector map as
+described <A HREF ="#atributes">below</A>).
+
+<DT><B>-m</B> <DD> Differentiate between genuine zero data values and
+null data in the output map. When the <EM>-f</EM> flag has not been
+set this ensures that all non-viewpoint cells are coded `-1'. When
+the <EM>-f</EM> flag has been set then all cells coded `0' in the
+<EM>target_mask</EM> (if there is one) are coded `-1'. See <A HREF =
+"#notes">below</A> for guidance when to set this flag.
+
+<DT><B>-n</B> <DD> In conjunction with <B>-m</B>: write NULL values instead of '-1' to output map.
+
+<DT><B>-o</B>
+<DD> Overwrite the output raster map if it already exists.
+
+<DT><B>-q</B>
+<DD> Be quiet about progress, confirming only the output file and
+sampling frequency. This is useful when <EM>r.cva</EM> is invoked from a
+script with standard output redirected to a disk file.
+
+<DT><B>-r</B>
+<DD> Allow replacement during random sampling (i.e. allow a given cell
+to be picked as a viewpoint more than once). Note that no cumulation
+occurs for replacement viewpoints; this ensures that the resulting
+cumulative viewshed will be as though the viewpoint was selected once
+only (the final output map will, of course, be different because
+replacement effectively results in the use of fewer viewpoints for a
+given sampling frequency).
+
+<DT><B>-s</B> <DD>Be silent about progress, in which case it is
+strongly recommended that you use <EM><A
+HREF="../../html/r.info.html">r.info</A></EM> to confirm what sampling frequency
+was actually achieved.
+
+<A NAME = "parameters"><H3>Parameters:</H3></A>
+
+<DT><B>input=</B><EM>name</EM>
+<DD>Name of a raster map layer containing elevation
+data, used as program input.
+
+<DT><B>target_mask=</B><EM>name</EM> <DD>Name of a binary raster map
+layer in which target cells within the areas of interest are assigned
+the category value `1' and all other cells are assigned the category
+value `0'. A cell is treated as a target when it is viewed from a
+viewpoint (note, therefore, that viewpoints are targets when viewed
+from other viewpoints).
+
+<DT><B>viewpt_mask=</B><EM>name</EM>
+<DD>Name of a binary raster map layer in which viewpoints within the areas
+of interest are assigned the category value `1' and all other cells
+are assigned the category value `0'.
+
+<DT><B>sites=</B><EM>name</EM>
+<DD>The name of a vector map containing the sites (vector points) that are to be treated as
+viewpoints. Sites falling outside the current region are ignored.
+
+<DT><B>obs_elev=</B><EM>value</EM>
+<DD>Height of the observer (in metres) above the elevation of the
+viewpoint. Default: 1.75 m
+
+<DT><B>target_elev=</B><EM>value</EM>
+<DD>Height of the object of interest (in metres) above the elevation
+of the target cell. Default: 0.0 m
+
+<DT><B>max_dist=</B><EM>value</EM>
+<DD>Maximum distance (in metres) from the viewing point within which
+the line-of-sight analysis will be performed. Options: 0-99999
+(stated in map units). Default: 100
+
+<DT><B>seed=</B><EM>value</EM>
+<DD>The seed for the random number generator used during random
+sampling. This option allows exact reproduction of a previous result.
+Options: 0-32767. Default: 1
+
+<DT><B>sample=</B><EM>value</EM>
+<DD>The sampling frequency as a percentage of the number of map cells
+in the current region. Options: 0.0-100.0. Default: 10.0
+
+<DT><B>type=</B><EM>value</EM> <DD>The type of sampling regime used to
+select viewpoints for analysis. See <A HREF = "#description">above</A>
+for details. Options: `all', `systematic', `random', `sites'.
+
+<DT><B>spot, offseta, offsetb, azimuth1, azimuth2, vert1, vert2, radius1, radius2</B>
+<DD> See <A HREF ="#atributes">below</A>) for a description of what they do.
+
+<DT><B>curvc=</B><EM>value</EM>
+<DD>Earth curvature correction will be applied if the distance from the observer
+point exceeds the given threshold. Set this to 0.0 to disable curvature correction
+(default). While correction slows down calculations, it is recommend for
+viewing distances of more than about 3 km. The correction uses a simplified
+formula which is accurate to 1 in 10,000 for distances up to 40 km:<p>
+<i>d(height)=dist^2/(2*Re)</i><p>
+With <i>Re</i> being the radius of the earth (ca. 6,356,766 m).
+</DL>
+
+
+
+
+<A NAME = "attributes"><H3>Vector map attributes:</H3>
+
+When using 'sites' type operation mode, <EM>r.cva</EM> checks for the
+presence of a number of attributes in the input vector map. The names
+of these attributes are all uppercase and match those used by the
+viewshed module of ESRI's ArcGIS Spatial Extension (trademarks of ESRI Corp.).
+They control <EM>r.cva</EM> operation in an analogous way.
+Each attribute, if present, must store a double value. The following
+attribute names and value ranges are recognised and interpreted:
+
+<dl>
+<dt>SPOT</dt>
+ <dd>This determines absolute observer height. If this attribute is given,
+ observer height will no longer be determined from the DEM.
+ <EM>r.cva</EM> can automatically raise this value to that of the DEM
+ point below the oberserver's position, if the value in SPOT would
+ cause the observer to get below ground level (use the <b>-h</b> flag).</dd>
+<dt>OFFSETA</dt>
+ <dd>Use this to add an offset to each oberver point's height. This is
+ just an alias for the <b>obs_elev</b> option and exists for
+ compatibility reasons.</dd>
+<dt>OFFSETB</dt>
+ <dd>The value stored in 'OFFSETB' will be added to each target point
+ at the moment it is considered for visibility. All blocking points
+ on the line-of-sight between observer and target will be left untouched.
+ This is different from the <b>target_elev</b> parameter. It can be
+ used to simulate visibility of targets with a certain height.</dd>
+<dt>AZIMUTH1, AZIMUTH2</dt>
+ <dd>These attributes must be given as values between 0.0 and 360.0 and
+ can be used to limit the horizontal scan range of the viewshed.
+ Limits are measured in degrees with 0.0 in the North and values
+ increasing clock-wise. Thus, if you wanted to limit visibility scans
+ to eastward directions only, you would set AZIMUTH1=1 and AZIMUTH2=180.
+ </dd>
+<dt>VERT1, VERT2</dt>
+ <dd>In analogy to 'AZIMUTH1' and 'AZIMUTH2', these attributes can be
+ used to limit the vertical range of the scan. The range is -90 to 0 for
+ 'VERT1' and 0 to 90 for 'VERT2'.
+ At -90 degrees, the observer looks straight down to the ground and at
+ 90 degress up into the sky.
+ </dd>
+<dt>RADIUS1, RADIUS2</dt>
+ <dd>Finally, two attributes exist to limit the distance of the scan.
+ 'RADIUS1' determines the minimum distance from the observer and 'RADIUS2'
+ the maximum distance. The latter is thus equivalent to using the
+ <b>max_dist</b> parameter and exists for compatibility reasons.
+ </dd>
+</dl>
+
+Interpretation of these attributes can be turned off by using the <b>-i</b> flag.
+<p>
+Note that you can use the same control values in operation modes
+other than 'sites', if you specify the equivalent parameters
+(see <A HREF = "#parameters">above</A>). The latter will also work
+as global settings in 'sites' mode. Attributes stored in the vector
+map table will override these global settings unless the <b>-i</b>
+flag is specified. The meaning of <b>OFFSETA, OFFSETB</b> and <b>RADIUS1</b> is
+the same as <b>obs_elev, target_elev</b> and <b>max_dist</b>. The latter still
+exist for compatibility with earlier versions of <em>r.cva</em>.
+
+
+<A NAME = "notes"><H2>NOTES</H2></A>
+
+<EM>r.cva</EM> works in the current geographic region, but ignores the
+current mask. Masking is achieved by specifying a target_mask and/or a
+viewpt_mask. For accurate results, the program must be run with the
+resolution of the geographic region set equal to the resolution of the
+data (see <EM><A HREF="../../html/g.region.html">g.region</A></EM>).
+<P>
+This program assumes all height and distance measurements in meters.
+Legacy units must be converted to meters by the user before running <EM>r.cva</EM>.
+<P>
+It is also assumed, that N-S and E-W resolution of the current region are the
+same. If they differ, distances calculated for limiting the viewshed will
+be imprecise and results affected accordingly.
+<P>
+The number of cells visible from a given viewpoint includes the
+viewpoint itself only if the viewpoint is a target cell.
+Consequently, if neither the viewpoint nor any other cell within
+<EM>max_dist</EM> of it is a target cell then the size of that
+viewpoint's viewshed will be zero. Set the <EM>-m</EM> flag to locate
+such viewpoints (they will show as cells coded `0' as opposed to
+`-1'). Similarly, set the <EM>-m</EM> flag when the <EM>-f</EM> flag has
+been set in order to differentiate cells which are genuinely not
+visible from any viewpoint as opposed to those that are not target
+cells.
+<P>
+<EM>r.cva</EM> always informs the user what sampling frequency was
+attempted and subsequently achieved. This information is printed to
+standard output and is also placed in the output map's history file
+(which can be viewed using <EM><A HREF="../../html/r.info.html">r.info</A></EM>).
+The sample attempted may differ slightly from that requested owing to
+the need to sample a whole number of cells. Any such discrepancy will
+be greatest in the case of small maps of coarse resolution. The
+sample achieved may be further reduced when a <EM>viewpt_mask</EM> has been
+specified. This will not occur with the option `random', because its
+algorithm repeatedly allocates viewpoints until the number that fall
+within mask cells coded `1' produces the frequency attempted. Such a
+reduction can, however, occur with the options `all', `systematic' and
+`sites', because their algorithms allocate viewpoints without
+reference to the mask and then simply ignore those that fall in mask
+cells coded `0'.
+<P>
+
+<H2>SEE ALSO</H2>
+
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>
+<br>
+<EM><A HREF="../../html/r.info.html">r.info</A></EM>
+<br>
+<EM><A HREF="../../html/g.region.html">g.region</A></EM>
+
+<P>
+
+Cumulative viewshed analysis was defined by D. Wheatley, 1995,
+`Cumulative viewshed analysis: a GIS-based method for investigating
+intervisibility, and its archaeological application', in G. R. Lock
+and Z. Stancic (eds.) <EM>Archaeology and Geographic Information
+Systems: A European Perspective</EM>, London: Taylor and Francis,
+pp. 171-186.
+<P>
+Lake et al. (M. W. Lake, P. E. Woodman and S. J. Mithen, 1998,
+`Tailoring GIS software for archaeological applications: an example
+concerning viewshed analysis', <EM>Journal of Archaeological
+Science</EM> 25: 27-38) discuss the first version of <EM>r.cva</EM>
+and provide some timing estimates. They also discuss the importance
+of eliminating the edge effect that occurs when calculating the number
+of cells visible from a viewpoint that is closer than
+<EM>max_dist</EM> to the edge of the map.
+<P>
+Fisher et al. (P. Fisher, C. Farrelly, A. Maddocks and C. Ruggles,
+1997, `Spatial Analysis of Visible Areas from the Bronze Age Cairns of
+Mull', <EM>Journal of Archaeological Science</EM> 24: 581--592)
+provide a useful review of some of the pitfalls that may trap unwary
+users of cumulative viewshed analysis. Most of these are easily
+avoided using <EM>r.cva</EM>.
+
+<H2>AUTHOR</H2>
+
+Mark Lake, Institute of Archaeology, University College London (the author).
+<P>
+With contributions by:
+Benjamin Ducke, Institute of Archaeology, University of Bamberg, Germany.
+
+<H2>ACKNOWLEDGEMENTS</H2>
+
+This version of <EM>r.cva</EM> was completed during the author's
+tenure of a Leverhulme Trust Special Research Fellowship at the
+Institute of Archaeology, University College London, U.K.
+<P>
+The first version of <EM>r.cva</EM> was written by the author as part of the
+MAGICAL Project directed by Dr. Steven Mithen (Dept. of Archaeology,
+University of Reading, U.K.). The MAGICAL Project was made possible
+by a Natural Environment Research Council award (NERC GR3/9540) to
+Dr. Mithen.
+<P>
+<EM>r.cva</EM> draws heavily on the code for <EM><A
+HREF="../../html/r.los.html">r.los</A></EM>, which was written by Kewan
+Q. Khawaja, Intelligent Engineering Systems Laboratory, M.I.T.
+<P>
+Adaptation to GRASS 6 vector and raster formats, additional functionality,
+documentation updates by Benjamin Ducke (2005).
+
+<H2>BUGS</H2>
+
+Missing values in site attributes are currently not handled.
+<p>
+
+Earth curvature correction has not been sufficiently tested. The 'offsetb' option (=target_elev)
+has not been sufficiently tested. Please feel free to test and report!
+<p>
+
+When <EM>max_dist</EM> is set to 100, the N-S and E-W resolutions are
+50, and a <EM>target_mask</EM> has been specified, the cells
+immediately NE, SE, SW and NW of the viewpoint may not be masked when
+then should be. Please report any others that you find to the author.
+<P>
+
+Note that <EM>r.cva</EM> fixes two bugs in the current version of
+<EM><A HREF="../../html/r.los.html">r.los</A></EM>: that <EM>obs_elev</EM> is
+truncated to its integer component and that the target mask
+(<EM>patt_map</EM> in <EM><A HREF="../../html/r.los.html">r.los</A></EM>) does
+not apply to the viewpoint or its eight immediate neighbours.
+<HR>
+<P><a href="index.html">Main index</a> - <a href="raster.html">raster index</a> - <a href="full_index.html">Full index</a></P>
+<P>© 2003-2012 <a href="http://grass.osgeo.org">GRASS Development Team</a></p>
+</body>
+</html>
Added: grass-addons/grass6/raster/r.cva/radians.h
===================================================================
--- grass-addons/grass6/raster/r.cva/radians.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/radians.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,16 @@
+/****************************************************************/
+/* */
+/* radians.h in ~/src/Glos */
+/* */
+/* This header file defines the syntactic sugar to be */
+/* used for some prominent radian measures. */
+/* */
+/****************************************************************/
+
+#define PI 3.141592654
+#define PIBYFOUR 0.785398163
+#define PIBYTWO 1.570796327
+#define TWOPI 6.283185308
+#define THREEPIBYTWO 4.71238898
+
+/****************************************************************/
Added: grass-addons/grass6/raster/r.cva/random_int.c
===================================================================
--- grass-addons/grass6/raster/r.cva/random_int.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/random_int.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,47 @@
+
+/****************************************************************
+ * random_int.c
+ *
+ * Return a random integer between a specified minimum
+ * and maximum value. Integer is a long and is drawn from
+ * a uniform population. Minimum and maximum values may
+ * be returned.
+ *
+ ****************************************************************/
+
+/* Written by Mark Lake on 1/3/96 to:-
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Return a random integer between a specified minimum
+ * and maximum value.
+ *
+ * Called by:-
+ *
+ * 1) random_sample()
+ *
+ * Calls:-
+ *
+ * 1) zufall()
+ */
+
+#include <string.h>
+#include <math.h>
+
+#include <grass/gis.h>
+
+#include "config.h"
+#include "zufall.h"
+
+long int random_int (long int min, long int max)
+{
+ double u[1];
+ char message[32];
+
+ strcpy (message, "Random number generator failed!");
+
+ if (zufall (1, u) != 0)
+ G_fatal_error (message);
+
+ return (long int) min + floor ((max - min + 1) * u[0]);
+}
Added: grass-addons/grass6/raster/r.cva/random_sample.c
===================================================================
--- grass-addons/grass6/raster/r.cva/random_sample.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/random_sample.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,273 @@
+/***********************************************************************/
+/*
+ random_sample.c
+
+ Updated by Mark Lake on 17/8/01 for GRASS 5.x
+
+ CONTAINS
+
+ 1) random_sample */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ random_sample
+
+ Called from:
+
+ 1) main main.c
+
+ Calls:
+
+ 1) line_of_sight line_of_site.c
+ 2) count_visible_cells count_visible_cells.c
+ 3) cumulate_visible_cells cumulate_visible_cells.c
+ 4) delete_lists delete_lists.c */
+
+/***********************************************************************/
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+#include "global_vars.h"
+#include "zufall.h"
+
+#include "init_segment_file.h"
+#include "line_of_sight.h"
+#include "count_visible_cells.h"
+#include "cumulate_visible_cells.h"
+#include "delete_lists.h"
+
+#define MAX_RANDOM_ALLOCATIONS 1000 /* number of attempts to pick cell
+ from which los has not yet been
+ conducted before it is assumed
+ that prog has gone into infinite
+ loop */
+
+void random_sample (int nrows, int ncols, SEGMENT *seg_in_p,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p, SEGMENT *seg_patt_p_v,
+ SEGMENT *seg_rnd_p,
+ double *attempted_sample, double *actual_sample,
+ int replace, int terse, int absolute, RASTER_MAP_TYPE data_type)
+{
+ int row_viewpt, col_viewpt;
+ void *value;
+ double viewpt_elev = 0.0;
+ struct point *heads[16];
+
+ extern long int random_int();
+ int loop_check, duplicate;
+ int row, col;
+ int null_value;
+ CELL mask;
+ CELL cell_no;
+ long int cells_in_map, cells_to_sample;
+ long int cells_sampled, cells_sampled_uniquely;
+ CELL data, already_analysed;
+ char message[128];
+ /* the following are used to store different raster map types */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+
+ /* Seed random number generator */
+ if ((seed < MIN_SEED) || (seed > MAX_SEED))
+ G_fatal_error ("Seed out of range ");
+
+ if (zufalli (seed) != 0) {
+ sprintf (message, "Failed to seed random number generator!");
+ G_fatal_error (message);
+ }
+
+
+ /* Initialise output segment file with appropriate data value */
+ init_segment_file (nrows, ncols, seg_out_2_p, seg_patt_p);
+
+
+ /* Initialise segment file used to prevent / report replacement */
+ data = background;
+ value = (char*) &data;
+ for (row = 0; row < nrows; row ++) {
+ for (col = 0; col < ncols; col ++) {
+ segment_put (seg_rnd_p, value, row, col);
+ }
+ }
+
+
+ /* Initialise variables */
+ already_analysed = background + 1;
+ cells_in_map = nrows * ncols;
+ if (absolute) /* Sample is no. of cells not % */
+ cells_to_sample = (long int) sample;
+ else
+ cells_to_sample = (long int) (sample / 100.0 * (double) cells_in_map);
+ cell_no = 0;
+ cells_sampled = 0;
+ cells_sampled_uniquely = 0;
+
+
+ /* Randomly select cells and perform line of sight analysis from each */
+ if (! terse)
+ fprintf (stdout, "\nPercentage complete:\n");
+
+ for (cells_sampled = 1; cells_sampled <= cells_to_sample; cells_sampled += 1) {
+ cell_no ++;
+
+ /* Do preparatory things to prevent long or infinite loop */
+ loop_check = 0;
+ do {
+ /* Check that number of attempts is not excessive */
+ loop_check ++;
+ if (loop_check > MAX_RANDOM_ALLOCATIONS) {
+ sprintf (message, "Failed to pick a viewpoint which has not\nyet been analysed.");
+ G_fatal_error (message);
+ }
+
+ /* Randomly allocate the viewpoint... */
+ row_viewpt = (int) random_int (0, nrows - 1);
+ col_viewpt = (int) random_int (0, ncols - 1);
+
+ /* check for NULL values in DEM */
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ null_value = 0;
+ segment_get (seg_in_p, value, row_viewpt, col_viewpt);
+
+ if ( data_type == CELL_TYPE ) {
+ viewpt_elev = c_value + obs_elev;
+ if (G_is_c_null_value (&c_value)) {
+ null_value = 1;
+ }
+ }
+ if ( data_type == FCELL_TYPE ) {
+ viewpt_elev = f_value + obs_elev;
+ if (G_is_f_null_value (&f_value)) {
+ null_value = 1;
+ }
+ }
+ if ( data_type == DCELL_TYPE ) {
+ viewpt_elev = d_value + obs_elev;
+ if (G_is_d_null_value (&d_value)) {
+ null_value = 1;
+ }
+ }
+
+ /* ...and check that this cell has not already been analysed */
+ value = (char*) &data;
+ segment_get (seg_rnd_p, value, row_viewpt, col_viewpt);
+ if (replace) {
+ if (data == already_analysed) {
+ data = background; /* Trick to break out out
+ of repeat loop */
+ duplicate = 1;
+ } else
+ duplicate = 0;
+ } else
+ duplicate = 0;
+
+ /* Check whether viewpoint is in area of interest */
+ if(patt_flag_v) {
+ value = (char*) &mask;
+ segment_get(seg_patt_p_v, value, row_viewpt, col_viewpt);
+ } else
+ mask = 1;
+
+ } while ((data == already_analysed) || (! mask) || (null_value) );
+
+
+ /* Don't actually reanalyse duplicates. This saves time.
+ More importantly, it prevents erroneous results
+ with cumulate (-f) option */
+ if (! duplicate) {
+ cells_sampled_uniquely ++;
+
+ /* Mark viewpoint in temporary file to
+ prevent duplication or allow number of
+ duplications to be counted */
+ value = (char*) &already_analysed;
+ segment_put (seg_rnd_p, value, row_viewpt, col_viewpt);
+
+
+ /* Read elevation of current view point */
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ null_value = 0;
+ segment_get (seg_in_p, value, row_viewpt, col_viewpt);
+
+ if ( data_type == CELL_TYPE ) {
+ viewpt_elev = c_value + obs_elev;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ viewpt_elev = f_value + obs_elev;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ viewpt_elev = d_value + obs_elev;
+ }
+
+ if (SPOT_SET) {
+ if (ADJUST) {
+ if (SPOT >= viewpt_elev) {
+ viewpt_elev = SPOT + obs_elev;
+ }
+ } else {
+ viewpt_elev = SPOT;
+ }
+ }
+
+
+ /* Do line of sight analysis from current viewpoint */
+ line_of_sight (cell_no, row_viewpt, col_viewpt, nrows, ncols, viewpt_elev,
+ heads, seg_in_p, seg_out_1_p, seg_patt_p, patt_flag, data_type);
+
+ /* Calculate viewshed */
+ if (! from_cells) {
+ count_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+ } else {
+ cumulate_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+ }
+ /* Free memory used in line-of-sight analysis */
+ delete_lists (heads);
+ }
+
+
+ /* Print percentage complete */
+ if (! terse) {
+ G_percent (cells_sampled, cells_to_sample, 1);
+ }
+
+ } /* END (loop through sampled cells) */
+
+
+ /* Calculate_sample */
+ if (absolute) {
+ *attempted_sample = cell_no;
+ *actual_sample = cells_sampled_uniquely;
+ } else {
+ *attempted_sample = 100.0 / (double) cells_in_map * (double) cell_no;
+ *actual_sample = 100.0 /
+ (double) cells_in_map * (double) cells_sampled_uniquely;
+ }
+}
Added: grass-addons/grass6/raster/r.cva/random_sample.h
===================================================================
--- grass-addons/grass6/raster/r.cva/random_sample.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/random_sample.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,13 @@
+#ifndef __RANDOM_SAMPLE_H__
+#define __RANDOM_SAMPLE_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void random_sample (int nrows, int ncols, SEGMENT *seg_in_p,
+ SEGMENT *seg_out_1_p, SEGMENT *seg_out_2_p,
+ SEGMENT *seg_patt_p, SEGMENT *seg_patt_p_v,
+ SEGMENT *seg_rnd_p,
+ double *attempted_sample, double *actual_sample,
+ int replace, int terse, int absolute, RASTER_MAP_TYPE data_type);
+#endif
Added: grass-addons/grass6/raster/r.cva/scan_all_cells.c
===================================================================
--- grass-addons/grass6/raster/r.cva/scan_all_cells.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/scan_all_cells.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,171 @@
+/***********************************************************************/
+/*
+ scan_all_cells.c
+
+ Updated by Mark Lake on 17/8/01 for GRASS 5.x
+
+ CONTAINS
+
+ 1) scan_all_cells */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ scan_all_cells
+
+ Called from:
+
+ 1) main main.c
+
+ Calls:
+
+ 1) line_of_sight line_of_site.c
+ 2) count_visible_cells count_visible_cells.c
+ 3) cumulate_visible_cells cumulate_visible_cells.c
+ 4) delete_lists delete_lists.c */
+
+/***********************************************************************/
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+#include "global_vars.h"
+
+#include "init_segment_file.h"
+#include "line_of_sight.h"
+#include "count_visible_cells.h"
+#include "cumulate_visible_cells.h"
+#include "delete_lists.h"
+
+void scan_all_cells (int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ int terse, RASTER_MAP_TYPE data_type)
+{
+ int row_viewpt, col_viewpt;
+ int mask;
+ int null_value;
+ int null_count;
+ void *value = NULL;
+ double cells_in_map;
+ CELL cell_no;
+ long int cells_analysed;
+ double viewpt_elev = 0.0;
+ struct point *heads[16];
+ /* the following are used to store different raster map types */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+
+ /* Initialise output segment file with appropriate data value */
+ init_segment_file (nrows, ncols, seg_out_2_p, seg_patt_p);
+
+ /* Initialise variables */
+ cell_no = 0;
+ cells_analysed = 0;
+
+ /* Loop through all cells and perform line of sight analysis from each */
+ if (! terse)
+ fprintf (stdout, "\nPercentage complete:\n");
+
+ cells_in_map = nrows * ncols;
+ null_count = 0;
+
+ for (row_viewpt = 0; row_viewpt < nrows; row_viewpt ++) {
+ for (col_viewpt = 0; col_viewpt < ncols; col_viewpt ++) {
+ cell_no ++;
+
+ if(patt_flag_v) {
+ value = (char*) &mask;
+ segment_get(seg_patt_p_v, value,
+ row_viewpt, col_viewpt);
+ } else
+ mask = 1;
+
+ /* Only include if viewpoint is in area of interest */
+ if (mask) {
+ cells_analysed ++;
+
+ /* Read elevation of current view point */
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ null_value = 0;
+ segment_get (seg_in_p, value, row_viewpt, col_viewpt);
+
+ if ( data_type == CELL_TYPE ) {
+ viewpt_elev = c_value + obs_elev;
+ if (G_is_c_null_value (&c_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+ if ( data_type == FCELL_TYPE ) {
+ viewpt_elev = f_value + obs_elev;
+ if (G_is_f_null_value (&f_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+ if ( data_type == DCELL_TYPE ) {
+ viewpt_elev = d_value + obs_elev;
+ if (G_is_d_null_value (&d_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+
+ if (!null_value) {
+ if (SPOT_SET) {
+ if (ADJUST) {
+ if (SPOT >= viewpt_elev) {
+ viewpt_elev = SPOT + obs_elev;
+ }
+ } else {
+ viewpt_elev = SPOT;
+ }
+ }
+
+
+ line_of_sight (cell_no, row_viewpt, col_viewpt,
+ nrows, ncols, viewpt_elev,
+ heads, seg_in_p, seg_out_1_p,
+ seg_patt_p, patt_flag, data_type);
+
+ if (! from_cells)
+ count_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+ else
+ cumulate_visible_cells (cell_no, row_viewpt, col_viewpt,
+ seg_out_1_p, seg_out_2_p, seg_patt_p, heads);
+
+ delete_lists (heads);
+ }
+ } /* END (process one cell in mask) */
+
+ /* Print progress report */
+ if (! terse) {
+ G_percent (cell_no, (int) cells_in_map, 2);
+ }
+ }
+
+ } /* END (loop through ALL cells) */
+
+ /* Calculate sample */
+ *attempted_sample = 100.0;
+ cells_analysed = cells_analysed - null_count;
+ *actual_sample = 100 / (double) cells_in_map * (double) cells_analysed;
+}
Added: grass-addons/grass6/raster/r.cva/scan_all_cells.h
===================================================================
--- grass-addons/grass6/raster/r.cva/scan_all_cells.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/scan_all_cells.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,15 @@
+#ifndef __SCAN_ALL_CELLS_H__
+#define __SCAN_ALL_CELLS_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+
+void scan_all_cells (int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ int terse, RASTER_MAP_TYPE data_type);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/segment3.c
===================================================================
--- grass-addons/grass6/raster/r.cva/segment3.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/segment3.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,123 @@
+/****************************************************************/
+/* */
+/* segment3.c */
+/* */
+/* This function picks up all the points in one segment */
+/* and performs los analysis on them. */
+/* */
+/****************************************************************/
+
+/* Modified from r.los by MWL on 14/12/95 to:
+ *
+ * Updated by Mark Lake on 17/8/01 for GRASS 5.x
+ *
+ * 1) Return the head of a list of cells in the current
+ * segment that are visible from the current viewpoint.
+ *
+ * Version 3.0, differs from v.1.0 as follows:
+ *
+ * 1) Calls make_list.c version 3;
+ * 2) Calls pts_elim.c version3;
+ * 3) Recieves cell_no and passes it to make_list.c.
+ * and pts_elim.c.
+ *
+ * Implicated files:
+ *
+ * 1) See main3.c.
+ *
+ * Modifications to r.los:
+ *
+ * 1) Recieves cell_no and passes it on to make_list.c
+ * and pts_elim.c.
+ */
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "config.h"
+#include "point.h"
+
+#define NULL_HERE 0
+
+#define NEXT_PT PRESENT_PT->next
+#define NEXT_PT_BACK_PTR PRESENT_PT->next->previous
+#define HEAD_BACK_PTR head->previous
+
+struct point *segment(int segment_no, int xmax, int ymax,
+ double slope_1, double slope_2,
+ int flip, int sign_on_y, int sign_on_x,
+ double viewpt_elev,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_p, SEGMENT *seg_patt_p,
+ int row_viewpt, int col_viewpt, int patt_flag,
+ CELL cell_no, RASTER_MAP_TYPE data_type)
+{
+ int lower_limit_y ,upper_limit_y,less,x,y;
+ int x_actual,y_actual,x_flip,y_flip;
+ struct point *head = NULL_HERE, *PRESENT_PT, *make_list();
+ int quadrant;
+ struct point *hidden_point_elimination();
+
+ /* decide which one of the four quadrants */
+ quadrant = 1 + (segment_no-1)/4;
+
+#ifdef DEBUG
+ printf ("\n Quadrant = %d", quadrant);
+#endif
+
+ if(slope_1 != 0){ less= ymax/slope_1 + 0.99;
+ xmax= (xmax<=less)?xmax:less;
+ }
+
+ /* outer loop over x coors for picking up points */
+
+ for(x=xmax ; x > 0; x--)
+ {
+
+ /* calculate limits for range of y for a single x */
+ lower_limit_y = x * slope_1 + 0.9;
+ upper_limit_y = x * slope_2 ;
+ upper_limit_y = (upper_limit_y<=ymax)?upper_limit_y:ymax;
+
+ /* loop over y range to pick up correct points */
+ for(y= upper_limit_y; y>= lower_limit_y ; y--)
+ {
+ /* calculate actual x, y that lie in current segment */
+
+ if(flip==0) { x_flip = x; y_flip = y;}
+ else { y_flip = x; x_flip = y;}
+
+ x_actual = sign_on_x*x_flip ;
+ y_actual = sign_on_y*y_flip;
+
+
+ /* add chosen point to the point list */
+ head = make_list(head,y_actual,x_actual,seg_in_p,seg_out_p,
+ viewpt_elev,quadrant,row_viewpt,col_viewpt,cell_no, data_type);
+ }
+ } /* end of outer loop */
+
+ if(head != NULL_HERE)
+ {
+ /* assign back pointers in linked list */
+ HEAD_BACK_PTR = NULL_HERE;
+ PRESENT_PT = head;
+ while(NEXT_PT != NULL_HERE)
+ {
+ NEXT_PT_BACK_PTR = PRESENT_PT;
+ PRESENT_PT = NEXT_PT;
+ }
+ head = hidden_point_elimination(head,viewpt_elev,
+ seg_in_p,seg_out_p,seg_patt_p,quadrant,sign_on_y,
+ sign_on_x,row_viewpt,col_viewpt,patt_flag,cell_no,
+ data_type);
+ }
+
+#ifdef DEBUG
+ printf ("\nSegment, head = %p", head); fflush (stdout);
+#endif
+
+ return(head);
+
+}
+
+/**************** END OF FUNCTION "SEGMENT" *********************/
Added: grass-addons/grass6/raster/r.cva/segment3.h
===================================================================
--- grass-addons/grass6/raster/r.cva/segment3.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/segment3.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,12 @@
+#ifndef __SEGMENT3_H__
+#define __SEGMENT3_H__
+
+struct point *segment(int segment_no, int xmax, int ymax,
+ double slope_1, double slope_2,
+ int flip, int sign_on_y, int sign_on_x,
+ double viewpt_elev,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_p, SEGMENT *seg_patt_p,
+ int row_viewpt, int col_viewpt, int patt_flag,
+ CELL cell_no, RASTER_MAP_TYPE data_type);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/segment_file.c
===================================================================
--- grass-addons/grass6/raster/r.cva/segment_file.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/segment_file.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,332 @@
+/* modified by Benjamin Ducke to work with different raster map types */
+/* May 2004 */
+
+/**************************************************************************/
+/*
+ Use fract to set the number of segments. On a machine with 32Kb RAM
+ it is possible that only 1 segment is required unless many maps
+ are to be opened simultaneously. Increasing the number of segments
+ seriously degrades performance */
+
+/**************************************************************************/
+
+#include <string.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <fcntl.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+char* Create_segment_file (int, int, RASTER_MAP_TYPE);
+int Init_segment_file (char*, SEGMENT*);
+
+
+/**************************************************************************/
+/* Public functions */
+/**************************************************************************/
+
+void Close_segmented_infile (int map_fd, char* seg_name, int seg_fd, SEGMENT* seg)
+{
+ segment_release (seg);
+ close (seg_fd);
+ unlink (seg_name);
+ G_close_cell (map_fd);
+}
+
+
+/* Create result map in CELL format */
+/* NULL cell handling added by Benjamin Ducke, May 2004 */
+void Close_segmented_outfile (char* map_name, int map_fd,
+ char* seg_name, int seg_fd, SEGMENT* seg,
+ CELL* cell_buf, int terse,
+ SEGMENT* elevation_seg, RASTER_MAP_TYPE data_type, int make_nulls)
+{
+ unsigned long row, nrows, col, ncols;
+ void *value = NULL;
+ char *null_flags;
+ /* the following are used to store different raster map types */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+
+ /* Find number of rows and columns in elevation map */
+ nrows = G_window_rows();
+ ncols = G_window_cols();
+
+ /* allocate memory for a NULL data row */
+ null_flags = G_calloc ((signed) ncols, sizeof (char));
+
+ /* Write pending updates by segment_put() to output map */
+ segment_flush(seg);
+
+ if (!terse) {
+ fprintf (stdout, "\nWriting raster map '%s': \n", map_name);
+ }
+
+ /* Convert output submatrices to full cell overlay */
+ for(row=0; row< nrows; row++) {
+ segment_get_row(seg, cell_buf, (signed) row);
+ /* check for NULL values in the corresponding row of the */
+ /* elevation input file. If there are any, set the CELLs */
+ /* in the output file to NULL, as well */
+
+ /* initialize null data row */
+ for (col=0; col<ncols; col++) {
+ null_flags[col] = 0;
+ }
+
+ /* convert all -1 cells to NULL, if user wants it so */
+ if ( make_nulls == 1 ) {
+ for (col=0; col<ncols; col++) {
+ if ( cell_buf[col] == -1 ) {
+ null_flags[col] = 1;
+ }
+ }
+ }
+
+ /* update NULL flags in row */
+ for (col=0; col<ncols; col++) {
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ segment_get (elevation_seg, value, (signed) row, (signed) col);
+
+ /* check for NULL value and record in null_flags row */
+ if ( data_type == CELL_TYPE ) {
+ if (G_is_c_null_value (&c_value) ) {
+ null_flags[col] = 1; /* signal NULL value in null_flags row */
+ }
+ }
+ if ( data_type == FCELL_TYPE ) {
+ if (G_is_f_null_value (&f_value) ) {
+ null_flags[col] = 1; /* signal NULL value in null_flags row */
+ }
+ }
+ if ( data_type == DCELL_TYPE ) {
+ if (G_is_d_null_value (&d_value) ) {
+ null_flags[col] = 1; /* signal NULL value in null_flags row */
+ }
+ }
+ }
+ /* set all NULL cells according to null_flags row */
+ G_insert_c_null_values (cell_buf, null_flags, (signed) ncols);
+
+ /* now, write this row to disk */
+ if(G_put_raster_row (map_fd, cell_buf, CELL_TYPE) < 0) {
+ G_fatal_error ("Failed to convert output submatrices ");
+ }
+
+ /* progress display */
+ if (! terse) {
+ G_percent ((signed) row, (signed) nrows-1,2);
+ }
+ }
+
+ G_free (null_flags);
+
+/* Close files */
+ segment_release (seg);
+ close (seg_fd);
+ unlink (seg_name);
+ G_close_cell (map_fd);
+}
+
+void Close_segmented_tmpfile (char *seg_name, int seg_fd, SEGMENT *seg)
+{
+ /* Close files */
+ segment_release (seg);
+ close (seg_fd);
+ unlink (seg_name);
+}
+
+/* this routine creates four segment files on disk, then stores */
+/* a raster map in them */
+/* store a raster map in a mem segment on disk */
+void Segment_infile (char* map_name, char* mapset, char* search_mapset, int* map_fd,
+ char* seg_name, int* seg_fd, SEGMENT* seg, void *cell_buf, int fract,
+ RASTER_MAP_TYPE data_type)
+{
+ int nrows, row;
+
+ strcpy (mapset, G_find_cell (map_name, search_mapset));
+ if (mapset == NULL)
+ G_fatal_error ("Can't find displayed data layer ");
+ if ((*map_fd = G_open_cell_old (map_name, mapset)) < 0)
+ G_fatal_error ("Can't open displayed data layer ");
+
+ nrows = G_window_rows ();
+
+ /* next line creates the actual segment, returns name of disk file */
+ strcpy (seg_name, Create_segment_file (nrows, fract, data_type));
+ *seg_fd = open (seg_name, O_RDWR);
+
+ if (*seg_fd < 0)
+ G_fatal_error ("Cannot open segment disk file ");
+
+ segment_init (seg, *seg_fd, 4);
+
+ /* store map data in segment file IGNORE MASK */
+ /* cell_buf must be pre-allocated by caller! */
+ for (row = 0; row < nrows; row++) {
+ if (G_get_raster_row_nomask (*map_fd, cell_buf, row, data_type) < 0) {
+ G_fatal_error ("Unable to read data layer into segment file ");
+ }
+ segment_put_row (seg, cell_buf, row); /* store raster map in segment file */
+ }
+
+}
+
+
+/**************************************************************************/
+
+void Segment_outfile (char* map_name, char* mapset, int* map_fd,
+ char* seg_name, int* seg_fd, SEGMENT* seg, char* prompt, int fract,
+ RASTER_MAP_TYPE data_type)
+{
+ int nrows;
+ char* mapset_address;
+
+ mapset_address = G_ask_cell_new (prompt, map_name);
+
+ if (mapset_address == NULL)
+ G_fatal_error ("Can't continue without filename ");
+
+ strcpy (mapset, mapset_address);
+
+ if ((*map_fd = G_open_cell_new (map_name)) < 0)
+ G_fatal_error ("Can't create output layer ");
+
+ nrows = G_window_rows ();
+ strcpy (seg_name, Create_segment_file (nrows, fract, data_type));
+ *seg_fd = open (seg_name, 2);
+
+ if (*seg_fd < 0)
+ G_fatal_error ("Can't open segment disk file ");
+
+ segment_init (seg, *seg_fd, 4);
+}
+
+
+/**************************************************************************/
+
+void Segment_named_outfile (char* map_name, char* mapset, int* map_fd,
+ char* seg_name, int* seg_fd, SEGMENT* seg, int overwrite,
+ int terse, int fract, RASTER_MAP_TYPE data_type)
+{
+ int nrows;
+ char* mapset_address;
+ char message [64];
+
+ mapset_address = G_find_cell (map_name, G_mapset ());
+ strcpy (mapset, G_mapset ());
+ if (mapset_address != NULL)
+ {
+ if (!overwrite)
+ {
+ sprintf (message, "Raster map '%s' exists ", map_name);
+ G_fatal_error (message);
+ }
+ else
+ {
+ if (!terse)
+ fprintf (stdout, "\nOverwriting raster map '%s' ", map_name);
+ }
+ }
+ else
+ {
+ if (!terse)
+ fprintf (stdout, "\nCreating raster map '%s' ", map_name);
+ }
+
+ if ((*map_fd = G_open_cell_new (map_name)) < 0)
+ G_fatal_error ("Can't create output layer ");
+
+ nrows = G_window_rows ();
+
+ strcpy (seg_name, Create_segment_file (nrows, fract, data_type));
+
+ *seg_fd = open (seg_name, O_RDWR);
+
+ if (*seg_fd < 0)
+ G_fatal_error ("Can't open segment disk file ");
+
+ segment_init (seg, *seg_fd, 4);
+}
+
+
+void Segment_tmpfile (char* seg_name, int* seg_fd, SEGMENT* seg, int fract,
+ RASTER_MAP_TYPE data_type)
+{
+ int nrows;
+
+ nrows = G_window_rows ();
+
+ strcpy (seg_name, Create_segment_file (nrows, fract, data_type));
+ *seg_fd = open (seg_name, O_RDWR);
+ if (*seg_fd < 0)
+ G_fatal_error ("Can't open segment disk file ");
+
+ /* get a handle to read segment file data values */
+ /* dimensions of segment and size of stored items are read */
+ /* automatically from the file header and stored in 'seg' */
+ segment_init (seg, *seg_fd, 4);
+}
+
+
+/**************************************************************************/
+/* Private functions */
+/**************************************************************************/
+
+char* Create_segment_file (int nrows, int fract, RASTER_MAP_TYPE data_type)
+{
+ int seg_fd;
+ char* seg_name;
+ int ncols;
+ int submatrix_rows, submatrix_cols;
+ int length_data_item = 0;
+ int seg_error;
+
+ ncols = G_window_cols ();
+ if (data_type == CELL_TYPE) length_data_item = sizeof (CELL);
+ if (data_type == FCELL_TYPE) length_data_item = sizeof (FCELL);
+ if (data_type == DCELL_TYPE) length_data_item = sizeof (DCELL);
+ submatrix_rows = nrows/fract + 1;
+ submatrix_cols = ncols/fract + 1;
+ seg_name = G_tempfile ();
+ seg_fd = creat (seg_name, 0666);
+
+ /* create a new segment file on disk */
+ /* size of data items (CELL, DCELL, FCELL ...) is stored in the file header */
+ seg_error = segment_format (seg_fd, nrows, ncols,
+ submatrix_rows, submatrix_cols, length_data_item);
+
+ if ( seg_error != 1 ) {
+ G_fatal_error ("Cannot create segment file(s).\nCheck available diskspace and access rights.");
+ }
+ close (seg_fd);
+
+ return seg_name;
+}
+
+
+/**************************************************************************/
+
+int Init_segment_file (char* seg_name, SEGMENT* seg)
+{
+ int seg_fd;
+
+ seg_fd = open (seg_name, 2);
+ if (seg_fd < 0)
+ G_fatal_error ("Cannot open segment disk file");
+ segment_init (seg, seg_fd, 4);
+ return seg_fd;
+}
Added: grass-addons/grass6/raster/r.cva/segment_file.h
===================================================================
--- grass-addons/grass6/raster/r.cva/segment_file.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/segment_file.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,27 @@
+#ifndef __SEGMENT_FILE_H__
+#define __SEGMENT_FILE_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void Segment_infile (char* map_name, char* mapset, char* search_mapset, int* map_fd,
+ char* seg_name, int* seg_fd, SEGMENT* seg, void *cell_buf, int fract,
+ RASTER_MAP_TYPE data_type);
+
+void Segment_named_outfile (char* map_name, char* mapset, int* map_fd,
+ char* seg_name, int* seg_fd, SEGMENT* seg, int overwrite,
+ int terse, int fract, RASTER_MAP_TYPE data_type);
+
+void Segment_tmpfile (char* seg_name, int* seg_fd, SEGMENT* seg, int fract,
+ RASTER_MAP_TYPE data_type);
+
+void Close_segmented_tmpfile (char *seg_name, int seg_fd, SEGMENT *seg);
+
+void Close_segmented_infile (int map_fd, char* seg_name, int seg_fd, SEGMENT* seg);
+
+void Close_segmented_outfile (char* map_name, int map_fd,
+ char* seg_name, int seg_fd, SEGMENT* seg,
+ CELL* cell_buf, int terse,
+ SEGMENT* elevation_seg, RASTER_MAP_TYPE data_type, int make_nulls);
+
+#endif
Added: grass-addons/grass6/raster/r.cva/systematic_sample.c
===================================================================
--- grass-addons/grass6/raster/r.cva/systematic_sample.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/systematic_sample.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,216 @@
+/***********************************************************************/
+/*
+ systematic_sample.c
+
+ Updated by Mark Lake on 17/8/01 for GRASS 5.x
+
+ CONTAINS
+
+ 1) systematic_sample */
+
+/***********************************************************************/
+
+/***********************************************************************/
+/*
+ systematic_sample
+
+ Called from:
+
+ 1) main main.c
+
+ Calls:
+
+ 1) line_of_sight line_of_site.c
+ 2) count_visible_cells count_visible_cells.c
+ 3) cumulate_visible_cells cumulate_visible_cells.c
+ 4) delete_lists delete_lists.c */
+
+/***********************************************************************/
+
+#include "config.h"
+#include <math.h>
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+#include "point.h"
+#include "global_vars.h"
+
+#include "init_segment_file.h"
+#include "line_of_sight.h"
+#include "count_visible_cells.h"
+#include "cumulate_visible_cells.h"
+#include "delete_lists.h"
+
+void systematic_sample (int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ int terse,
+ int absolute,
+ RASTER_MAP_TYPE data_type)
+{
+ int row_viewpt, col_viewpt;
+ int mask;
+ int null_value;
+ int null_count;
+ void *value = NULL;
+ CELL cell_no, actual_cells;
+ double viewpt_elev = 0.0;
+ struct point *heads[16];
+
+ double row_increment, col_increment;
+ double no_cells_to_sample, actual_cells_to_sample, cells_in_map;
+ double side_of_square_sample;
+ double ratio_of_sides, rows_in_sample, cols_in_sample;
+ double sample_row, sample_col;
+
+ /* the following are used to store different raster map types */
+ CELL c_value;
+ FCELL f_value;
+ DCELL d_value;
+
+
+
+ /* Initialise output segment file with appropriate data value */
+ init_segment_file (nrows, ncols, seg_out_2_p, seg_patt_p);
+
+ /* Calculate the row and col spacings */
+ if (absolute) /* Sample is no. of cells not % */
+ no_cells_to_sample = (long int) sample;
+ else
+ no_cells_to_sample = sample / 100.0 * (double) nrows * (double) ncols;
+ side_of_square_sample = sqrt (no_cells_to_sample);
+ ratio_of_sides = (double) nrows / (double) ncols;
+ rows_in_sample = ratio_of_sides * side_of_square_sample;
+ if (fmod (rows_in_sample, (double) 1) < 0.5)
+ rows_in_sample = floor (rows_in_sample);
+ else
+ rows_in_sample = ceil (rows_in_sample);
+ cols_in_sample = (1 / ratio_of_sides) * side_of_square_sample;
+ if (fmod (cols_in_sample, (double) 1) < 0.5)
+ cols_in_sample = floor (cols_in_sample);
+ else
+ cols_in_sample = ceil (cols_in_sample);
+ row_increment = nrows / rows_in_sample;
+ col_increment = ncols / cols_in_sample;
+ actual_cells_to_sample = rows_in_sample * cols_in_sample;
+ cells_in_map = nrows * ncols;
+ null_count = 0;
+
+
+ /* Loop through sample of cells and perform
+ line of sight analysis from each */
+
+ if (! terse)
+ fprintf (stdout, "\nPercentage complete:\n");
+
+ actual_cells = 0;
+ cell_no = 0;
+
+ sample_row = ((double) nrows - ((rows_in_sample - 1) * row_increment)) / 2;
+ while (sample_row < nrows) {
+ sample_col = ((double) ncols - ((cols_in_sample - 1) * col_increment)) / 2;
+
+ while (sample_col < ncols) {
+ row_viewpt = floor (sample_row);
+ col_viewpt = floor (sample_col);
+
+ cell_no ++;
+
+ if (patt_flag_v) {
+ value = (char*) &mask;
+ segment_get(seg_patt_p_v, value,
+ row_viewpt, col_viewpt);
+ } else
+ mask = 1;
+
+ /* Only include if viewpoint is in area of interest */
+ if (mask) {
+ actual_cells ++;
+
+ /* Read elevation of current view point */
+ /* value = (char*) &read_viewpt_elev; */
+ if ( data_type == CELL_TYPE ) {
+ value = (CELL*) &c_value;
+ }
+ if ( data_type == FCELL_TYPE ) {
+ value = (FCELL*) &f_value;
+ }
+ if ( data_type == DCELL_TYPE ) {
+ value = (DCELL*) &d_value;
+ }
+
+ null_value = 0;
+ segment_get (seg_in_p, value, row_viewpt, col_viewpt);
+
+ if ( data_type == CELL_TYPE ) {
+ viewpt_elev = c_value + obs_elev;
+ if (G_is_c_null_value (&c_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+ if ( data_type == FCELL_TYPE ) {
+ viewpt_elev = f_value + obs_elev;
+ if (G_is_f_null_value (&f_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+ if ( data_type == DCELL_TYPE ) {
+ viewpt_elev = d_value + obs_elev;
+ if (G_is_d_null_value (&d_value)) {
+ null_value = 1;
+ null_count ++;
+ }
+ }
+
+ if (!null_value) {
+ if (SPOT_SET) {
+ if (ADJUST) {
+ if (SPOT >= viewpt_elev) {
+ viewpt_elev = SPOT + obs_elev;
+ }
+ } else {
+ viewpt_elev = SPOT;
+ }
+ }
+
+ line_of_sight (cell_no, row_viewpt, col_viewpt,
+ nrows, ncols, viewpt_elev,
+ heads, seg_in_p, seg_out_1_p,
+ seg_patt_p, patt_flag, data_type);
+
+ if (! from_cells)
+ count_visible_cells (cell_no, row_viewpt, col_viewpt, seg_out_1_p, seg_out_2_p,
+ seg_patt_p, heads);
+ else
+ cumulate_visible_cells (cell_no, row_viewpt, col_viewpt, seg_out_1_p, seg_out_2_p,
+ seg_patt_p, heads);
+ delete_lists (heads);
+ }
+
+ /* Print progress report */
+ if (! terse) {
+ G_percent (actual_cells, (int) actual_cells_to_sample, 1);
+ }
+ }
+ sample_col += col_increment;
+ } /* END (while (sample_col < ncols) */
+
+ sample_row += row_increment;
+ } /* END ( while (sample_row < nrows) ) */
+
+
+ /* Calculate sample */
+ actual_cells = actual_cells - null_count;
+ if (absolute) {
+ *attempted_sample = actual_cells_to_sample;
+ *actual_sample = actual_cells;
+ } else {
+ *attempted_sample = 100.0 / cells_in_map * actual_cells_to_sample;
+ *actual_sample = 100.0 / cells_in_map * actual_cells;
+ }
+}
Added: grass-addons/grass6/raster/r.cva/systematic_sample.h
===================================================================
--- grass-addons/grass6/raster/r.cva/systematic_sample.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/systematic_sample.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,17 @@
+#ifndef __SYSTEMATIC_SAMPLE_H__
+#define __SYSTEMATIC_SAMPLE_H__
+
+#include <grass/gis.h>
+#include <grass/segment.h>
+
+void systematic_sample (int nrows, int ncols,
+ SEGMENT *seg_in_p, SEGMENT *seg_out_1_p,
+ SEGMENT *seg_out_2_p, SEGMENT *seg_patt_p,
+ SEGMENT *seg_patt_p_v,
+ double *attempted_sample, double *actual_sample,
+ int terse,
+ int absolute,
+ RASTER_MAP_TYPE data_type);
+
+
+#endif
Added: grass-addons/grass6/raster/r.cva/version
===================================================================
--- grass-addons/grass6/raster/r.cva/version (rev 0)
+++ grass-addons/grass6/raster/r.cva/version 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1 @@
+0.9.5
Added: grass-addons/grass6/raster/r.cva/zufall.c
===================================================================
--- grass-addons/grass6/raster/r.cva/zufall.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/zufall.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,183 @@
+#include<stdio.h>
+#include<math.h>
+#include"zufall.h"
+
+/*-
+ * portable lagged Fibonacci series uniform random number generator with
+ * "lags" -273 und -607:
+ *
+ * t = u(i-273)+buff(i-607) (floating pt.)
+ * u(i) = t - float(int(t))
+ *
+ * W.P. Petersen, IPS, ETH Zuerich, 19 Mar. 92
+ */
+
+struct klotz0 klotz0_1;
+struct klotz1 klotz1_1;
+
+int zufall (int n, double *a)
+{
+ static int buffsz = 607;
+ int left, aptr, bptr, aptr0, i, k, q, nn, vl, qq, k273, k607;
+ double t;
+ extern struct klotz0 klotz0_1;
+
+ --a;
+ aptr = 0;
+ nn = n;
+
+L1:
+ if (nn <= 0)
+ return 0;
+
+ q = (nn - 1) / 607; /* factor nn = q*607 + r */
+ left = buffsz - klotz0_1.ptr;
+
+ if (q <= 1)
+ { /* only one or fewer full segments */
+ if (nn < left)
+ {
+ for (i = 1; i <= nn; ++i)
+ a[i + aptr] = klotz0_1.buff[klotz0_1.ptr + i - 1];
+ klotz0_1.ptr += nn;
+ return 0;
+ }
+ else
+ {
+ for (i = 1; i <= left; ++i)
+ a[i + aptr] = klotz0_1.buff[klotz0_1.ptr + i - 1];
+ klotz0_1.ptr = 0;
+ aptr += left;
+ nn -= left;
+ /* buff -> buff case */
+ vl = 273;
+ k273 = 334;
+ k607 = 0;
+ for (k = 1; k <= 3; ++k)
+ {
+ /* dir$ ivdep */
+ /* vdir nodep */
+ /* VOCL LOOP, TEMP(t), NOVREC(buff) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = klotz0_1.buff[k273 + i - 1] + klotz0_1.buff[k607 + i - 1];
+ klotz0_1.buff[k607 + i - 1] = t - (float) ((int) t);
+ }
+ k607 += vl;
+ k273 += vl;
+ vl = 167;
+ if (k == 1)
+ {
+ k273 = 0;
+ }
+ }
+ goto L1;
+ }
+ }
+ else
+ {
+ /* more than 1 full segment */
+ for (i = 1; i <= left; ++i)
+ a[i + aptr] = klotz0_1.buff[klotz0_1.ptr + i - 1];
+ nn -= left;
+ klotz0_1.ptr = 0;
+ aptr += left;
+
+ /* buff -> a(aptr0) */
+
+ vl = 273;
+ k273 = 334;
+ k607 = 0;
+ for (k = 1; k <= 3; ++k)
+ {
+ if (k == 1)
+ {
+ /* VOCL LOOP, TEMP(t) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = klotz0_1.buff[k273 + i - 1] + klotz0_1.buff[k607 + i - 1];
+ a[aptr + i] = t - (float) ((int) t);
+ }
+ k273 = aptr;
+ k607 += vl;
+ aptr += vl;
+ vl = 167;
+ }
+ else
+ {
+ /* dir$ ivdep */
+ /* vdir nodep */
+ /* VOCL LOOP, TEMP(t) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = a[k273 + i] + klotz0_1.buff[k607 + i - 1];
+ a[aptr + i] = t - (float) ((int) t);
+ }
+ k607 += vl;
+ k273 += vl;
+ aptr += vl;
+ }
+ }
+ nn += -607;
+
+ /* a(aptr-607) -> a(aptr) for last of the q-1 segments */
+
+ aptr0 = aptr - 607;
+ vl = 607;
+
+ /* vdir novector */
+ for (qq = 1; qq <= q - 2; ++qq)
+ {
+ k273 = aptr0 + 334;
+ /* dir$ ivdep */
+ /* vdir nodep */
+ /* VOCL LOOP, TEMP(t), NOVREC(a) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = a[k273 + i] + a[aptr0 + i];
+ a[aptr + i] = t - (float) ((int) t);
+ }
+ nn += -607;
+ aptr += vl;
+ aptr0 += vl;
+ }
+
+ /* a(aptr0) -> buff, last segment before residual */
+
+ vl = 273;
+ k273 = aptr0 + 334;
+ k607 = aptr0;
+ bptr = 0;
+ for (k = 1; k <= 3; ++k)
+ {
+ if (k == 1)
+ {
+ /* VOCL LOOP, TEMP(t) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = a[k273 + i] + a[k607 + i];
+ klotz0_1.buff[bptr + i - 1] = t - (double) ((int) t);
+ }
+ k273 = 0;
+ k607 += vl;
+ bptr += vl;
+ vl = 167;
+ }
+ else
+ {
+ /* dir$ ivdep */
+ /* vdir nodep */
+ /* VOCL LOOP, TEMP(t), NOVREC(buff) */
+ for (i = 1; i <= vl; ++i)
+ {
+ t = klotz0_1.buff[k273 + i - 1] + a[k607 + i];
+ klotz0_1.buff[bptr + i - 1] = t - (float) ((int) t);
+ }
+ k607 += vl;
+ k273 += vl;
+ bptr += vl;
+ }
+ }
+ goto L1;
+ }
+} /* zufall */
Added: grass-addons/grass6/raster/r.cva/zufall.h
===================================================================
--- grass-addons/grass6/raster/r.cva/zufall.h (rev 0)
+++ grass-addons/grass6/raster/r.cva/zufall.h 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,31 @@
+#include "config.h"
+
+#define min(a,b) ((a) <= (b) ? (a) : (b))
+
+/* clock cycle for machine */
+/* #define CYCLE 2.9e-9 cycle for SX-3: */
+/* #define CYCLE 6.0e-9 cycle for Y-MP: */
+/* #define CYCLE 3.2e-9 cycle for VP2200 */
+
+# define MIN_SEED 0
+#define MAX_SEED 31328
+
+struct klotz0 {
+ double buff[607];
+ int ptr;
+};
+
+struct klotz1 {
+ double xbuff[1024];
+ int first, xptr;
+};
+
+int normal00();
+int normalen();
+int normalrs();
+int normalsv();
+int fische();
+int zufall();
+int zufalli();
+int zufallrs();
+int zufallsv();
Added: grass-addons/grass6/raster/r.cva/zufalli.c
===================================================================
--- grass-addons/grass6/raster/r.cva/zufalli.c (rev 0)
+++ grass-addons/grass6/raster/r.cva/zufalli.c 2012-03-31 07:58:03 UTC (rev 51216)
@@ -0,0 +1,46 @@
+#include "zufall.h"
+int zufalli(int seed)
+{
+ static int kl = 9373;
+ static int ij = 1802;
+
+ static int i, j, k, l, m;
+ static double s, t;
+ static int ii, jj;
+ extern struct klotz0 klotz0_1;
+
+/* generates initial seed buffer by linear congruential */
+/* method. Taken from Marsaglia, FSU report FSU-SCRI-87-50 */
+/* variable seed should be 0 < seed <31328 */
+
+ if (seed != 0) {
+ ij = seed;
+ }
+
+ i = ij / 177 % 177 + 2;
+ j = ij % 177 + 2;
+ k = kl / 169 % 178 + 1;
+ l = kl % 169;
+
+ for (ii = 1; ii <= 607; ++ii) {
+ s = (float)0.;
+ t = (float).5;
+ for (jj = 1; jj <= 24; ++jj) {
+ m = i * j % 179 * k % 179;
+ i = j;
+ j = k;
+ k = m;
+ l = (l * 53 + 1) % 169;
+ if (l * m % 64 >= 32) {
+ s += t;
+ }
+ t *= (float).5;
+/* L2: */
+ }
+
+ /* printf("DIAG: s %g\n",s); */
+ klotz0_1.buff[ii - 1] = s;
+/* L1: */
+ }
+ return 0;
+} /* zufalli_ */
More information about the grass-commit
mailing list