Skip to content

Commit cf119e5

Browse files
authored
LOS: Add (placeholder) intersection location parameters (#9585)
Implementation will come later
1 parent 309d4e0 commit cf119e5

File tree

7 files changed

+122
-64
lines changed

7 files changed

+122
-64
lines changed

alg/gdal_alg.h

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -357,11 +357,10 @@ GDALDatasetH CPL_DLL GDALViewshedGenerate(
357357
void *pProgressArg, GDALViewshedOutputType heightMode,
358358
CSLConstList papszExtraOptions);
359359

360-
bool CPL_DLL GDALIsLineOfSightVisible(const GDALRasterBandH, const int xA,
361-
const int yA, const double zA,
362-
const int xB, const int yB,
363-
const double zB,
364-
CSLConstList papszOptions);
360+
bool CPL_DLL GDALIsLineOfSightVisible(
361+
const GDALRasterBandH, const int xA, const int yA, const double zA,
362+
const int xB, const int yB, const double zB, int *pnxTerrainIntersection,
363+
int *pnyTerrainIntersection, CSLConstList papszOptions);
365364

366365
/************************************************************************/
367366
/* Rasterizer API - geometries burned into GDAL raster. */

alg/los.cpp

Lines changed: 24 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -162,21 +162,31 @@ static bool IsAboveTerrain(const GDALRasterBandH hBand, const int x,
162162
* For example, datasets referenced against geographic coordinate at high latitudes may have issues.
163163
*
164164
* @param hBand The band to read the DEM data from. This must NOT be null.
165-
*
165+
*
166166
* @param xA The X location (raster column) of the first point to check on the raster.
167167
*
168168
* @param yA The Y location (raster row) of the first point to check on the raster.
169-
*
169+
*
170170
* @param zA The Z location (height) of the first point to check.
171-
*
171+
*
172172
* @param xB The X location (raster column) of the second point to check on the raster.
173173
*
174174
* @param yB The Y location (raster row) of the second point to check on the raster.
175-
*
175+
*
176176
* @param zB The Z location (height) of the second point to check.
177-
*
177+
*
178+
* @param[out] pnxTerrainIntersection The X location where the LOS line
179+
* intersects with terrain, or nullptr if it does not intersect
180+
* terrain. Not implemented currently (*pnxTerrainIntersection
181+
* is set to -1)
182+
*
183+
* @param[out] pnyTerrainIntersection The Y location where the LOS line
184+
* intersects with terrain, or nullptr if it does not intersect
185+
* terrain. Not implemented currently (*pnyTerrainIntersection
186+
* is set to -1)
187+
*
178188
* @param papszOptions Options for the line of sight algorithm (currently ignored).
179-
*
189+
*
180190
* @return True if the two points are within Line of Sight.
181191
*
182192
* @since GDAL 3.9
@@ -185,10 +195,18 @@ static bool IsAboveTerrain(const GDALRasterBandH hBand, const int x,
185195
bool GDALIsLineOfSightVisible(const GDALRasterBandH hBand, const int xA,
186196
const int yA, const double zA, const int xB,
187197
const int yB, const double zB,
198+
int *pnxTerrainIntersection,
199+
int *pnyTerrainIntersection,
188200
CPL_UNUSED CSLConstList papszOptions)
189201
{
190202
VALIDATE_POINTER1(hBand, "GDALIsLineOfSightVisible", false);
191203

204+
if (pnxTerrainIntersection)
205+
*pnxTerrainIntersection = -1;
206+
207+
if (pnyTerrainIntersection)
208+
*pnyTerrainIntersection = -1;
209+
192210
// Perform a preliminary check of the start and end points.
193211
if (!IsAboveTerrain(hBand, xA, yA, zA))
194212
{

autotest/alg/los.py

Lines changed: 8 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -37,10 +37,14 @@ def test_los_basic():
3737

3838
mem_ds = gdal.GetDriverByName("MEM").Create("", 2, 1)
3939

40-
assert gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 1, 0, 1)
41-
assert gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 0, 0, 1)
42-
assert not gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, -1, 1, 0, 1)
43-
assert not gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 1, 0, -1)
40+
success, x, y = gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 1, 0, 1)
41+
assert success
42+
assert x == -1
43+
assert y == -1
44+
45+
assert gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 0, 0, 1)[0]
46+
assert not gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, -1, 1, 0, 1)[0]
47+
assert not gdal.IsLineOfSightVisible(mem_ds.GetRasterBand(1), 0, 0, 1, 1, 0, -1)[0]
4448

