地図のことば

2013年2月15日 (金)

『基準点成果等閲覧サービス』(国土地理院ホーム)のサイト

■「国土地理院が設置した各種基準点の成果等閲覧サービス」のサイト
・『基準点成果等閲覧サービス
http://sokuseikagis1.gsi.go.jp/

●閲覧例の概略メモ
1) 配点図から閲覧対象を検索
Haiten_zu_web

2) 基準点検索
Kijyunten_kensaku_b_web

3) 基準点一覧の表示
Kijyunten_jyouho
赤枠部をクリックすると
以下の成果情報が表示されるようです.

4) 基準点基本情報/基準点成果情報の表示

Seika_j_b

・閲覧概要
基準点名  折木沢
三等三角点 TR35240608801
北緯  35°14′11″.6438
東経 140°06′28″.3583
標高(m) 190.00
ジオイド高(m)  34.28

基準点成果等閲覧サービス閲覧機能ユーザマニュアル
http://sokuseikagis1.gsi.go.jp/SysMsg/help/manual.pdf

●参考サイト
標高がわかるWeb地図(右クリックで標高が表示されます)
http://saigai.gsi.go.jp/2012demwork/checkheight/index.html
本システムで得られる標高値について
http://saigai.gsi.go.jp/2012demwork/checkheight/readme.pdf

直接,緯度経度を入力して標高を表示する仕様にはなっていないようです.
できるとよいと思いました.
標高結果はCSVファイルとして出力可能.

【折木沢周辺の検索例】
35.236512    140.108222    180    10m(種別)
標高精度は5.0m以内のようです.

●参考サイト 
標高精度は各サイトで確認しましょう.
位置情報データベース
http://www.hucc.hokudai.ac.jp/~x10795/
携帯版位置情報データベース
http://www.hucc.hokudai.ac.jp/~x10795/kindex.html
緯度・経度・標高・水深を求める(2010/6/12)
http://www.hucc.hokudai.ac.jp/~x10795/Latlonele.html
緯度・経度・標高・水深リストの作成(2010/7/6)
http://www.hucc.hokudai.ac.jp/~x10795/Latlonele2.html

「緯度・経度・標高・水深リストの作成(2010/7/6)」で
上記基準点(TR35240608801)の緯度経度を計算した結果は160mとなるようです.

LatLng2Height:::緯度経度から標高算出API
http://lab.uribou.net/ll2h/


2013年1月29日 (火)

『UTM座標値と経度緯度の換算』などに関するサイト

■『UTM座標値と経度緯度の換算』などに関するサイトを探してみました.

以下に,学習に役立つ・利便性の高いと思われるサイトをメモさせて頂きました.
サイトを公開されている方々・関係機関に感謝いたします.
他にもたくさんあると思います.
貴重なサイトから多くの情報と学習の支援を頂いています.

現地測量で得た測点の地磁気補正と測地系変換
http://www2.ipcku.kansai-u.ac.jp/~moto/GISContents/47.htm#5
情報を得るのに役立ちました.感謝

How to Use the Spreadsheet for Converting UTM to Latitude and Longitude (Or Vice Versa)
http://www.uwgb.edu/dutchs/UsefulData/HowUseExcel.HTM

UTM
http://homepage3.nifty.com/Nowral/21_UTM/21_UTM.html
Utm2

平面直角座標 ⇔ 緯度・経度 変換
http://www.n-survey.com/online/xybl.htm

MSP GEOTRANS 3.2 (Geographic Translator)
http://earth-info.nga.mil/GandG/geotrans/

Transverse Mercator Calculator
http://www.dmap.co.uk/ll2tm.htm
Tmc

緯度経度座標をUTM座標に変換する
http://asp.ncm-git.co.jp/QuickConvert/BL2UTM.aspx

Convert Coordinates - Calculate a position in a variety of formats.
http://www.earthpoint.us/Convert.aspx
Google Earthに各種グリッドを重ねて表示してくれます.
学習支援ソフトとしてもすばらしい.
このようなソフトが利用できるとは.感謝!

Convert_coord

Earth_point_b

■参考サイト
Universal Transverse Mercator (UTM), the Military Grid Reference System (MGRS), and the Universal Polar Stereographic (UPS)
http://earth-info.nga.mil/GandG/coordsys/grids/universal_grid_system.html#zz8

Reference Guide to the Universal Grid System (1.19 MB)
http://earth-info.nga.mil/GandG/coordsys/images/utm_mgrs_images/universal_grid_basics_20070319.pdf

・GRS Grid Zone Designators (GZD) Graphic  (4.72 MB)
http://earth-info.nga.mil/GandG/coordsys/images/utm_mgrs_images/MGRS_GZD.pdf

