id	summary	reporter	owner	description	type	status	component	version	severity	resolution	keywords	cc	stage	has_patch	needs_docs	needs_tests	needs_better_patch	easy	ui_ux
30267	GeoDjango does not create correct PostGIS index for 3D geometry fields	Casper van der Wel	nobody	"When I create a model with a 3D geometryfield, an index USING gist (geom gist_geometry_ops_nd) is created automatically by geodjango in my PostGIS database. However, the spatial lookups 'intersects' and 'bboverlaps' do not use this index, which I can confirm by explaining the queries. Below some code to reproduce this issue.

I am not sure if this is a bug or if I am reporting a non-intended use case for 3D geometryfields. It appears to me as a pretty straightforward use case.

The model:

{{{
from django.contrib.gis.db import models
from django.db import connections
from django.db.models.query import QuerySet

class QuerySetExplain(QuerySet):
    def explain(self):
        cursor = connections[self.db].cursor()
        query, params = self.query.sql_with_params()
        cursor.execute('EXPLAIN %s' % query, params)
        return '\n'.join(r[0] for r in cursor.fetchall())

class PolyModel(models.Model):
    objects = QuerySetExplain.as_manager()
    geometry = models.PolygonField(srid=28992, dim=3)
}}}


After creating and running migrations, I explain the following query:

{{{
>>> from django.contrib.gis.geos import Polygon
>>> from my_app.models import PolyModel
>>> bbox = Polygon.from_bbox((0, 0, 1, 1))
>>> qs = PolyModel.objects.filter(geometry__intersects=bbox)
>>> print(qs.explain())

Seq Scan on my_app_polymodel  (cost=0.00..2964.00 rows=1 width=241)
  Filter: ((geometry && '0103000020E61000000100000005000000000000000000000000000000000000000000000000000000000000000000F03F000000000000F03F000000000000F03F000000000000F03F000000000000000000000000000000000000000000000000'::geometry) AND _st_intersects(geometry, '0103000020E61000000100000005000000000000000000000000000000000000000000000000000000000000000000F03F000000000000F03F000000000000F03F000000000000F03F000000000000000000000000000000000000000000000000'::geometry))
}}}

This says that is used a sequential scan. **It did not use the index**
The solution for me is creating the indices ""USING gist (geometry)"" with raw SQL. 
Looking at PostGIS manuals, it seems you have to use the operator &&& to actually use an ND index. I can confirm that in the sql shell.

Versions: I am using django 1.11.20, PostGRES 9.3 with PostGIS 2.4. I also reproduced the issue in raw SQL using PostGRES 10. 

After glancing the django source code, it appears to me that it persists in master. But I did not run any tests.

A quickfix would be always creating 2 indexes on 3D geometry fields. But that does seem a waste of resources. Using &&& instead of && as the 'bboverlaps' lookup would probably a better fix, but I am not such a PostGIS expert that I can oversee the changes necessary in the django ORM.

Thanks in advance for looking at this bug report.

"	Cleanup/optimization	new	GIS	1.11	Normal		geodjango, postgis	Ülgen Sarıkavak	Accepted	0	0	0	0	0	0