4549
with pytest.raises(Exception, match="Received a NULL pointer"):
4650
gdal.IsLineOfSightVisible(None, 0, 0, 0, 0, 0, 0)

autotest/cpp/test_alg.cpp

Lines changed: 46 additions & 46 deletions
Original file line numberDiff line numberDiff line change
@@ -301,14 +301,14 @@ TEST_F(test_alg, GDALIsLineOfSightVisible_single_point_dataset)
301301
ASSERT_TRUE(poDS->RasterIO(GF_Write, 0, 0, 1, 1, &val, 1, 1, GDT_Int8, 1,
302302
nullptr, 0, 0, 0, nullptr) == CE_None);
303303
// Both points below terrain
304-
EXPECT_FALSE(
305-
GDALIsLineOfSightVisible(pBand, 0, 0, 0.0, 0, 0, 0.0, nullptr));
304+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 0, 0, 0.0, 0, 0, 0.0, nullptr,
305+
nullptr, nullptr));
306306
// One point below terrain
307-
EXPECT_FALSE(
308-
GDALIsLineOfSightVisible(pBand, 0, 0, 0.0, 0, 0, 43.0, nullptr));
307+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 0, 0, 0.0, 0, 0, 43.0, nullptr,
308+
nullptr, nullptr));
309309
// Both points above terrain
310-
EXPECT_TRUE(
311-
GDALIsLineOfSightVisible(pBand, 0, 0, 44.0, 0, 0, 43.0, nullptr));
310+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 0, 0, 44.0, 0, 0, 43.0, nullptr,
311+
nullptr, nullptr));
312312
}
313313

314314
// Test GDALIsLineOfSightVisible() with 10x10 default dataset
@@ -331,25 +331,25 @@ TEST_F(test_alg, GDALIsLineOfSightVisible_default_square_dataset)
331331
const int y2 = 2;
332332

333333
// Both points are above terrain.
334-
EXPECT_TRUE(
335-
GDALIsLineOfSightVisible(pBand, x1, y1, 1.0, x2, y2, 1.0, nullptr));
334+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, x1, y1, 1.0, x2, y2, 1.0,
335+
nullptr, nullptr, nullptr));
336336
// Flip the order, same result.
337-
EXPECT_TRUE(
338-
GDALIsLineOfSightVisible(pBand, x2, y2, 1.0, x1, y1, 1.0, nullptr));
337+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, x2, y2, 1.0, x1, y1, 1.0,
338+
nullptr, nullptr, nullptr));
339339

340340
// One point is below terrain.
341-
EXPECT_FALSE(
342-
GDALIsLineOfSightVisible(pBand, x1, y1, -1.0, x2, y2, 1.0, nullptr));
341+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, x1, y1, -1.0, x2, y2, 1.0,
342+
nullptr, nullptr, nullptr));
343343
// Flip the order, same result.
344-
EXPECT_FALSE(
345-
GDALIsLineOfSightVisible(pBand, x2, y2, -1.0, x1, y1, 1.0, nullptr));
344+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, x2, y2, -1.0, x1, y1, 1.0,
345+
nullptr, nullptr, nullptr));
346346

347347
// Both points are below terrain.
348-
EXPECT_FALSE(
349-
GDALIsLineOfSightVisible(pBand, x1, y1, -1.0, x2, y2, -1.0, nullptr));
348+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, x1, y1, -1.0, x2, y2, -1.0,
349+
nullptr, nullptr, nullptr));
350350
// Flip the order, same result.
351-
EXPECT_FALSE(
352-
GDALIsLineOfSightVisible(pBand, x2, y2, -1.0, x1, y1, -1.0, nullptr));
351+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, x2, y2, -1.0, x1, y1, -1.0,
352+
nullptr, nullptr, nullptr));
353353
}
354354

