Clamp BBox on round-earth before passing to ST_IntersectsRect* in Qgs…

…HanaFeatureIterator
qgis · Dec 7, 2020 · f287fdd · f287fdd
1 parent 444276e
commit f287fdd
Show file tree

Hide file tree

Showing 4 changed files with 100 additions and 19 deletions.
diff --git a/src/providers/hana/qgshanafeatureiterator.cpp b/src/providers/hana/qgshanafeatureiterator.cpp
@@ -27,19 +27,68 @@
 #include "qgsmessagelog.h"
 #include "qgssettings.h"
 
+#include "ogr_srs_api.h"
+
+namespace
+{
+  double getAngularUnits( const QgsCoordinateReferenceSystem &crs )
+  {
+    OGRSpatialReferenceH hCRS = OSRNewSpatialReference( nullptr );
+    int errcode = OSRImportFromProj4( hCRS, crs.toProj().toUtf8() );
+    if ( errcode != OGRERR_NONE )
+    {
+      if ( hCRS )
+        OSRRelease( hCRS );
+      throw QgsHanaException( "Unable to parse a spatial reference system" );
+    }
+
+    char *angularUnits = nullptr;
+    double factor = OSRGetAngularUnits( hCRS, &angularUnits );
+    OSRRelease( hCRS );
+    return factor;
+  }
+
+  QgsRectangle clampBBOX( QgsRectangle bbox, const QgsCoordinateReferenceSystem &crs, double allowedExcessFactor )
+  {
+    // In geographic CRS', HANA will reject any points outside the "normalized"
+    // range, which is (in radian) [-PI;PI] for longitude and [-PI/2;PI/2] for
+    // latitude. As QGIS seems to expect that larger bounding boxes should
+    // be allowed and should not wrap-around, we clamp the bounding boxes for
+    // geographic CRS' here.
+    if ( !crs.isGeographic() )
+      return bbox;
+
+    double factor = getAngularUnits( crs );
+
+    double minx = -M_PI / factor;
+    double maxx = M_PI / factor;
+    double spanx = maxx - minx;
+    minx -= allowedExcessFactor * spanx;
+    maxx += allowedExcessFactor * spanx;
+
+    double miny = -0.5 * M_PI / factor;
+    double maxy = 0.5 * M_PI / factor;
+    double spany = maxy - miny;
+    miny -= allowedExcessFactor * spany;
+    maxy += allowedExcessFactor * spany;
+
+    return bbox.intersect( QgsRectangle( minx, miny, maxx, maxy ) );
+  }
+}
+
 QgsHanaFeatureIterator::QgsHanaFeatureIterator(
   QgsHanaFeatureSource *source,
   bool ownSource,
   const QgsFeatureRequest &request )
   : QgsAbstractFeatureIteratorFromSource<QgsHanaFeatureSource>( source, ownSource, request )
+  , mDatabaseVersion( source->mDatabaseVersion )
   , mConnection( source->mUri )
   , mSrsExtent( source->mSrsExtent )
   , mFidColumn( source->mFidColumn )
 {
-  mClosed = true;
-
   if ( mConnection.isNull() )
   {
+    mClosed = true;
     iteratorClosed();
     return;
   }
@@ -53,20 +102,19 @@ QgsHanaFeatureIterator::QgsHanaFeatureIterator(
   }
   catch ( QgsCsException & )
   {
-    iteratorClosed();
+    close();
     return;
   }
 
   try
   {
     mSqlQuery = buildSqlQuery( request );
-    mClosed = false;
 
     rewind();
   }
   catch ( const QgsHanaException & )
   {
-    iteratorClosed();
+    close();
   }
 }
 
@@ -83,8 +131,9 @@ bool QgsHanaFeatureIterator::rewind()
   mResultSet.reset();
   if ( !( mFilterRect.isNull() || mFilterRect.isEmpty() ) && mSource->isSpatial() && mHasGeometryColumn )
   {
-    QString ll = QStringLiteral( "POINT(%1 %2)" ).arg( QString::number( mFilterRect.xMinimum() ),  QString::number( mFilterRect.yMinimum() ) );
-    QString ur = QStringLiteral( "POINT(%1 %2)" ).arg( QString::number( mFilterRect.xMaximum() ),  QString::number( mFilterRect.yMaximum() ) );
+    QgsRectangle filterRect = getFilterRect();
+    QString ll = QStringLiteral( "POINT(%1 %2)" ).arg( QString::number( filterRect.xMinimum() ),  QString::number( filterRect.yMinimum() ) );
+    QString ur = QStringLiteral( "POINT(%1 %2)" ).arg( QString::number( filterRect.xMaximum() ),  QString::number( filterRect.yMaximum() ) );
     mResultSet = mConnection->executeQuery( mSqlQuery, { ll, mSource->mSrid, ur, mSource->mSrid } );
   }
   else
