2

I am trying to do a spatial join on two tables, the join runs successfully but there's no data in the new table. what am i doing wrong?

CREATE TABLE us_countyparts_2008_with_count_v5 AS(
  SELECT cp.statefp, cp.countyfp, cp.cd110fp, sum(bd.pop00)
  FROM countyparts2008_with_cd AS cp
  JOIN us_block_2000_with_data AS bd
    ON ST_Contains(cp.geom, bd.geom)
  GROUP BY cp.statefp, cp.countyfp, cp.cd110fp
)

I did IsValid on both and they both returned as valid geometries

This is some of the data from countyparts enter image description here

and this is some data from us_blocks enter image description here

I'm pretty sure the error is something to do with the extent this is the extent for the us_block file

"BOX(-2356113.74289801 -783964.164622752,1044933.24146018 845925.228214008)"

and this is the extent for the county file "BOX(-124.041983 24.3963080000001,-69.3723069999999 48.540437)"

If this is the problem then how do I fix the extent for the us_block file

  • Its difficult to say without detail on the schema. What does SELECT *, ST_AsEWKT(geom) FROM countyparts2008_with_cd LIMIT 2 show? What does SELECT *, ST_AsEWKT(geom) FROM us_block_2000_with_data LIMIT 2 show? Do you possibly have any invalid geometries in your data? – BradHards Oct 29 '16 at 5:30
  • 1
    And if you run only inner select it returns expected data? It seems good to me. Maybe you should try another spatial operator (ST_Intersects)? – DavidP Oct 29 '16 at 5:36
  • @BradHards they both return multipolygons – captnvitman Oct 29 '16 at 13:47
  • and I've tried intersects and a few others and they all return an empty table @DavidP – captnvitman Oct 29 '16 at 13:47
  • 2
    Please show the data - edit the question (click edit below the question) to show the top few rows. Also, please tell us if you checked for invalid geometries. We can't really help you if you don't help us. – BradHards Oct 29 '16 at 22:51
2

The query you created should work without any issues, as shown by doing the following simple query which uses known geometries.

-- Use a CTE to create some known data
WITH countyparts2008_with_cd AS (
  SELECT *
  FROM (VALUES
    (1,1,1,'srid=4326;MULTIPOLYGON(((0 0, 10 0, 10 10, 0 10, 0 0)))'::Geometry)
    )A(statefp, countyfp, cd110fp, geom)
  ),
  us_block_2000_with_data AS (
  SELECT *
  FROM (VALUES
    (10,'srid=4326;MULTIPOLYGON(((0 0, 5 0, 5 5, 0 5, 0 0)))'::Geometry)
    ,(10,'srid=4326;MULTIPOLYGON(((5 5, 10 5, 10 10, 5 10, 5 5)))'::Geometry)
    )A(pop00,geom)
  )
-- The select statement from the create table query
  SELECT cp.statefp, cp.countyfp, cp.cd110fp, sum(bd.pop00)
  FROM countyparts2008_with_cd AS cp
  JOIN us_block_2000_with_data AS bd
    ON ST_Contains(cp.geom, bd.geom)
  GROUP BY cp.statefp, cp.countyfp, cp.cd110fp

The issue you are having is that the us_block table geometries are not in the same projection as the country file, which is in WGS84. I do not know which projection the us_block table is in.

There are a couple of questions that cover reprojecting the geometries in the table, here and here. You will need to determine what the projection is.

You could also use ST_Transform within your query. You will need to get the SRID set correctly against the geometries in the us_block table.

CREATE TABLE us_countyparts_2008_with_count_v5 AS(
  SELECT cp.statefp, cp.countyfp, cp.cd110fp, sum(bd.pop00)
  FROM countyparts2008_with_cd AS cp
  JOIN us_block_2000_with_data AS bd
    ON ST_Contains(cp.geom, ST_Transform(bd.geom,4326))
  GROUP BY cp.statefp, cp.countyfp, cp.cd110fp
)
  • Now that I've looked at it more I think it has something to do with the extent. just edited my original post with more info – captnvitman Oct 31 '16 at 2:00
  • @captnvitman I have altered the answer to suit. You will need to determine what the projection the us_block table is in and alter the SRID in the table to suit. – MickyT Oct 31 '16 at 8:54
1

Did you try to join the tables with a CROSS JOIN and a WHERE?

Something like:

CREATE TABLE us_countyparts_2008_with_count_v5 AS(
  SELECT cp.statefp, cp.countyfp, cp.cd110fp, sum(bd.pop00)
  FROM countyparts2008_with_cd AS cp CROSS JOIN us_block_2000_with_data AS bd
  WHERE ST_CONTAINS(cp.geom, bd.geom)
  GROUP BY cp.statefp, cp.countyfp, cp.cd110fp
)

I learned it that way in school. I have never seen, in my limited experience with PostGIS, spatial joins the way you do it. This syntax is very interesting though and I will try it. CROSS JOIN operations are very time consuming for the server. Specially when you cross join hundred of thousands of rows from each relation!!

  • In general this is correct, but when you look on the table of the author the column cd110fp might not be part of the group by. The first 10 rows indicate that the rows and therfore the geometries might be a product of state and county only. But only the author might tell if this is the case. – Matte Oct 30 '16 at 17:57
  • so what im doing here is finding out what counties are split by congressional district boundaries then using census blocks to find out the populations of each of those county splits, so i need to group by statefp, county fips, and cd number – captnvitman Oct 30 '16 at 18:13
  • I would think it is something to do with the census blocks, because a similar spatial join works when i am working with just counties and congressional districts, but there's nothing to indicate that there is something wrong with the shapefiles I imported – captnvitman Oct 30 '16 at 18:27
0

I would try to change the spatial predicate to ST_Intersects (just for sure). ST_Contains(cp.geom, ST_Buffer(bd.geom, -0.01)) might also help, I suppose the two datasets don't match on the edges.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.