Point in poly dataset mapping

[Seth Villanueva]

Hello,

I am trying to map points in one dataset to the polygons defined in a second dataset in such a way that the result is a dataset of all point with the polygon in which they are contained. I have attempted to use the command geoinpoly to accomplish this, but the command is failing to place virtually all points into any polygon. I have been unable to find any alternative commands, and so my question is if there is a known way to accomplish this in Stata or if it would be relatively easier to work in R? I will attach examples of each dataset below for reference.



Polygon data: (California County data)

* Example generated by -dataex-. To install: ssc install dataex
clear
input str2 STATEFP str3 COUNTYFP str11 INTPTLAT str12 INTPTLON int _ID double(_X _Y)
"06" "091" "+39.5769252" "-120.5219926"    9                   .                  .
"06" "067" "+38.4500161" "-121.3404408"  326                   .                  .
"06" "083" "+34.5370572" "-120.0399729"  330                   .                  .
"06" "009" "+38.1838996" "-120.5614415"  347                   .                  .
"06" "111" "+34.3587415" "-119.1331432"  395                   .                  .
"06" "037" "+34.1963983" "-118.2618616"  398                   .                  .
"06" "097" "+38.5251824" "-122.9261095"  456                   .                  .
"06" "031" "+36.0724780" "-119.8155301"  539                   .                  .
"06" "073" "+33.0236041" "-116.7761174"  540                   .                  .
"06" "061" "+39.0620323" "-120.7227181"  553                   .                  .
"06" "075" "+37.7272391" "-123.0322294"  637                   .                  .
"06" "041" "+38.0518169" "-122.7459738"  676                   .                  .
"06" "043" "+37.5700295" "-119.9128613"  687                   .                  .
"06" "035" "+40.7210862" "-120.6299518"  712                   .                  .
"06" "055" "+38.5070999" "-122.3259045"  744                   .                  .
"06" "089" "+40.7605142" "-122.0435556"  780                   .                  .
"06" "053" "+36.2401070" "-121.3155723"  804                   .                  .
"06" "105" "+40.6478033" "-123.1144760"  955                   .                  .
"06" "045" "+39.4323876" "-123.4428811"  966                   .                  .
"06" "027" "+36.5619770" "-117.4039269"  982                   .                  .
"06" "051" "+37.9158363" "-118.8751668" 1053                   .                  .
"06" "109" "+38.0214344" "-119.9647335" 1054                   .                  .
"06" "095" "+38.2672255" "-121.9395940" 1059                   .                  .
"06" "071" "+34.8572198" "-116.1811967" 1111                   .                  .
"06" "013" "+37.9194790" "-121.9515431" 1130                   .                  .
"06" "003" "+38.6217831" "-119.7983522" 1317                   .                  .
"06" "017" "+38.7855320" "-120.5343981" 1374                   .                  .
"06" "113" "+38.6795954" "-121.9024412" 1410                   .                  .
"06" "115" "+39.2662469" "-121.3407809" 1467                   .                  .
"06" "069" "+36.6107024" "-121.0852960" 1481                   .                  .
"06" "023" "+40.7066554" "-123.9261757" 1611                   .                  .
"06" "065" "+33.7298275" "-116.0022389" 1624                   .                  .
"06" "029" "+35.3466288" "-118.7295064" 1633                   .                  .
"06" "011" "+39.1777385" "-122.2375629" 1693                   .                  .
"06" "015" "+41.7639580" "-124.0036140" 1756                   .                  .
"06" "049" "+41.5929185" "-120.7183704" 1808                   .                  .
"06" "019" "+36.7610058" "-119.6550193" 1826                   .                  .
"06" "039" "+37.2098213" "-119.7498023" 1937                   .                  .
"06" "085" "+37.2207774" "-121.6906224" 1946                   .                  .
"06" "103" "+40.1261573" "-122.2322737" 1954                   .                  .
"06" "077" "+37.9349815" "-121.2722440" 2051                   .                  .
"06" "001" "+37.6471385" "-121.9124880" 2098                   .                  .
"06" "057" "+39.2975082" "-120.7713429" 2175                   .                  .
"06" "007" "+39.6659588" "-121.6019188" 2237                   .                  .
"06" "047" "+37.1948063" "-120.7228019" 2312                   .                  .
"06" "107" "+36.2288339" "-118.7810551" 2364                   .                  .
"06" "099" "+37.5623162" "-121.0028311" 2435                   .                  .
"06" "059" "+33.6756872" "-117.7772068" 2553                   .                  .
"06" "025" "+33.0408143" "-115.3554001" 2562                   .                  .
"06" "101" "+39.0361898" "-121.7039397" 2605                   .                  .
"06" "005" "+38.4435501" "-120.6538563" 2724                   .                  .
"06" "033" "+39.0948019" "-122.7467569" 2811                   .                  .
"06" "063" "+39.9922953" "-120.8243709" 2825                   .                  .
"06" "081" "+37.4146725" "-122.3715457" 2882                   .                  .
"06" "093" "+41.5879861" "-122.5332868" 2978                   .                  .
"06" "087" "+37.0124883" "-122.0072050" 3025                   .                  .
"06" "021" "+39.6025462" "-122.4016998" 3140                   .                  .
"06" "079" "+35.3852242" "-120.4475445" 3180                   .                  .
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655585           39.69356
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655524 39.693616999999996
"06" "091" "+39.5769252" "-120.5219926"    9 -120.65540899999999          39.693842
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655252          39.693906
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654887 39.693937999999996
"06" "091" "+39.5769252" "-120.5219926"    9          -120.65473          39.694115
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654708            39.6943
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654955          39.695003
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655035 39.695575999999996
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654888          39.695826
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654615 39.696228999999995
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654718          39.696479
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654473 39.697373999999996
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654794          39.697851
"06" "091" "+39.5769252" "-120.5219926"    9           -120.6554 39.699546999999995
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655419          39.699918
"06" "091" "+39.5769252" "-120.5219926"    9         -120.655207          39.700491
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654924          39.700821
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654453           39.70136
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654378          39.701626
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654563 39.702110999999995
"06" "091" "+39.5769252" "-120.5219926"    9 -120.65454199999999 39.702279999999995
"06" "091" "+39.5769252" "-120.5219926"    9          -120.65411          39.703094
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654098          39.703441
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654281 39.704288999999996
"06" "091" "+39.5769252" "-120.5219926"    9          -120.65453          39.704515
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654818          39.705041
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654901          39.705404
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654628          39.705694
"06" "091" "+39.5769252" "-120.5219926"    9 -120.65436700000001          39.705822
"06" "091" "+39.5769252" "-120.5219926"    9         -120.654124 39.706233999999995
"06" "091" "+39.5769252" "-120.5219926"    9 -120.65422699999999          39.706629
"06" "091" "+39.5769252" "-120.5219926"    9         -120.646465            39.7066
"06" "091" "+39.5769252" "-120.5219926"    9         -120.646162          39.706607
"06" "091" "+39.5769252" "-120.5219926"    9         -120.645708 39.706596999999995
"06" "091" "+39.5769252" "-120.5219926"    9         -120.643744          39.706589
"06" "091" "+39.5769252" "-120.5219926"    9         -120.642713          39.706575
"06" "091" "+39.5769252" "-120.5219926"    9 -120.64250200000001          39.706585
"06" "091" "+39.5769252" "-120.5219926"    9         -120.636785          39.706565
"06" "091" "+39.5769252" "-120.5219926"    9         -120.636731          39.706565
"06" "091" "+39.5769252" "-120.5219926"    9         -120.626064          39.706526
"06" "091" "+39.5769252" "-120.5219926"    9         -120.619899          39.706542
end

