overlay_nearest() and other overlay_* functions return incorrect results when used with the same layer (target layer = source layer) in some situations

[agiudiceandrea]

What is the bug or the crash?

The overlay_nearest() function, used with a target layer which is the same as the source layer and such layer is in a Shapefile or GeoJSON or FlatGeobuf format, incorrectly returns always the same current geometry instead of the nearest geometry.

When the layer is in a Memory or GeoPackage or Spatialite format, then the function correctly returns the nearest geometry to the current geometry.

The issue doesn't occur with the overlay_equals() and overlay_intersects() functions when the cache parameter is not specified or is set to false (the default value).
The issue does also occur with the overlay_equals() and overlay_intersects() functions when the cache parameter is set to true.

Steps to reproduce the issue

1a) create a new scratch (Memory) point layer
1b) add a single point to the memory layer
1c) use the expression overlay_nearest(@layer,$geometry) in the field calculator for such layer: the preview will correctly show [] (i.e. no nearest geometry returned)

2a) export the Memory layer to a Shapefile (or GeoJSON or FlatGeobuf) layer and add it to the map
2b) use the same expression in the field calculator for the Shapefile layer: the preview will incorrectly show [ <geometry: Point> ] (i.e. the expression incorrectly returns the very some point as the nearest point)

3a) export the Memory layer to a GeoPackage or a Spatialite layer and add it to the map
3b) use the same expression in the field calculator for the GeoPackage or a Spatialite layer: the preview will correctly show [] (i.e. no nearest geometry returned)

Versions

QGIS 3.16.0 and QGIS 3.16.5 on Windows 7 (from OSGeo4W v1)
QGIS 3.16.16, QGIS 3.22.3 and QGIS 3.23.0-Master on Windows 10 (from OSGeo4W v2)

Supported QGIS version

• I'm running a supported QGIS version according to the roadmap.

New profile

• I tried with a new QGIS profile

Additional context

The overlay_nearest() function and the other overlay_* functions have been ported in core QGIS from the refFunction plugin with PR #38405 by @olivierdalang, @m-kuhn and @enricofer.

The functions code has logic in place to allegedly handle the "some layer" (target layer = source layer) situation in the right way:

QGIS/src/core/expression/qgsexpressionfunction.cpp

Lines 6842 to 6875 in d5d3ed6

	QList<QgsFeatureId> fidsList;
	if ( isNearestFunc )
	{
	fidsList = spatialIndex.nearestNeighbor( geometry, sameLayers ? limit + 1 : limit, max_distance );
	}
	else
	{
	fidsList = spatialIndex.intersects( intDomain );
	}
	QListIterator<QgsFeatureId> i( fidsList );
	while ( i.hasNext() )
	{
	QgsFeatureId fId2 = i.next();
	if ( sameLayers && feat.id() == fId2 )
	continue;
	features.append( cachedTarget->getFeature( fId2 ) );
	}
	}
	else
	{
	// If the cache (local spatial index) is not enabled, we directly
	// get the features from the target layer
	request.setFilterRect( intDomain );
	QgsFeatureIterator fit = targetLayer->getFeatures( request );
	QgsFeature feat2;
	while ( fit.nextFeature( feat2 ) )
	{
	if ( sameLayers && feat.id() == feat2.id() )
	continue;
	features.append( feat2 );
	}
	}

The function, in this case, should drop the current feature from the result.

See also olivierdalang@751ad05

[agiudiceandrea]

Additionally, the overlay_nearest(), when used with the same layer (target layer = source layer), doesn't respect the limit:=-1 option.
In this case, the functions behave as the parameter limit was set to 0 if the layer is in a Memory or GeoPackage or Spatialite formt, while it behave as the parameter limit was set to 1 if the layer is in a Shapefile or GeoJSON or FlatGeobuf format.

[pigreco]

I confirm

OSGeo4W v2 win 10 64 bit

[pigreco]

https://gis.stackexchange.com/questions/424128/qgis-expression-when-using-overlay-nearest-sometimes-a-feature-considers-itself