■60進数⇔10進数の変換コード例(R) 不具合は修正してください
# 引数を文字列で与える."129:30:35.187" 分離記号は : - /の3種類
Scale60to10<- function(sexagesimal60){
  sep<- c(":","-","/")
  nc<- nchar(sexagesimal60)
  tmp<- strsplit(sexagesimal60,"")[[1]]
  n<- length(sep)
  for (i in 1:nc){
    for (j in 1:n) {
      if (tmp[i]==sep[j] ){
        kugiri<- sep[j]
        break
      }
    }
  }
  Tmp<- strsplit(sexagesimal60,kugiri)[[1]]
  Out<- as.numeric(Tmp[1])+as.numeric(Tmp[2])/60+as.numeric(Tmp[3])/60/60
  return(Out)
}

#10進数を60進数に変換
Scale10to60<- function(decimal10){
  sign<- sign(decimal10)
  decimal10<- abs(decimal10)
  sc1<- trunc(decimal10)
  tmp1<- decimal10-sc1
  sc2<- trunc(tmp1*60)
  tmp2<- tmp1*60-sc2
  sc3<- tmp2*60
  sc1<- sign*sc1
  #Tmp<-paste(sc1,"°",sc2,"′",sc3,"″",sep="")
  Tmp<-paste(sc1,"/",sc2,"/",sc3,sep="")
  return(noquote(Tmp))#noquoteで返す
}

●使用例
木更津高専の緯度/経度  35.384502/139.955692
を例にして

>  Scale10to60(35.384502)
[1] 35/23/4.20719999999164
>  Scale10to60(139.955692)
[1] 139/57/20.4911999999968

>  Scale60to10("35/23/4.20719999999164")
[1] 35.3845
>  Scale60to10("139/57/20.4911999999968")
[1] 139.9557

■Rで座標系の変換の例(不具合は修正してください)
Rで座標系の変換
http://spatiohack.blogspot.jp/2010/11/r_15.html
PBSmapping: Mapping Fisheries Data and Spatial Analysis Tools
http://cran.r-project.org/web/packages/PBSmapping/index.html
PBSmapping.pdf
http://cran.r-project.org/web/packages/PBSmapping/PBSmapping.pdf
Introduction to PBSmapping
http://cran.r-project.org/web/packages/PBSmapping/vignettes/PBSmappingIntro.pdf
Package ‘rgdal’January 18, 2013
http://cran.r-project.org/web/packages/rgdal/rgdal.pdf
spTransform longlat to utm
http://comments.gmane.org/gmane.comp.lang.r.general/277732

●緯度経度をUTM座標に変換
>  library(rgdal)
>  xy <- cbind(c(139.955692), c(35.384502))#複数の時はc(139.955,**,**), c(35.384,**,**)
>  project(xy, "+proj=utm +zone=54 ellps=WGS84")
       [,1]    [,2]
[1,] 405151 3916185
>  utm=project(xy, "+proj=utm +zone=54 ellps=WGS84")
>  cat(format.default(utm,nsmall=2),"\n")
405151.34 3916184.79

●UTM座標を緯度経度に変換
>  xy <- cbind(c(405151.34), c(3916184.79))
>  project(xy, "+proj=utm +zone=54 ellps=WGS84",inv=TRUE)
         [,1]    [,2]
[1,] 139.9557 35.3845

●緯度経度を平面直角座標に変換
>  pj <- CRS("+init=epsg:2451")
>  sp.xy <- SpatialPoints(matrix(c(139.955692,35.384502),nrow=1,byrow=TRUE),
  proj4string = CRS("+proj=longlat +datum=WGS84"))
>  spTransform(sp.xy, CRS=pj)
SpatialPoints:
     coords.x1 coords.x2
[1,]  11116.29 -68277.84
Coordinate Reference System (CRS) arguments: +init=epsg:2451 +proj=tmerc
+lat_0=36 +lon_0=139.8333333333333 +k=0.9999 +x_0=0 +y_0=0 +ellps=GRS80
+towgs84=0,0,0,0,0,0,0 +units=m +no_defs

●平面直角座標を緯度経度に変換
【例1】
>  library(rgdal)
>  pj <- CRS("+init=epsg:2451")
>  sp.xy <- SpatialPoints(matrix(c(11116.29,-68277.84),nrow=1,byrow=TRUE),proj4string = pj)
>  pj2 <- CRS("+proj=longlat +datum=WGS84")
>  spTransform(sp.xy, CRS=pj2)
SpatialPoints:
     coords.x1 coords.x2
[1,]  139.9557   35.3845
Coordinate Reference System (CRS) arguments: +proj=longlat
+datum=WGS84 +ellps=WGS84 +towgs84=0,0,0