Point Data: (CA Wind Turbine locations)

* Example generated by -dataex-. To install: ssc install dataex
clear
input double(_X _Y) long ID byte WPC
-124.19923070720745  41.99504410912535   1 2
-124.15194436766573 41.996515579862084   2 2
-124.13626966740948  41.99659619708002   3 2
-124.13267554213068  41.99682075496593   4 4
-124.13027943019631  41.99697038591703   5 5
-124.13247492619587 41.995033388336076   6 2
-124.13007888025852  41.99518300996258   7 3
-124.12648476875978  41.99540733094536   8 4
-124.12368872850877  41.99414184406832   9 3
 -124.0238445627746 41.994566112181914  10 2
-124.01665472944521  41.99500802096399  11 2
 -124.0106630367236  41.99537584461126  12 1
-124.02125048096023  41.99292566574144  13 2
-124.00054545030719  41.99479652422508  14 3
 -123.9888956072007  41.99490992597796  15 3
-124.09044341852564  41.98595128053652  16 2
-123.92258802611877 41.995347661932705  17 3
 -123.9201908848489  41.99549290238514  18 4
 -124.1800354509185  41.98422734170503  19 2
-124.17159167955208 41.978184525656836  20 3
-124.16919647079071  41.97833491505535  21 2
-124.08774855395916 41.983418670437295  22 2
-123.93497482546229  41.99387633679682  23 1
-123.92239330025338  41.99355948510758  24 2
-123.91846301999209  41.99439758206417  25 3
-123.90321620545144  41.99471994599033  26 3
-124.17139010638246  41.97639730587627  27 2
-124.16899496353301  41.97654768590305  28 3
-124.16540220647177 41.976773144590055  29 2
-123.91512421317672  41.99327974790947  30 2
 -123.8526781905957 41.995962593932774  31 2