@@ -95,6 +144,7 @@ bool QgsHanaFeatureIterator::rewind()
 
 bool QgsHanaFeatureIterator::close()
 {
+
   if ( mClosed )
     return false;
 
@@ -114,7 +164,7 @@ bool QgsHanaFeatureIterator::fetchFeature( QgsFeature &feature )
 {
   feature.setValid( false );
 
-  if ( mClosed )
+  if ( mClosed || !mResultSet )
     return false;
 
   if ( !mResultSet->next() )
@@ -164,7 +214,6 @@ bool QgsHanaFeatureIterator::fetchFeature( QgsFeature &feature )
   feature.setFields( mSource->mFields ); // allow name-based attribute lookups
   geometryToDestinationCrs( feature, mTransform );
 
-
   return true;
 }
 
@@ -176,16 +225,31 @@ bool QgsHanaFeatureIterator::nextFeatureFilterExpression( QgsFeature &feature )
     return fetchFeature( feature );
 }
 
-QString QgsHanaFeatureIterator::getBBOXFilter( const QVersionNumber &dbVersion ) const
+QString QgsHanaFeatureIterator::getBBOXFilter( ) const
 {
-  if ( dbVersion.majorVersion() == 1 )
+  if ( mDatabaseVersion.majorVersion() == 1 )
     return QStringLiteral( "%1.ST_SRID(%2).ST_IntersectsRect(ST_GeomFromText(?, ?), ST_GeomFromText(?, ?)) = 1" )
            .arg( QgsHanaUtils::quotedIdentifier( mSource->mGeometryColumn ), QString::number( mSource->mSrid ) );
   else
     return QStringLiteral( "%1.ST_IntersectsRectPlanar(ST_GeomFromText(?, ?), ST_GeomFromText(?, ?)) = 1" )
            .arg( QgsHanaUtils::quotedIdentifier( mSource->mGeometryColumn ) );
 }
 