【例2】  (【例1】と同じ内容)
>  library(rgdal)
>  sp.xy <- SpatialPoints(matrix(c(11116.29,-68277.84),nrow=1,byrow=TRUE),proj4string = CRS("+init=epsg:2451"))
>  spTransform(sp.xy, CRS=CRS("+proj=longlat +datum=WGS84"))
SpatialPoints:
     coords.x1 coords.x2
[1,]  139.9557   35.3845
Coordinate Reference System (CRS) arguments:
+proj=longlat +datum=WGS84 +ellps=WGS84
+towgs84=0,0,0

>  CRSargs(CRS("+init=epsg:2451"))
[1] "+init=epsg:2451 +proj=tmerc +lat_0=36 +lon_0=139.8333333333333 +k=0.9999 +x_0=0 +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs"

・参考サイト
平面直角座標 ⇔ 緯度・経度変換 with Google Map
http://www.n-survey.com/online/gmap.htm
Henkan

■GMTによる方法
●緯度経度をUTM座標に変換
・cygwin環境
$ echo 139.955692 35.384502 0 | mapproject  -Ju54R/1.0 -F -C
405151.34090197 3916184.79375918        0

・コマンドプロンプト(MSDOS) 環境
C:\GMTwork>echo 139.955692 35.384502 0 | mapproject  -Ju54R/1.0 -F -C
405151.34090197 3916184.79375918        0

書式指定(gmtset D_FORMAT %.4f)
C:\GMTwork>gmtset D_FORMAT %.4f
C:\GMTwork>echo 139.955692 35.384502 0 | mapproject  -Ju54R/1.0 -F -C
405151.3409     3916184.7938    0

●UTM座標を緯度経度に変換(cygwin環境)
$ gmtset D_FORMAT %.6f
$ echo 405151.34090197 3916184.79375918 0 | mapproject  -Ju54R/1.0 -F -C -I
139.955692      35.384502       0

●参考サイト
MAPPROJECT
http://gmt.soest.hawaii.edu/gmt/html/man/mapproject.html
mapproject - 2 次元座標のフォワード及びインバースの図法変換 
http://gmt.shin-gen.jp/GMT4.1.4/mapproject.html
GMTDEFAULTS
http://gmt.soest.hawaii.edu/gmt/html/man/gmtdefaults.html
GMTDEFAULTS
http://gmt.shin-gen.jp/GMT4.1.4/gmtdefaults.html
測地系 (datum)
http://ryuiki.agbi.tsukuba.ac.jp/~nishida/lecture/GIS/datum.html
地図投影 (map projection)
http://ryuiki.agbi.tsukuba.ac.jp/~nishida/lecture/GIS/proj.html
地理座標系  EPSGコードに関する情報
http://www.geopacific.org/opensourcegis/573074065ea76a197cfb
PROJ.4 - Cartographic Projections Library
http://trac.osgeo.org/proj/
ここで,proj446_win32_bin.zipが入手可能のようです.
README.TXTに説明があります.

C:\>proj.exe  -I +proj=utm +zone=54 +ellps=WGS84
405151.34 3916184.79                   (input)
139d57'20.491"E 35d23'4.207"N     (output)
<cntl-Z>
C:\>
C:\>proj.exe  +proj=utm +zone=54 +ellps=WGS84
139d57'20.491"E 35d23'4.207"N        (input)
405151.34       3916184.79                (output)
139d57'20.491"  35d23'4.207"           (input)
405151.34       3916184.79                (output)

<cntl-Z>
C:\>

■度分秒の経度・緯度データファイルを読み,度に変換してファイル出力するコード例.
関数Scale60to10を使用する

dd=read.csv("keido_ido_data.txt",sep=",",header=FALSE)
allmrow=c()
for (i in 1:nrow(dd)){
  vrow=c()
  for (j in 1:2){
     tmp= dd[i,j]
     henkan=Scale60to10(as.character(tmp))
     vrow=c(vrow,henkan)
  }
  mrow=matrix(vrow,nrow=1)
  allmrow=rbind(allmrow,mrow)
}
write.table(sprintf('%.6f,%.6f',allmrow[,1],allmrow[,2]),"out_test.csv",row.names=FALSE,col.names=FALSE,quote = FALSE)
#write.table(formatC(allmrow,format="f",digits=6),"out_test.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)

●使用例
・入力ファイル(keido_ido_data.txt)
140/10/23,35/33/16
140/5/54,35/41/38

出力ファイル(out_test.csv)
140.173056,35.554444
140.098333,35.693889