-124.08505392622261  41.98088600379725  32 2
-124.01558359673137 41.989674116846444  33 3
-123.85008796574716 41.994317926344095  34 2
-124.16436147205165  41.97395862063676  35 3
 -124.1051813688967  41.98971215857214  36 2
-123.90151357520836  41.99173717137083  37 2
-123.82112528621323 41.994250770926854  38 3
-124.08142824861955  41.97721102939065  39 2
-123.86787525697741  41.99145004658092  40 3
-123.74899948361264  41.99673387875674  41 1
-123.88057614005508  41.98618583761484  42 3
-123.80616496066004  41.98974289272031  43 3
-123.76352453306922  41.99263488309432  44 2
-123.75753874154594 41.993829126324485  45 3
-123.75581344599868   41.9927311221101  46 4
-123.75341556488739  41.99287287209101  47 3
-123.70545433389998 41.995697432551175  48 2
-124.13785041433972 41.967699133213586  49 2
-124.08555489725545   41.9788095133928  50 3
 -123.7556234069209  41.99094231882569  51 3
-123.75322559200758   41.9910840599659  52 4
-123.75216683876506  41.99474666787404  53 2
-124.14166485294949  41.96826027329252  54 3
-123.93861036728332 41.983991321535484  55 3
-123.75303563635404  41.98929524658054  56 2
-123.74344446934528  41.98986167879306  57 2
-123.66190845871776  41.99464418922891  58 3
 -124.0440363300365 41.969923535748876  59 2
-123.81406484250962  41.98747207726343  60 1
-123.79143994901887  41.99089610352341  61 2
-123.74924913375074  41.98771890540841  62 2
-123.65692434179236   41.9931346099319  63 2
 -123.6545259582115  41.99327429508221  64 3
-123.65319922601188  41.99443154542198  65 4
-123.64825489299496  41.99483934024041  66 3
-124.10252250078325  41.96450182826098  67 3
-124.04024513452715  41.96835735847832  68 2
-123.88449894945937 41.977849548074694  69 2
-123.88395813742933  41.98028211144183  70 3
-123.85639414379649  41.98193962568957  71 2
-123.78898271558198 41.987001496567935  72 3
-123.74306513921164  41.98628397758737  73 2
-123.74574752106672   41.9888256785592  74 3
-123.70110286298826  41.98875244100158  75 2
-123.67352480042507  41.99036627827407  76 3
 -123.6702863665122  41.99235540055622  77 2
-123.50201108674051  42.00025266527839  78 2
 -123.8764255812738  41.98695671910493  79 2
-123.78358290291071  41.98538584032165  80 2
-123.77763928255803   41.9824383978933  81 2
-123.77404320398804  41.98265165815697  82 3
-123.77257191084425  41.98393889089112  83 2
-123.73084499428728  41.98784435768682  84 2
-123.64526960810396  41.99321294741084  85 2
-123.45263846712704   42.0012503975518  86 2
-124.15632706435248  41.98057849903591  87 1
-124.09134610781233 41.961596436558324  88 2
-124.07241283164576  41.97361932296619  89 2
 -124.0471543653671  41.96541140858483  90 2
-123.94661110938358 41.972289739930645  91 2
-123.85457279006548  41.97836324578878  92 3
 -123.8506549546869  41.97718391419475  93 2
-123.77145537565835   41.9810050137775  94 3
-123.73549311630308 41.983130467796094  95 2
-123.68535363721296  41.99067469747613  96 2
-123.51984541121412  41.99766163217692  97 2
-124.49070010063689  41.95307624098682  98 7
-124.06839423016922  41.96121916946119  99 2
-124.04204700283034 41.962846507899485 100 2
end

