6

R での O'Reilly の Data Mashupsをインスピレーションとして使用して、ユタ州ソルトレイク郡のシェープファイルにいくつかの住所をプロットしようとしています。

私はデータフレームのジオテーブルを持っています:

> geoTable
         address        Y         X EID
1    130 E 300 S 40.76271 -111.8872   1
2    875 E 900 S 40.74992 -111.8660   2
3   2200 S 700 E 40.72298 -111.8714   3
4    702 E 100 S 40.76705 -111.8707   4
5 177 East 200 S 40.76518 -111.8859   5
6    702 3rd ave 40.77264 -111.8683   6
7   2175 S 900 E 40.72372 -111.8652   7
8   803 E 2100 S 40.72556 -111.8680   8

そして、それを eventData オブジェクトに強制しました:

> addressEvents<-as.EventData(geoTable,projection=NA)
> addressEvents
         address        Y         X EID
1    130 E 300 S 40.76271 -111.8872   1
2    875 E 900 S 40.74992 -111.8660   2
3   2200 S 700 E 40.72298 -111.8714   3
4    702 E 100 S 40.76705 -111.8707   4
5 177 East 200 S 40.76518 -111.8859   5
6    702 3rd ave 40.77264 -111.8683   6
7   2175 S 900 E 40.72372 -111.8652   7
8   803 E 2100 S 40.72556 -111.8680   8

プロットに必要なものはすべて揃っているように見えますが、機能していません。シェープファイルをロードしてプロットすると

addPoints(addressEvents,col="red",cex=.5)

空のシェープファイルを見たままです。さらに、eventData オブジェクトに対して findPolys を実行しようとすると、NULL が返されます。

> findPolys(addressEvents,myShapeFile)
NULL

どうすればこれを機能させることができますか？O'Reilly のチュートリアルを問題なく完了することができましたが、ここでどこが間違っているのかを理解するのに苦労しています。それがシェープファイルなのか、データ フレームなのか、その他のものなのかわかりません。

データとシェープファイルをインポートするために使用するコマンドは次のとおりです

slc<-read.table('~/utah.txt',sep=',',header=TRUE,strip.white=TRUE,stringsAsFactors=FALSE)

myShapeFile<-importShapefile("/Users/neil/Downloads/SGID93_DEMOGRAPHIC_CensusTracts2000/SGID93_DEMOGRAPHIC_CensusTracts2000",readDBF=TRUE)
4

2 に答える 2

5

これらの関連する質問、特に Eduardo の回答も参照してください。

于 2009-09-26T19:03:01.707 に答える
4

PBSmapping は、粗雑なヒューリスティックを使用して .prj ファイルから投影を計算しているようです。(ヘルプ (importShapefile) を参照)。私は個人的にprjファイル内のすべてのものを理解していませんが、このWebサイトwww.spatialreference.orgを使用すると、マップが一致すると思います

http://www.spatialreference.org/ref/epsg/26912/

新しい形状ファイルを取得するたびに、この Web サイトでその投影システムを見つけて、proj4 文字列を探します。この場合は、「+proj=utm +zone=12 +ellps=GRS80 +datum=NAD83 +units=m +」です。 no_defs"

(私が言ったように、私は PBSmapping を知りませんが、次のように maptools を使用してこれを読むことができます)

library(maptools)
sf=readShapeSpatial("SGID93_DEMOGRAPHIC_CensusTracts2000.shp",proj4string=CRS("+proj=utm +zone=12 +ellps=GRS80 +datum=NAD83 +units=m +no_defs"))

次に、次を使用してlatlongに変換します

library(rgdal)

sftransformed=spTransform(sf,CRS("+proj=longlat"))

と

plot(sftransformed,axes=T)

軸上に適切な単位を含むプロットを提供します。

PBSmapping が proj4 文字列を理解するかどうかわからない? 正直いらないらしい。

于 2009-09-26T13:27:33.177 に答える