[GDAL/OGR In CSharp] C#のOGRを使って作図する(Point)

GDAL/OGRのC#(.net)用ラッパーを使ってシェープファイルを作図してみました。
C#での解説では見つけられなかったのですが、Pythonでの解説1)http://invisibleroads.com/tutorials/gdal-shapefile-points-save.htmlを見ると、下記のような感じでデータを持っているようです。

  • Point,Lineなどのジオメトリ(Geometry)は地物(Feature)に所属している。
  • 地物はレイヤ(Layer)に所属している
  • レイヤはデータソース(DataSource)に所属している
  • データソースはドライバ(Driver)のデータフォーマットによって保存されている。

 Driver
  Datasource
   Layer
    Feature
     Geometry(Point、Line…)

ということで、試行錯誤の末にコードを書いてみました。
単に図形を登録するだけでなく、付属のDBFファイルに保存される属性値も文字列と実数(浮動小数点)のカラムを作って保存してあります。

using OGR = OSGeo.OGR;
using OSR = OSGeo.OSR;

/*中略*/
const string sGdalDir = @"C:\Program Files\GDAL"; //GDALのインストールディレクトリ
const string sShp = @"D:\hogehoge\1_pt.shp";
const string sCsv = @"D:\hogehoge\1_pt.csv";

//投影法:OGC WKTフォーマット
const string sPrj = "PROJCS[\"JGD2000 / UTM zone 54N\",GEOGCS[\"JGD2000\",DATUM[\"Japanese_Geodetic_Datum_2000\",SPHEROID[\"GRS 1980\",6378137,298.257222101,AUTHORITY[\"EPSG\",\"7019\"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY[\"EPSG\",\"6612\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AUTHORITY[\"EPSG\",\"4612\"]],PROJECTION[\"Transverse_Mercator\"],PARAMETER[\"latitude_of_origin\",0],PARAMETER[\"central_meridian\",141],PARAMETER[\"scale_factor\",0.9996],PARAMETER[\"false_easting\",500000],PARAMETER[\"false_northing\",0],UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],AXIS[\"Easting\",EAST],AXIS[\"Northing\",NORTH],AUTHORITY[\"EPSG\",\"3100\"]]";

static void Main(string[] args)
{
    #region "OGRの準備(環境変数Pathの設定→Ogr.RegisterAll())"
    string sVal = Environment.GetEnvironmentVariable("PATH") + ";" +
    sGdalDir + ";" +
    Path.Combine(sGdalDir, "csharp") + ";" +
    Path.Combine(sGdalDir, "gdal-data") + ";" +
    Path.Combine(sGdalDir, "gdalplugins") + ";" +
    Path.Combine(sGdalDir, "projlib");
    Environment.SetEnvironmentVariable("PATH", sVal);
    OGR.Ogr.RegisterAll();
    #endregion

    //Shpを作成する
    //Driver → DataSource → Layer
    OGR.Driver drv = OGR.Ogr.GetDriverByName("ESRI Shapefile");
    using (OGR.DataSource ds = drv.CreateDataSource(sShp, null))
    using (OGR.Layer layer = ds.CreateLayer("layer", new OSR.SpatialReference(sPrj), OGR.wkbGeometryType.wkbPoint, null))
    {
        OGR.FeatureDefn fetDef = layer.GetLayerDefn(); //レイヤと地物を結びつける?

        //属性テーブルにフィールド追加
        OGR.FieldDefn fdStr = new OGR.FieldDefn("fStr", OGR.FieldType.OFTString); //文字列型
        layer.CreateField(fdStr, 1);

        OGR.FieldDefn fldDef = new OGR.FieldDefn("fX", OGR.FieldType.OFTReal); //実数型
        layer.CreateField(fldDef, 1);

        for (long l = 0; l < 100; l++) //地物を100個作成(Pointは1地物に1ジオメトリ:Pointのみ)
        {
            //Pointを作成
            double dX = l % 10; double dY = l / 10; //座標
            OGR.Geometry pt = new OGR.Geometry(OGR.wkbGeometryType.wkbPoint);
            pt.SetPoint(0, dX, dY, 0d);

            //地物を作成→上記のPointを登録&属性値の設定→レイヤに登録
            using (OGR.Feature f = new OGR.Feature(fetDef))
            {
                //属性テーブルのカラムに値を割り当て
                f.SetField("fStr", string.Format("ID:{0}", l));
                f.SetField("fX", dX);
                
                //地物のジオメトリを設定
                f.SetGeometry(pt);
                f.SetFID(l);
                
                //レイヤに地物を登録
                layer.CreateFeature(f);
            }
        }
    }
}

オマケ:OGR.Ogr.RegisterAll()で「System.TypeInitializationException」が出たときは、構成マネージャーでx68かx64に変更(64bit版GDALを入れていればx64、32bit版ならx86)すると幸せになれるでしょう。

つづき:[GDAL/OGR In CSharp] C#のOGRを使って作図する(Polygon)

脚注

脚注
1 http://invisibleroads.com/tutorials/gdal-shapefile-points-save.html