Please let me know if there is anymore information I could provide that would be useful to know, and thank you very much in advance for any and all help.

[Robert Picard]

You do not show what code you tried so I'm guessing here that you tried to match points to the county centroids instead of the polygons.

I copied the second dataex example and saved the data in "points.dta". I found a US-wide county shapefile from this page:

https://catalog.data.gov/dataset/tig...onal-shapefile

and downloaded the shapefile using this link:

https://www2.census.gov/geo/tiger/TI..._us_county.zip

Once unzipped to the current directory, I used the following code to convert the shapefile to Stata datasets and then match each point to a county:

* downloaded from https://www2.census.gov/geo/tiger/TIGER2017/COUNTY/tl_2017_us_county.zip
shp2dta using "tl_2017_us_county.shp", data("shp.dta")  coor("shp_coor.dta") replace

* load second dataex example
use "points.dta", clear
rename ID point_id

* match lat/lon to polygons in the county shapefile
geoinpoly _Y _X using "shp_coor.dta"

* recover county data by merging on the matched shapefile identifier
merge m:1 _ID using "shp.dta", keep(master match) nogen keepusing(STATEFP COUNTYFP NAME)

* show 5 random examples, -listsome- if from SSC
sort point_id
set seed 3123
listsome, max(5) random

Here's where 5 random points were matched:

. listsome, max(5) random

     +---------------------------------------------------------------------------------+
     |         _X          _Y   point_id   WPC    _ID   STATEFP   COUNTYFP        NAME |
     |---------------------------------------------------------------------------------|
 11. | -124.01665   41.995008         11     2   1756        06        015   Del Norte |
 46. | -123.75581   41.992731         46     4   1756        06        015   Del Norte |
 56. | -123.75304   41.989295         56     2   1756        06        015   Del Norte |
 72. | -123.78898   41.987001         72     3   1756        06        015   Del Norte |
 99. | -124.06839   41.961219         99     2   1756        06        015   Del Norte |
     +---------------------------------------------------------------------------------+

.

[Seth Villanueva]

Hey Robert, thanks for the response.
I think part of my problem was in the sorting of my polygon data, with all of a polygon's points needing to appear together with their empty coordinate entry coming first.

Here is my data and code that I failed to provide before:
county data: https://www.census.gov/geographies/m...file.2017.html
wind data: https://www.nrel.gov/gis/data-wind.html

specifically I am working with CA wind data

//All counties spatial data
shp2dta using tl_2017_us_county.shp, database(shp2dta_counties_database) coordinates(shp2dta_counties_coord) replace
use shp2dta_counties_database.dta
merge 1:m _ID using shp2dta_counties_coord.dta
drop LSAD CLASSFP MTFCC CSAFP CBSAFP METDIVFP FUNCSTAT _merge
save counties_GIS.dta, replace


//California

spshape2dta california_50m_wind, replace
use california_50m_wind.dta
rename _CX _X
rename _CY _Y
rename _ID point_id
save california_wind_GIS.dta, replace
//CA only counties
use counties_GIS
keep if STATEFP == "06"
//getting empty coordinate pairs to the start of each polygon data block
gen empty = 1
replace empty = 0 if _X == . & _Y == .
sort _ID empty
drop empty
save CA_counties_GIS.dta, replace

//placing wind points into county polygons
use california_wind_GIS
geoinpoly _Y _X using "CA_counties_GIS.dta"

I am choosing to use the spshape2dta command rather than shp2dta, however I tried to implement your method using shp2dta as well and got the same result.

Following the execution of this code I get the message: "(note: dropped 65 shapes that are not polygons)"
And the data then looks like such:

* Example generated by -dataex-. To install: ssc install dataex
clear
input long point_id double(_X _Y) long ID byte WPC long _ID
  1 -124.19923070720745  41.99504410912535   1 2 .
  2 -124.15194436766573 41.996515579862084   2 2 .
  3 -124.13626966740948  41.99659619708002   3 2 .
  4 -124.13267554213068  41.99682075496593   4 4 .
  5 -124.13027943019631  41.99697038591703   5 5 .
  6 -124.13247492619587 41.995033388336076   6 2 .
  7 -124.13007888025852  41.99518300996258   7 3 .
  8 -124.12648476875978  41.99540733094536   8 4 .
  9 -124.12368872850877  41.99414184406832   9 3 .
 10  -124.0238445627746 41.994566112181914  10 2 .
 11 -124.01665472944521  41.99500802096399  11 2 .
 12  -124.0106630367236  41.99537584461126  12 1 .
 13 -124.02125048096023  41.99292566574144  13 2 .
 14 -124.00054545030719  41.99479652422508  14 3 .
 15  -123.9888956072007  41.99490992597796  15 3 .
 16 -124.09044341852564  41.98595128053652  16 2 .
 17 -123.92258802611877 41.995347661932705  17 3 .
 18  -123.9201908848489  41.99549290238514  18 4 .
 19  -124.1800354509185  41.98422734170503  19 2 .
 20 -124.17159167955208 41.978184525656836  20 3 .
 21 -124.16919647079071  41.97833491505535  21 2 .
 22 -124.08774855395916 41.983418670437295  22 2 .
 23 -123.93497482546229  41.99387633679682  23 1 .
 24 -123.92239330025338  41.99355948510758  24 2 .
 25 -123.91846301999209  41.99439758206417  25 3 .
 26 -123.90321620545144  41.99471994599033  26 3 .
 27 -124.17139010638246  41.97639730587627  27 2 .
 28 -124.16899496353301  41.97654768590305  28 3 .
 29 -124.16540220647177 41.976773144590055  29 2 .
 30 -123.91512421317672  41.99327974790947  30 2 .
 31  -123.8526781905957 41.995962593932774  31 2 .
 32 -124.08505392622261  41.98088600379725  32 2 .
 33 -124.01558359673137 41.989674116846444  33 3 .
 34 -123.85008796574716 41.994317926344095  34 2 .
 35 -124.16436147205165  41.97395862063676  35 3 .
 36  -124.1051813688967  41.98971215857214  36 2 .
 37 -123.90151357520836  41.99173717137083  37 2 .
 38 -123.82112528621323 41.994250770926854  38 3 .
 39 -124.08142824861955  41.97721102939065  39 2 .
 40 -123.86787525697741  41.99145004658092  40 3 .
 41 -123.74899948361264  41.99673387875674  41 1 .
 42 -123.88057614005508  41.98618583761484  42 3 .
 43 -123.80616496066004  41.98974289272031  43 3 .
 44 -123.76352453306922  41.99263488309432  44 2 .
 45 -123.75753874154594 41.993829126324485  45 3 .
 46 -123.75581344599868   41.9927311221101  46 4 .
 47 -123.75341556488739  41.99287287209101  47 3 .
 48 -123.70545433389998 41.995697432551175  48 2 .
 49 -124.13785041433972 41.967699133213586  49 2 .
 50 -124.08555489725545   41.9788095133928  50 3 .
 51  -123.7556234069209  41.99094231882569  51 3 .
 52 -123.75322559200758   41.9910840599659  52 4 .
 53 -123.75216683876506  41.99474666787404  53 2 .
 54 -124.14166485294949  41.96826027329252  54 3 .
 55 -123.93861036728332 41.983991321535484  55 3 .
 56 -123.75303563635404  41.98929524658054  56 2 .
 57 -123.74344446934528  41.98986167879306  57 2 .
 58 -123.66190845871776  41.99464418922891  58 3 .
 59  -124.0440363300365 41.969923535748876  59 2 .
 60 -123.81406484250962  41.98747207726343  60 1 .
 61 -123.79143994901887  41.99089610352341  61 2 .
 62 -123.74924913375074  41.98771890540841  62 2 .
 63 -123.65692434179236   41.9931346099319  63 2 .
 64  -123.6545259582115  41.99327429508221  64 3 .
 65 -123.65319922601188  41.99443154542198  65 4 .
 66 -123.64825489299496  41.99483934024041  66 3 .
 67 -124.10252250078325  41.96450182826098  67 3 .
 68 -124.04024513452715  41.96835735847832  68 2 .
 69 -123.88449894945937 41.977849548074694  69 2 .
 70 -123.88395813742933  41.98028211144183  70 3 .
 71 -123.85639414379649  41.98193962568957  71 2 .
 72 -123.78898271558198 41.987001496567935  72 3 .
 73 -123.74306513921164  41.98628397758737  73 2 .
 74 -123.74574752106672   41.9888256785592  74 3 .
 75 -123.70110286298826  41.98875244100158  75 2 .
 76 -123.67352480042507  41.99036627827407  76 3 .
 77  -123.6702863665122  41.99235540055622  77 2 .
 78 -123.50201108674051  42.00025266527839  78 2 .
 79  -123.8764255812738  41.98695671910493  79 2 .
 80 -123.78358290291071  41.98538584032165  80 2 .
 81 -123.77763928255803   41.9824383978933  81 2 .
 82 -123.77404320398804  41.98265165815697  82 3 .
 83 -123.77257191084425  41.98393889089112  83 2 .
 84 -123.73084499428728  41.98784435768682  84 2 .
 85 -123.64526960810396  41.99321294741084  85 2 .
 86 -123.45263846712704   42.0012503975518  86 2 .
 87 -124.15632706435248  41.98057849903591  87 1 .
 88 -124.09134610781233 41.961596436558324  88 2 .
 89 -124.07241283164576  41.97361932296619  89 2 .
 90  -124.0471543653671  41.96541140858483  90 2 .
 91 -123.94661110938358 41.972289739930645  91 2 .
 92 -123.85457279006548  41.97836324578878  92 3 .
 93  -123.8506549546869  41.97718391419475  93 2 .
 94 -123.77145537565835   41.9810050137775  94 3 .
 95 -123.73549311630308 41.983130467796094  95 2 .
 96 -123.68535363721296  41.99067469747613  96 2 .
 97 -123.51984541121412  41.99766163217692  97 2 .
 98 -124.49070010063689  41.95307624098682  98 7 .
 99 -124.06839423016922  41.96121916946119  99 2 .