355355
// Test GDALIsLineOfSightVisible() through a mountain (not a unit test)
@@ -397,57 +397,57 @@ TEST_F(test_alg, GDALIsLineOfSightVisible_through_mountain)
397397
// Both points are just above terrain, with terrain between.
398398
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, iMesaTopX, iMesaTopY, 222,
399399
iMesaBottomX, iMesaBottomY, 199,
400-
nullptr));
400+
nullptr, nullptr, nullptr));
401401
// Flip the order, same result.
402402
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, iMesaBottomX, iMesaBottomY,
403403
199, iMesaTopX, iMesaTopY, 222,
404-
nullptr));
404+
nullptr, nullptr, nullptr));
405405

406406
// Both points above terrain.
407407
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, iMesaTopX, iMesaTopY, 322,
408408
iMesaBottomX, iMesaBottomY, 322,
409-
nullptr));
409+
nullptr, nullptr, nullptr));
410410

411411
// Both points below terrain.
412412
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, iMesaTopX, iMesaTopY, 0,
413413
iMesaBottomX, iMesaBottomY, 0,
414-
nullptr));
414+
nullptr, nullptr, nullptr));
415415

416416
// Test negative slope bresenham diagonals across the whole raster.
417417
// Both high above terrain.
418-
EXPECT_TRUE(
419-
GDALIsLineOfSightVisible(pBand, 0, 0, 460, 120, 120, 460, nullptr));
418+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 0, 0, 460, 120, 120, 460,
419+
nullptr, nullptr, nullptr));
420420
// Both heights are 1m above in the corners, but middle terrain violates LOS.
421-
EXPECT_FALSE(
422-
GDALIsLineOfSightVisible(pBand, 0, 0, 295, 120, 120, 183, nullptr));
421+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 0, 0, 295, 120, 120, 183,
422+
nullptr, nullptr, nullptr));
423423

424424
// Test positive slope bresenham diagnoals across the whole raster.
425425
// Both high above terrain.
426-
EXPECT_TRUE(
427-
GDALIsLineOfSightVisible(pBand, 0, 120, 460, 120, 0, 460, nullptr));
426+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 0, 120, 460, 120, 0, 460,
427+
nullptr, nullptr, nullptr));
428428
// Both heights are 1m above in the corners, but middle terrain violates LOS.
429-
EXPECT_FALSE(
430-
GDALIsLineOfSightVisible(pBand, 0, 120, 203, 120, 0, 247, nullptr));
429+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 0, 120, 203, 120, 0, 247,
430+
nullptr, nullptr, nullptr));
431431

432432
// Vertical line tests with hill between two points, in both directions.
433-
EXPECT_FALSE(
434-
GDALIsLineOfSightVisible(pBand, 83, 111, 154, 83, 117, 198, nullptr));
435-
EXPECT_FALSE(
436-
GDALIsLineOfSightVisible(pBand, 83, 117, 198, 83, 111, 154, nullptr));
437-
EXPECT_TRUE(
438-
GDALIsLineOfSightVisible(pBand, 83, 111, 460, 83, 117, 460, nullptr));
439-
EXPECT_TRUE(
440-
GDALIsLineOfSightVisible(pBand, 83, 117, 460, 83, 111, 460, nullptr));
433+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 83, 111, 154, 83, 117, 198,
434+
nullptr, nullptr, nullptr));
435+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 83, 117, 198, 83, 111, 154,
436+
nullptr, nullptr, nullptr));
437+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 83, 111, 460, 83, 117, 460,
438+
nullptr, nullptr, nullptr));
439+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 83, 117, 460, 83, 111, 460,
440+
nullptr, nullptr, nullptr));
441441

442442
// Horizontal line tests with hill between two points, in both directions.
443-
EXPECT_FALSE(
444-
GDALIsLineOfSightVisible(pBand, 75, 115, 192, 89, 115, 191, nullptr));
445-
EXPECT_FALSE(
446-
GDALIsLineOfSightVisible(pBand, 89, 115, 191, 75, 115, 192, nullptr));
447-
EXPECT_TRUE(
448-
GDALIsLineOfSightVisible(pBand, 75, 115, 460, 89, 115, 460, nullptr));
449-
EXPECT_TRUE(
450-
GDALIsLineOfSightVisible(pBand, 89, 115, 460, 75, 115, 460, nullptr));
443+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 75, 115, 192, 89, 115, 191,
444+
nullptr, nullptr, nullptr));
445+
EXPECT_FALSE(GDALIsLineOfSightVisible(pBand, 89, 115, 191, 75, 115, 192,
446+
nullptr, nullptr, nullptr));
447+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 75, 115, 460, 89, 115, 460,
448+
nullptr, nullptr, nullptr));
449+
EXPECT_TRUE(GDALIsLineOfSightVisible(pBand, 89, 115, 460, 75, 115, 460,
450+
nullptr, nullptr, nullptr));
451451
}
452452

