Add GDALViewshedGenerate binding

errors.go

@@ -154,6 +154,7 @@ func ErrLogger(fn ErrorHandler) interface {
 	GridOption
 	NearblackOption
 	DemOption
+	ViewshedOption
 	SetGCPsOption
 	GCPsToGeoTransformOption
 	RegisterPluginOption
@@ -388,6 +389,9 @@ func (ec errorCallback) setNearblackOpt(o *nearBlackOpts) {
 func (ec errorCallback) setDemOpt(o *demOpts) {
 	o.errorHandler = ec.fn
 }
+func (ec errorCallback) setViewshedOpt(o *viewshedOpts) {
+	o.errorHandler = ec.fn
+}
 func (ec errorCallback) setSetGCPsOpt(o *setGCPsOpts) {
 	o.errorHandler = ec.fn
 }

godal.cpp

@@ -1905,6 +1905,22 @@ GDALDatasetH godalDem(cctx *ctx, const char *pszDest, const char *pszProcessing,
 	return ret;
 }
 
+GDALDatasetH godalViewshedGenerate(cctx *ctx, GDALRasterBandH bnd, const char *pszDriverName, const char *pszTargetRasterName, const char **papszCreationOptions, double dfObserverX,
+    					   double dfObserverY, double dfObserverHeight, double dfTargetHeight, double dfVisibleVal, double dfInvisibleVal, double dfOutOfRangeVal, 
+						   double dfNoDataVal, double dfCurvCoeff, GUInt32 eMode, double dfMaxDistance, GUInt32 heightMode) {
+	godalWrap(ctx);
+	GDALDatasetH ret = nullptr;
+#if GDAL_VERSION_NUM >= GDAL_COMPUTE_VERSION(3, 1, 0)
+	ret = GDALViewshedGenerate(bnd, pszDriverName, pszTargetRasterName, papszCreationOptions, dfObserverX, dfObserverY, dfObserverHeight, dfTargetHeight, dfVisibleVal, dfInvisibleVal, dfOutOfRangeVal, 
+	dfNoDataVal, dfCurvCoeff, GDALViewshedMode(eMode), dfMaxDistance, nullptr, nullptr, GDALViewshedOutputType(heightMode), nullptr);
+#endif
+	if(ret == nullptr) {
+		forceError(ctx);
+	}
+	godalUnwrap();
+	return ret;
+}
+
 OGRSpatialReferenceH godalGetGCPSpatialRef(GDALDatasetH hSrcDS) {
 	return GDALGetGCPSpatialRef(hSrcDS);
 }

godal.go

@@ -4060,6 +4060,88 @@ func (ds *Dataset) Dem(destPath, processingMode string, colorFilename string, sw
 	return &Dataset{majorObject{C.GDALMajorObjectH(dsRet)}}, nil
 }
 
+// ViewshedMode is the "cell height calculation mode" for the viewshed process
+//
+// Source: https://github.com/OSGeo/gdal/blob/master/alg/viewshed/viewshed_types.h
+type ViewshedMode uint32
+
+const (
+	// MDiagonal is the "diagonal mode"
+	MDiagonal ViewshedMode = iota + 1
+	// MEdge is the "edge mode"
+	MEdge
+	// MMax is the "maximum value produced by Diagonal and Edge mode"
+	MMax
+	// MMin is the "minimum value produced by Diagonal and Edge mode"
+	MMin
+)
+
+// ViewshedOutputType sets the return type and information represented by the returned data
+//
+// Source: https://gdal.org/en/stable/programs/gdal_viewshed.html
+//
+// NOTE: "Cumulative (ACCUM)" mode not currently supported, as it's not available in the `GDALViewshedGenerate` function
+// (it's only used in the command line invocation of `viewshed`)
+type ViewshedOutputType uint32
+
+const (
+	// Normal returns a raster of type Byte containing visible locations
+	Normal ViewshedOutputType = iota + 1
+	// MinTargetHeightFromDem return a raster of type Float64 containing the minimum target height for target to be visible from the DEM surface
+	MinTargetHeightFromDem
+	// MinTargetHeightFromGround return a raster of type Float64 containing the minimum target height for target to be visible from ground level
+	MinTargetHeightFromGround
+)
+
+// Viewshed (binding for GDALViewshedGenerate), creates a viewshed from a raster DEM, these parameters (mostly) map to to parameters for GDALViewshedGenerate
+//
+// for more information see: https://gdal.org/en/stable/api/gdal_alg.html#_CPPv420GDALViewshedGenerate15GDALRasterBandHPKcPKc12CSLConstListddddddddd16GDALViewshedModed16GDALProgressFuncPv22GDALViewshedOutputType12CSLConstList
+//
+// Creations options can be set through options with:
+//
+//	Viewshed(bnd, "mem", "none", ..., CreationOption("TILED=YES","BLOCKXSIZE=256"))
+func Viewshed(targetBand Band, driverName *DriverName, targetRasterName string, observerX float64, observerY float64, observerHeight float64, targetHeight float64,
+	visibleVal float64, invisibleVal float64, outOfRangeVal float64, noDataVal float64, curveCoeff float64, mode ViewshedMode, maxDistance float64,
+	heightMode ViewshedOutputType, opts ...ViewshedOption) (*Dataset, error) {
+
+	// TODO: Should I put a 'warning' in the documentation for `Viewshed` instead of disallowing the last two configurations?
+	if !CheckMinVersion(3, 1, 0) {
+		return nil, errors.New("failed to run, 'viewshed' not supported on GDAL versions < 3.1.0")
+	} else if !CheckMinVersion(3, 4, 2) {
+		return nil, errors.New("cannot run 'viewshed' with GDAL version <= 3.4.1, as some tests produce invalid results under these conditions")
+	} else if !CheckMinVersion(3, 10, 0) && heightMode == MinTargetHeightFromDem {
+		return nil, errors.New("height mode CANNOT be `MinTargetHeightFromDem` when running a GDAL version < 3.10, as some tests produce invalid results under these conditions")
+	}
+
+	// Allow `driverName` to be null and handle it here to match parameter/behaviour of GDALViewshedGenerate
+	defaultDriverName := GTiff
+	if driverName == nil {
+		driverName = &defaultDriverName
+	}
+
+	viewshedOpts := viewshedOpts{}
+	for _, opt := range opts {
+		opt.setViewshedOpt(&viewshedOpts)
+	}
+
+	copts := sliceToCStringArray(viewshedOpts.creation)
+	defer copts.free()
+	driver := unsafe.Pointer(C.CString(string(*driverName)))
+	defer C.free(unsafe.Pointer(driver))
+	targetRaster := unsafe.Pointer(C.CString(targetRasterName))
+	defer C.free(unsafe.Pointer(targetRaster))
+
+	cgc := createCGOContext(nil, viewshedOpts.errorHandler)
+	dsRet := C.godalViewshedGenerate(cgc.cPointer(), targetBand.handle(), (*C.char)(driver), (*C.char)(targetRaster), copts.cPointer(), C.double(observerX),
+		C.double(observerY), C.double(observerHeight), C.double(targetHeight), C.double(visibleVal), C.double(invisibleVal), C.double(outOfRangeVal),
+		C.double(noDataVal), C.double(curveCoeff), C.uint(mode), C.double(maxDistance), C.uint(heightMode))
+	if err := cgc.close(); err != nil {
+		return nil, err
+	}
+
+	return &Dataset{majorObject{C.GDALMajorObjectH(dsRet)}}, nil
+}
+
 // Nearblack runs the library version of nearblack
 //
 // See the nearblack doc page to determine the valid flags/opts that can be set in switches.

godal.h

@@ -171,6 +171,10 @@ extern "C" {
 	GDALDatasetH godalGrid(cctx *ctx, const char *pszDest, GDALDatasetH hSrcDS, char **switches);
 	GDALDatasetH godalNearblack(cctx *ctx, const char *pszDest, GDALDatasetH hDstDS, GDALDatasetH hSrcDS, char **switches);
 	GDALDatasetH godalDem(cctx *ctx, const char *pszDest, const char *pszProcessing, const char *pszColorFilename, GDALDatasetH hSrcDS, char **switches);
+	GDALDatasetH godalViewshedGenerate(cctx *ctx, GDALRasterBandH bnd, const char *pszDriverName, const char *pszTargetRasterName, const char **papszCreationOptions,
+							   double dfObserverX, double dfObserverY, double dfObserverHeight, double dfTargetHeight, double dfVisibleVal, 
+							   double dfInvisibleVal, double dfOutOfRangeVal,  double dfNoDataVal, double dfCurvCoeff, GUInt32 eMode, 
+							   double dfMaxDistance, GUInt32 heightMode);
 
 	typedef struct {
 		const GDAL_GCP *gcpList;

godal_test.go

@@ -3851,7 +3851,7 @@
 
 func (dl *debugLogger) L(ec ErrorCategory, code int, msg string) error {
 	if ec >= CE_Warning {
 		return fmt.Errorf(msg)
 	}
 	if ec == CE_Debug {
 		dl.logs += ",GOTESTDEBUG:" + msg
@@ -4440,6 +4440,179 @@
 	assert.Error(t, err)
 }
 
+// Test Ported from: https://github.com/OSGeo/gdal/blob/6cdae8b8f7d09ecf67e24959e984d2e7bbe3ee62/autotest/cpp/test_viewshed.cpp#L98
+func TestViewshedSimpleHeight(t *testing.T) {
+	if !CheckMinVersion(3, 1, 0) {
+		_, err := Viewshed(Band{}, nil, "none", 1, 1, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal)
+		assert.EqualError(t, err, "failed to run, 'viewshed' not supported on GDAL versions < 3.1.0")
+		return
+	} else if !CheckMinVersion(3, 4, 2) {
+		_, err := Viewshed(Band{}, nil, "none", 1, 1, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal)
+		assert.EqualError(t, err, "cannot run 'viewshed' with GDAL version <= 3.4.1, as some tests produce invalid results under these conditions")
+		return
+	}
+
+	// setup common to all scopes
+	var (
+		ehc    = eh()
+		driver = DriverName("MEM")
+
+		identity = [6]float64{0, 1, 0, 0, 0, 1}
+
+		xLen = 5
+		yLen = 5
+		in   = []int8{
+			-1, 0, 1, 0, -1,
+			-1, 2, 0, 4, -1,
+			-1, 1, 0, -1, -1,
+			0, 3, 0, 2, 0,
+			-1, 0, 0, 3, -1,
+		}
+		observable = []float64{
+			4, 2, 0, 4, 8,
+			3, 2, 0, 4, 3,
+			2, 1, 0, -1, -2,
+			4, 3, 0, 2, 1,
+			6, 3, 0, 2, 4,
+		}
+	)
+	vrtDs, err := Create(Memory, "", 1, Int8, xLen, yLen)
+	if err != nil {
+		t.Error(err)
+		return
+	}
+	err = vrtDs.SetGeoTransform(identity)
+	if err != nil {
+		t.Error(err)
+		return
+	}
+	err = vrtDs.Bands()[0].IO(IOWrite, 0, 0, in, xLen, yLen)
+	if err != nil {
+		t.Error(err)
+		return
+	}
+
+	// from cpp scope 1: normal
+	{
+		rds, err := Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal, ErrLogger(ehc.ErrorHandler))
+		if err != nil {
+			t.Error(err)
+			return
+		}
+		defer rds.Close()
+
+		out := make([]int8, xLen*yLen)
+		err = rds.Bands()[0].IO(IORead, 0, 0, out, xLen, yLen)
+		if err != nil {
+			t.Error(err)
+			return
+		}
+
+		expected := make([]int8, xLen*yLen)
+		for i := 0; i < len(in); i++ {
+			var v int8 = 0
+			if in[i] >= int8(observable[i]) {
+				v = 127
+			}
+			expected[i] = v
+		}
+		assert.Equal(t, expected, out)
+	}
+
+	// from cpp scope 2: dem
+	// NOTE: This test fails in releases older than 3.10
+	if CheckMinVersion(3, 10, 0) {
+		rds, err := Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, MinTargetHeightFromDem)
+		if err != nil {
+			t.Error(err)
+			return
+		}
+		defer rds.Close()
+
+		dem := make([]float64, xLen*yLen)
+		err = rds.Bands()[0].IO(IORead, 0, 0, dem, xLen, yLen)
+		if err != nil {
+			t.Error(err)
+			return
+		}
+
+		expected := make([]float64, xLen*yLen)
+		copy(expected, observable)
+		for i := 0; i < len(expected); i++ {
+			expected[i] = math.Max(0.0, expected[i])
+		}
+		assert.Equal(t, expected, dem)
+	} else {
+		_, err := Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, MinTargetHeightFromDem)
+		assert.EqualError(t, err, "height mode CANNOT be `MinTargetHeightFromDem` when running a GDAL version < 3.10, as some tests produce invalid results under these conditions")
+	}
+
+	// from cpp scope 3: ground
+	{
+		rds, err := Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, MinTargetHeightFromGround)
+		if err != nil {
+			t.Error(err)
+			return
+		}
+		defer rds.Close()
+
+		ground := make([]float64, xLen*yLen)
+		err = rds.Bands()[0].IO(IORead, 0, 0, ground, xLen, yLen)
+		if err != nil {
+			t.Error(err)
+			return
+		}
+
+		expected := make([]float64, xLen*yLen)
+		for i := 0; i < len(expected); i++ {
+			expected[i] = math.Max(0.0, observable[i]-float64(in[i]))
+		}
+		assert.Equal(t, expected, ground)
+	}
+}
+
+func TestViewshedCreationOptions(t *testing.T) {
+	if !CheckMinVersion(3, 1, 0) {
+		_, err := Viewshed(Band{}, nil, "none", 1, 1, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal)
+		assert.EqualError(t, err, "failed to run, 'viewshed' not supported on GDAL versions < 3.1.0")
+		return
+	} else if !CheckMinVersion(3, 4, 2) {
+		_, err := Viewshed(Band{}, nil, "none", 1, 1, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal)
+		assert.EqualError(t, err, "cannot run 'viewshed' with GDAL version <= 3.4.1, as some tests produce invalid results under these conditions")
+		return
+	}
+
+	var (
+		driver  = GTiff
+		tmpname = tempfile()
+	)
+	defer os.Remove(tmpname)
+	vrtDs, err := Create(driver, tmpname, 1, Int8, 20, 20, CreationOption("TILED=YES", "BLOCKXSIZE=128", "BLOCKYSIZE=128"))
+	if err != nil {
+		t.Error(err)
+		return
+	}
+	identity := [6]float64{0, 1, 0, 0, 0, 1}
+	err = vrtDs.SetGeoTransform(identity)
+	if err != nil {
+		t.Error(err)
+		return
+	}
+
+	// Invalid - with error logger
+	ehc := eh()
+	_, err = Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal, CreationOption("INVALID_OPT=BAR"), ErrLogger(ehc.ErrorHandler))
+	assert.Error(t, err)
+
+	// Invalid - no error logger
+	_, err = Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal, CreationOption("INVALID_OPT=BAR"))
+	assert.Error(t, err)
+
+	// Valid
+	_, err = Viewshed(vrtDs.Bands()[0], &driver, "none", 2, 2, 0, 0, 255, 0, 0, -1, 0, MEdge, 0, Normal, CreationOption("TILED=YES", "BLOCKXSIZE=128", "BLOCKYSIZE=128"))
+	assert.NoError(t, err)
+}
+
 func TestNearblackBlack(t *testing.T) {
 	// 1. Create an image, linearly interpolated, from black (on the left) to white (on the right), using `Grid()`
 	var (
@@ -4489,7 +4662,7 @@
 	defer func() { _ = VSIUnlink(fname) }()
 	defer gridDs.Close()
 	originalColors := make([]byte, outXSize*outYSize)
 	gridDs.Read(0, 0, originalColors, outXSize, outYSize)
 
 	// 2. Put the Dataset generated above, through the `Nearblack` function, to set pixels near BLACK to BLACK
 	argsNbString := "-near 10 -nb 0"
@@ -4502,7 +4675,7 @@
 	defer func() { _ = VSIUnlink(fname2) }()
 	defer nbDs.Close()
 	nearblackColors := make([]byte, outXSize*outYSize)
 	nbDs.Read(0, 0, nearblackColors, outXSize, outYSize)
 
 	// 3. Test on all rows that pixels where abs(0 - pixelValue) <= 10, are set to black (0)
 	for i := 0; i < outYSize; i++ {
@@ -4561,7 +4734,7 @@
 	defer func() { _ = VSIUnlink(fname) }()
 	defer gridDs.Close()
 	originalColors := make([]byte, outXSize*outYSize)
 	gridDs.Read(0, 0, originalColors, outXSize, outYSize)
 
 	// 2. Put the Dataset generated above, through the `Nearblack` function, to set pixels near WHITE to WHITE
 	argsNbString := "-near 10 -nb 0 -white"
@@ -4574,7 +4747,7 @@
 	defer func() { _ = VSIUnlink(fname2) }()
 	defer nbDs.Close()
 	nearblackColors := make([]byte, outXSize*outYSize)
 	nbDs.Read(0, 0, nearblackColors, outXSize, outYSize)
 
 	// 3. Test on all rows that pixels where abs(255 - pixelValue) <= 10, are set to white (255)
 	for i := 0; i < outYSize; i++ {
@@ -4660,7 +4833,7 @@
 	defer func() { _ = VSIUnlink(fname) }()
 	defer gridDs.Close()
 	originalColors := make([]byte, outXSize*outYSize)
 	gridDs.Read(0, 0, originalColors, outXSize, outYSize)
 
 	// 2. Put the Dataset generated above, through the `Nearblack` function, to set pixels near BLACK to BLACK
 	nbDs, err := Create(Memory, "nbDs", 1, Byte, outXSize, outYSize)
@@ -4676,7 +4849,7 @@
 		return
 	}
 	nearblackColors := make([]byte, outXSize*outYSize)
 	nbDs.Read(0, 0, nearblackColors, outXSize, outYSize)
 
 	// 3. Test on all rows that pixels where abs(0 - pixelValue) <= 10, are set to black (0)
 	for i := 0; i < outYSize; i++ {
@@ -5088,7 +5261,7 @@
 			return
 		}
 	}
 	gridDs.SetProjection("epsg:32632")
 	defer func() { _ = VSIUnlink(fname) }()
 	defer gridDs.Close()
 
@@ -5223,7 +5396,7 @@
 			return
 		}
 	}
 	gridDs.SetProjection("epsg:32632")
 	defer func() { _ = VSIUnlink(fname) }()
 	defer gridDs.Close()
 

options.go

@@ -894,6 +894,7 @@ func CreationOption(opts ...string) interface {
 	DatasetVectorTranslateOption
 	GMLExportOption
 	RasterizeOption
+	ViewshedOption
 } {
 	return creationOpt{opts}
 }
@@ -916,6 +917,9 @@ func (co creationOpt) setGMLExportOpt(gmlo *gmlExportOpts) {
 func (co creationOpt) setRasterizeOpt(o *rasterizeOpts) {
 	o.create = append(o.create, co.creation...)
 }
+func (co creationOpt) setViewshedOpt(dc *viewshedOpts) {
+	dc.creation = append(dc.creation, co.creation...)
+}
 
 type configOpt struct {
 	config []string
@@ -1276,6 +1280,16 @@ type DemOption interface {
 	setDemOpt(demOpt *demOpts)
 }
 
+type viewshedOpts struct {
+	creation     []string
+	errorHandler ErrorHandler
+}
+
+// ViewshedOption is an option that can be passed to Viewshed()
+type ViewshedOption interface {
+	setViewshedOpt(viewshedOpt *viewshedOpts)
+}
+
 type setGCPsOpts struct {
 	errorHandler ErrorHandler
 	projString   string