100 -124.04204700283034 41.962846507899485 100 2 .
end

With the obvious problem being that _ID has no observations.

I'm not sure where I'm going wrong, it seems like my method properly mirrors your own.
Thanks again for the help.

[Robert Picard]

A shapefile has two components, the database dataset with one observation per shape and the coordinates dataset where you have a series of points that, when connected, create each polygon shape. The order of the coordinates matters. Since you used merge on the county coordinates, you corrupted the shapefile since the data is reordered to perform the merge.

The California wind data is presented by 200-1000m cells and the shapefile coordinates contain the vertices of the polygon that describe each cell. The spshape2dta command automatically includes shape centroids in the shapefile database. So all you have to do is to find the county where each cell centroid is located:

* load wind data centroids (200-1000m cell area) 
use "california_50m_wind.dta", clear
assert _ID == ID
drop _ID
rename ID cell_id

* match each cell centroid to a county polygon
geoinpoly _CY _CX using "shp2dta_counties_coord.dta"
gen match = !mi(_ID)
tab WPC match

* recover the state/county info from the county database
merge m:1 _ID using "shp2dta_counties_database.dta", keep(master match) keepusing(STATEFP COUNTYFP NAME) nogen
isid cell_id, sort

* exclude centroids located in other states
list if STATEFP != "06" & match
replace match = 0 if STATEFP != "06"
tab COUNTYFP WPC if match

[Seth Villanueva]

A shapefile has two components, the database dataset with one observation per shape and the coordinates dataset where you have a series of points that, when connected, create each polygon shape. The order of the coordinates matters. Since you used merge on the county coordinates, you corrupted the shapefile since the data is reordered to perform the merge.

The California wind data is presented by 200-1000m cells and the shapefile coordinates contain the vertices of the polygon that describe each cell. The spshape2dta command automatically includes shape centroids in the shapefile database. So all you have to do is to find the county where each cell centroid is located:

* load wind data centroids (200-1000m cell area)
use "california_50m_wind.dta", clear
assert _ID == ID
drop _ID
rename ID cell_id

* match each cell centroid to a county polygon
geoinpoly _CY _CX using "shp2dta_counties_coord.dta"
gen match = !mi(_ID)
tab WPC match

* recover the state/county info from the county database
merge m:1 _ID using "shp2dta_counties_database.dta", keep(master match) keepusing(STATEFP COUNTYFP NAME) nogen
isid cell_id, sort

* exclude centroids located in other states
list if STATEFP != "06" & match
replace match = 0 if STATEFP != "06"
tab COUNTYFP WPC if match

Thank you Robert, this works like a charm. Thank you very much for pointing out where I was mistaken and also for providing a correction.

[Paola Mallia]

Hi Robert,
This post is extremely useful.
Following all the instructions above, I still get this error code:

Y (latitude) must be between -90 and 90 in coordinates file
r(198);

That I do not get if use the option "noproj" in the geoinpoly command, but in that case I get zero matches (_ID all missing, like Seth).
How to overcome this issue?

thanks