+QgsRectangle QgsHanaFeatureIterator::getFilterRect() const
+{
+  const QgsCoordinateReferenceSystem &crs = mSource->mCrs;
+  if ( !crs.isGeographic() )
+    return mFilterRect;
+
+  if ( mDatabaseVersion.majorVersion() > 1 )
+    return clampBBOX( mFilterRect, crs, 0.0 );
+
+  int srid = QgsHanaUtils::toPlanarSRID( mSource->mSrid );
+  if ( srid == mSource->mSrid )
+    return mFilterRect;
+  return clampBBOX( mFilterRect, crs, 0.5 );
+}
+
 bool QgsHanaFeatureIterator::prepareOrderBy( const QList<QgsFeatureRequest::OrderByClause> &orderBys )
 {
   Q_UNUSED( orderBys )
@@ -295,7 +359,7 @@ QString QgsHanaFeatureIterator::buildSqlQuery( const QgsFeatureRequest &request
   QStringList sqlFilter;
   // Set spatial filter
   if ( mSource->isSpatial() && mHasGeometryColumn && !( filterRect.isNull() || filterRect.isEmpty() ) )
-    sqlFilter.push_back( getBBOXFilter( QgsHanaUtils::toHANAVersion( mConnection->getDatabaseVersion() ) ) );
+    sqlFilter.push_back( getBBOXFilter() );
 
   if ( !mSource->mQueryWhereClause.isEmpty() )
     sqlFilter.push_back( mSource->mQueryWhereClause );
@@ -385,7 +449,8 @@ QString QgsHanaFeatureIterator::buildSqlQuery( const QgsFeatureRequest &request
 }
 
 QgsHanaFeatureSource::QgsHanaFeatureSource( const QgsHanaProvider *p )
-  : mUri( p->mUri )
+  : mDatabaseVersion( p->mDatabaseVersion )
+  , mUri( p->mUri )
   , mSchemaName( p->mSchemaName )
   , mTableName( p->mTableName )
   , mFidColumn( p->mFidColumn )

diff --git a/src/providers/hana/qgshanafeatureiterator.h b/src/providers/hana/qgshanafeatureiterator.h
@@ -36,6 +36,7 @@ class QgsHanaFeatureSource : public QgsAbstractFeatureSource
     bool isSpatial() const { return !mGeometryColumn.isEmpty() && mGeometryType != QgsWkbTypes::Unknown; }
 
   private:
+    QVersionNumber mDatabaseVersion;
     QgsDataSourceUri mUri;
     QString mSchemaName;
     QString mTableName;
@@ -74,16 +75,18 @@ class QgsHanaFeatureIterator : public QgsAbstractFeatureIteratorFromSource<QgsHa
     bool prepareOrderBy( const QList<QgsFeatureRequest::OrderByClause> &orderBys ) override;
 
     QString buildSqlQuery( const QgsFeatureRequest &request );
-    QString getBBOXFilter( const QVersionNumber &dbVersion ) const;
+    QString getBBOXFilter() const;
+    QgsRectangle getFilterRect() const;
 
   private:
+    const QVersionNumber mDatabaseVersion;
     QgsHanaConnectionRef mConnection;
     QgsHanaResultSetRef mResultSet;
     QString mSqlQuery = QString( "" );
     QgsRectangle mFilterRect;
-    QgsRectangle mSrsExtent;
+    const QgsRectangle mSrsExtent;
     QgsAttributeList mAttributesToFetch;
-    QString mFidColumn;
+    const QString mFidColumn;
     QgsCoordinateTransform mTransform;
     bool mHasAttributes = false;
     bool mHasGeometryColumn = false;

diff --git a/src/providers/hana/qgshanaprovider.cpp b/src/providers/hana/qgshanaprovider.cpp
@@ -65,12 +65,14 @@ namespace
 
   void createCoordinateSystem( QgsHanaConnectionRef &conn, const QgsCoordinateReferenceSystem &srs )
   {
-    OGRSpatialReferenceH hCRS = nullptr;
-    hCRS = OSRNewSpatialReference( nullptr );
+    OGRSpatialReferenceH hCRS = OSRNewSpatialReference( nullptr );
     int errcode = OSRImportFromProj4( hCRS, srs.toProj().toUtf8() );
-
     if ( errcode != OGRERR_NONE )
+    {
+      if ( hCRS )
+        OSRRelease( hCRS );
       throw QgsHanaException( "Unable to parse a spatial reference system" );
+    }
 
     QgsCoordinateReferenceSystem srsWGS84;
     srsWGS84.createFromString( "EPSG:4326" );
@@ -111,6 +113,8 @@ namespace
                         QgsHanaUtils::quotedString( srs.toWkt() ),
                         QgsHanaUtils::quotedString( srs.toProj() ) );
 
+    OSRRelease( hCRS );
+
     QString errorMessage;
     conn->execute( sql, &errorMessage );
 

diff --git a/tests/src/python/test_provider_hana.py b/tests/src/python/test_provider_hana.py
@@ -27,6 +27,7 @@
     QgsFeatureRequest,
     QgsFeature,
     QgsProviderRegistry,
+    QgsRectangle,
     QgsSettings)
 from qgis.testing import start_app, unittest
 from test_hana_utils import QgsHanaProviderUtils
@@ -293,6 +294,14 @@ def testBinaryTypeEdit(self):
         expected = {1: QByteArray(b'bbbvx'), 2: QByteArray(b'dddd')}
         self.assertEqual(values, expected)
 
+    def testFilterRectOutsideSrsExtent(self):
+        """Test filterRect which partially lies outside of the srs extent"""
+        self.source.setSubsetString(None)
+        extent = QgsRectangle(-103, 46, -25, 97)
+        result = set([f[self.pk_name] for f in self.source.getFeatures(QgsFeatureRequest().setFilterRect(extent))])
+        expected = set([1, 2, 4, 5])
+        assert set(expected) == result, f'Expected {expected} and got {result} when testing setFilterRect {extent}'
+
     def testEncodeDecodeUri(self):
         """Test HANA encode/decode URI"""
         md = QgsProviderRegistry.instance().providerMetadata('hana')