・次の二つの文は同じ出力結果となる
write.table(sprintf('%.6f,%.6f',allmrow[,1],allmrow[,2]),"out_test.csv",row.names=FALSE,col.names=FALSE,quote = FALSE)
write.table(formatC(allmrow,format="f",digits=6),"out_test.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)
digits=6は小数点以下の桁数

※write.tableには次の書式もあるようですが,このdigitsは出力桁数を意味.
write.table(format(allmrow,digits=8),"out_test.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)

●write.tableを以下に書き換えると,
write.table(allmrow,"out_test.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)
・出力は次のようになる.
140.173055555556,35.5544444444444
140.098333333333,35.6938888888889

■度の経度・緯度データファイルを読み,度分秒に変換してファイル出力するコード例.
Scale10to60を使う(一部変更【1行追加】)

#10進数を60進数に変換
Scale10to60<- function(decimal10){
  sign<- sign(decimal10)
  decimal10<- abs(decimal10)
  sc1<- trunc(decimal10)
  tmp1<- decimal10-sc1
  sc2<- trunc(tmp1*60)
  tmp2<- tmp1*60-sc2
  sc3<- tmp2*60
  sc1<- sign*sc1
  #Tmp<-paste(sc1,"°",sc2,"′",sc3,"″",sep="")
  sc3=round(sc3,digits=7)#秒の小数点以下の桁数設定のため追加
  Tmp<-paste(sc1,"/",sc2,"/",sc3,sep="")
  return(noquote(Tmp))#noquoteで返す
}

dd=read.csv("indata.txt",sep=",",header=FALSE)
allmrow=c()
for (i in 1:nrow(dd)){
  vrow=c()
  for (j in 1:2){
     tmp= dd[i,j]
     henkan=Scale10to60(tmp)
     vrow=c(vrow,henkan)
  }
  mrow=matrix(vrow,nrow=1)
  allmrow=rbind(allmrow,mrow)
}

write.table(allmrow,"outdata.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)

●indata.txt
140.173055555556,35.5544444444444
140.098333333333,35.6938888888889

●outdata.csv
140/10/23,35/33/16
140/5/54,35/41/38

●秒の小数点以下の桁数設定のため追加行をコメントアウトしたときの出力
140/10/23.0000000016162,35/33/15.9999999998388
140/5/53.9999999987572,35/41/38.0000000000388

■緯度経度(度)をファイルから読み,UTM座標に変換してファイルに出力するコード例
library(rgdal)
dd=read.csv("chiba_longlat.csv",sep=",",header=FALSE)
allout=c()
for (i in 1:nrow(dd)){
  tmp=project(as.matrix(dd[i,]), "+proj=utm +zone=54 ellps=WGS84")
  allout=rbind(allout,tmp)
}
write.table(sprintf('%.1f,%.1f',allout[,1],allout[,2]),"chiba_utm.csv",row.names=FALSE,col.names=FALSE,quote = FALSE)
#write.table(formatC(allout,format="f",digits=1),"chiba_utm.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE) 
#write.table(allout,"chiba_utm.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)

●chiba_longlat.csv 【入力ファイル】
140.173056,35.554444
140.098333,35.693889
140.067500,35.649444
140.117500,35.636944
140.143056,35.664722
140.097778,35.653889

●chiba_utm.csv      【出力ファイル】
425051.7,3934846.1
418420.9,3950371.7
415584.4,3945468.2
420098.4,3944040.0
422439.1,3947100.5
418330.0,3945935.6

・次の二つの文は同じ出力結果となる
write.table(sprintf('%.1f,%.1f',allout[,1],allout[,2]),"chiba_utm.csv",row.names=FALSE,col.names=FALSE,quote = FALSE)
write.table(formatC(allout,format="f",digits=1),"chiba_utm.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE) 

●write.table(allout,"chiba_utm.csv",sep=",",row.names=FALSE,col.names=FALSE,quote = FALSE)のときの出力
425051.676178514,3934846.12675814
418420.918788606,3950371.67801163
415584.410123964,3945468.19202016
420098.389487822,3944040.00033986
422439.147528616,3947100.45228504
418329.959188607,3945935.61832967

■参考サイト
Formatting Using C-style Formats
http://stuff.mit.edu/afs/sipb/project/r-project/arch/i386_rhel3/lib/R/library/base/html/formatc.html
Use C-style String Formatting Commands
http://stuff.mit.edu/afs/sipb/project/r-project/arch/i386_rhel3/lib/R/library/base/html/sprintf.html
Encode in a Common Format
http://stuff.mit.edu/afs/sipb/project/r-project/arch/i386_rhel3/lib/R/library/base/html/format.html


2012年12月29日 (土)

平面直角座標系

■平面直角座標系
わかりやすい平面直角座標系
http://vldb.gsi.go.jp/sokuchi/patchjgd/download/Help/jpc/jpc.htm

平面直角座標系(平成十四年国土交通省告示第九号)
http://www.gsi.go.jp/LAW/heimencho.html



無料ブログはココログ
2017年2月
      1 2 3 4
5 6 7 8 9 10 11
12 13 14 15 16 17 18
19 20 21 22 23 24 25
26 27 28