[postgis-users] Using SRID with meters, but still not getting UNITS in meters with ST_DWithin

Mon Apr 11 20:00:42 PDT 2016

Go back to basics, use a geography column from the start:

ALTER TABLE tpostalcoordinate ADD COLUMN geog geography(point, 4326);
UPDATE tpostalcoordinate SET geog =
geography(ST_SetSRID(ST_MakePoint(longitude, latitude), 4326));
CREATE INDEX tpostalcoordinate_gix ON tpostalcoordinate USING GIST (geog);
SELECT a.*
  FROM tpostalcoordinate a
  JOIN  tpostalcoordinate  b
  ON ST_DWithin(a.geog, b.geog, 10000)
  WHERE
    a.countrycode2tcountry = 'US' AND
    b.countrycode2tcountry = 'US' AND
    b.postalcode = '94404';

A geocentric system is one where every point is referenced relative to
the center of the earth, so they all have an x,y,z. It's not something
you want to use.

P

On Mon, Apr 11, 2016 at 12:31 PM, Michael Moore <michaeljmoore at gmail.com> wrote:
>
>
> On Mon, Apr 11, 2016 at 11:50 AM, Paul Ramsey <pramsey at cleverelephant.ca>
> wrote:
>>
>> Before I can even begin to address what's going on with the
>> index/query side, I have to ask about the coordinate system:
>>
>> You are using a geocentric system for postal code points?
>>
>> https://epsg.io/4896
>>
>> This seems very odd indeed. If you were flying a satellite or drilling
>> a deep well, maybe you'd use geocentric coordinates, but simple postal
>> code queries? So, perhaps there's a core error to deal with first: why
>> are you in EPSG:4896?
>>
>> P
>>
>>
>> On Mon, Apr 11, 2016 at 11:21 AM, Michael Moore <michaeljmoore at gmail.com>
>> wrote:
>> > I am trying to find all zip codes withing a given range of current zip
>> > code.
>> > For example, User is at zip code 95076 and want to know all zip codes
>> > within
>> > a 30 mile range. This will always be a short distance, nothing over 50
>> > miles. The source zip code is any place in the USA.  My input table has
>> > a
>> > latitude field, a longitude field and a POINT field for each zip code.
>> > I want to supply in input in meters, not degrees. I do NOT want to cast
>> > to
>> > geography because this prevents the index from being used.
>> > Here is a working query that uses the index, but ST_DWithin  is not
>> > understanding ".05" as meters:
>> >
>> >
>> > lcd1_dev=>   select t2.city, t2.postalcode,
>> > ST_Distance(t1.geo_position,t2.geo_position)
>> > lcd1_dev->   from      tpostalcoordinate t1
>> > lcd1_dev->   left join tpostalcoordinate t2 on
>> > ST_DWithin(t1.geo_position,t2.geo_position , .05)
>> > lcd1_dev->   where t1.postalcode = '94404'and t1.countrycode2tcountry =
>> > 'US'
>> > and t2.countrycode2tcountry= 'US';
>> >      city     | postalcode |     st_distance
>> > --------------+------------+---------------------
>> >  Redwood City | 94065      |  0.0273766323714193
>> >  San Mateo    | 94408      | 0.00504738546179631
>> >  Belmont      | 94002      |  0.0440065904155286
>> >  San Mateo    | 94404      |                   0
>> >  San Mateo    | 94403      |  0.0370314731005901
>> >  San Mateo    | 94407      |  0.0416118372581607
>> >
>> > This shows that I am using SRID 4896
>> > lcd1_dev=> select * from geometry_columns where f_table_name =
>> > 'tpostalcoordinate';
>> >  f_table_catalog | f_table_schema |   f_table_name    |
>> > f_geometry_column |
>> > coord_dimension | srid | type
>> >
>> > -----------------+----------------+-------------------+-------------------+-----------------+------+-------
>> >  lcd1_dev        | qsn_app        | tpostalcoordinate | geo_position
>> > |
>> > 2 | 4896 | POINT
>> >
>> > 4896 is UNIT "metre" as show here:
>> >  4896 | EPSG      |      4896 |
>> >
>> > GEOCCS["ITRF2005",DATUM["International_Terrestrial_Reference_Frame_2005",SPHEROID["GRS
>> > 1980",6378137,298.257222101,AUTHORIT
>> >
>> > Y["EPSG","7019"]],AUTHORITY["EPSG","6896"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Geocentric
>> > X",OTH
>> > ER],AXIS["Geocentric Y",OTHER],AXIS["Geocentric
>> > Z",NORTH],AUTHORITY["EPSG","4896"]] | +proj=geocent +ellps=GRS80
>> > +units=m
>> > +no_de
>> >
>> > The following (from EXPLAIN) proves that the index is being used in the
>> > working query shown above:
>> >    ->  Index Scan using tpostalcoordinate_pk on tpostalcoordinate t1
>> > (cost=0.42..8.45 rows=1 width=32)
>> >
>> > What do I need to do to get this in meters without losing my index
>> > access?
>> >
>> > TIA,
>> > Mike
>> >
>> > _______________________________________________
>> > postgis-users mailing list
>> > postgis-users at lists.osgeo.org
>> > http://lists.osgeo.org/mailman/listinfo/postgis-users
>> _______________________________________________
>> postgis-users mailing list
>> postgis-users at lists.osgeo.org
>> http://lists.osgeo.org/mailman/listinfo/postgis-users
>
>
> Here is what we do to create the column on the table:
>
>
>
> 1.    SELECT AddGeometryColumn('qsn_app','tpostalcoordinate',
> 'geo_position',4326,'POINT',2);
>
>
>
> 2.    UPDATE tpostalcoordinate pc set pc.geo_position =
> ST_SetSRID(ST_Point(pc.longitude, pc.latitude), 4326);
>
>
>
> 3.    CREATE INDEX tpostal_geo_position_idx ON tpostalcoordinate USING
> gist(geo_position);
>
>
> I am on EPSG:4896 because I tried about 10 others that were also meters and
> that's just the one I ended up on .  There is one thing you said that leads
> me to believe there is something I don't understand: "...geocentric
> system...". My understanding is that there are two datatypes to choose from
> ; geometry and geography. I have chosen geometry because 1. it should be
> faster due to simpler calculations and 2. It should be faster because it can
> use the index.  Beyond that I thought that I would then need to choose a
> proper SRID. The two criteria for choosing an SRID for me are 1. the units
> are meters and 2. it is based on a datum that is appropriate for the USA.
> Now the fact that you said: "geocentric system" makes me think that maybe
> there is something else I need to be looking at when choosing a proper SRID.
> A little more on the application: The user wants to know all of the schools
> that are within a specific distance of his current location, his zip code.
> Precision is not important. We get the user's zip code from the screen and
> we have a list of schools and their zip codes.
>
>
> Thanks,
>
> Mike
>
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/postgis-users