453453
} // namespace

swig/include/Operations.i

Lines changed: 17 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -639,6 +639,21 @@ GDALDatasetShadow *ViewshedGenerate( GDALRasterBandShadow *srcBand,
639639
/* IsLineOfSightVisible() */
640640
/************************************************************************/
641641

642+
#ifdef SWIGPYTHON
643+
%feature( "kwargs" ) IsLineOfSightVisible;
644+
%apply Pointer NONNULL {GDALRasterBandShadow *band};
645+
%inline %{
646+
void IsLineOfSightVisible(GDALRasterBandShadow *band,
647+
int xA, int yA, double zA,
648+
int xB, int yB, double zB,
649+
bool *pbVisible, int *pnXIntersection, int *pnYIntersection,
650+
char** options = NULL)
651+
{
652+
*pbVisible = GDALIsLineOfSightVisible(band, xA, yA, zA, xB, yB, zB, pnXIntersection, pnYIntersection, options);
653+
}
654+
%}
655+
%clear GDALRasterBandShadow *band;
656+
#else
642657
#ifndef SWIGJAVA
643658
%feature( "kwargs" ) IsLineOfSightVisible;
644659
#endif
@@ -649,10 +664,11 @@ bool IsLineOfSightVisible(GDALRasterBandShadow *band,
649664
int xB, int yB, double zB,
650665
char** options = NULL)
651666
{
652-
return GDALIsLineOfSightVisible(band, xA, yA, zA, xB, yB, zB, options);
667+
return GDALIsLineOfSightVisible(band, xA, yA, zA, xB, yB, zB, NULL, NULL, options);
653668
}
654669
%}
655670
%clear GDALRasterBandShadow *band;
671+
#endif
656672

657673
/************************************************************************/
658674
/* AutoCreateWarpedVRT() */

swig/include/python/docs/gdal_operations_docs.i

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,8 @@ options : dict/list, optional
2929
3030
Returns
3131
-------
32-
bool
33-
True if the two points are within Line of Sight.
32+
(bool, int, int)
33+
First value of tuple is True if the two points are within Line of Sight.
34+
Second value is the X location where the LOS line intersects with terrain (will be set in the future, currently set to -1).
35+
Third value is the Y location where the LOS line intersects with terrain (will be set in the future, currently set to -1).
3436
";

swig/include/python/typemaps_python.i

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3449,3 +3449,22 @@ OBJECT_LIST_INPUT(GDALMDArrayHS);
34493449
/* %typemap(freearg) (int nUsages, GDALRATFieldUsage *paeUsages)*/
34503450
CPLFree( $2 );
34513451
}
3452+
3453+
3454+
%typemap(in,numinputs=0) (bool *pbVisible, int *pnXIntersection, int *pnYIntersection) ( bool visible = 0, int nxintersection = 0, int nyintersection = 0 )
3455+
{
3456+
/* %typemap(in) (bool *pbVisible, int *pnXIntersection, int *pnYIntersection) */
3457+
$1 = &visible;
3458+
$2 = &nxintersection;
3459+
$3 = &nyintersection;
3460+
}
3461+
3462+
%typemap(argout) (bool *pbVisible, int *pnXIntersection, int *pnYIntersection)
3463+
{
3464+
/* %typemap(argout) (bool *pbVisible, int *pnXIntersection, int *pnYIntersection) */
3465+
PyObject *r = PyTuple_New( 3 );
3466+
PyTuple_SetItem( r, 0, PyBool_FromLong(*$1) );
3467+
PyTuple_SetItem( r, 1, PyLong_FromLong(*$2) );
3468+
PyTuple_SetItem( r, 2, PyLong_FromLong(*$3) );
3469+
$result = t_output_helper($result,r);
3470+
}

0 commit comments

Comments
 (0)