gdal根据 中心点生成 缓冲区生成shape

引用ESRI.ArcGIS.Client.DLL和GDAL.CSharp.DLL
1、准备中心点坐标
DataTable tbdata = dbHelper.ExecuteDatatable(“tb”, sqlWhere, true);
if (tbdata != null && tbdata.Rows.Count > 0)
{
foreach (DataRow dr in tbdata.Rows)
{
string XZB = dr[“xzb”].ToString();
string YZB = dr[“YZB”].ToString();
string SHAPE = string.Empty;
//如果中心点坐标是球面坐标系则转成2000平面大地坐标
if (!string.IsNullOrEmpty(XZB) && Convert.ToDouble(XZB) < 360 && !string.IsNullOrEmpty(YZB) && Convert.ToDouble(YZB) < 360)
{
ESRI.ArcGIS.Client.Geometry.MapPoint mp = new ESRI.ArcGIS.Client.Geometry.MapPoint
{
X = Convert.ToDouble(XZB),
Y = Convert.ToDouble(YZB)
};
ESRI.ArcGIS.Client.Geometry.MapPoint mpNew = CoordinateTrans.Wgs84ToXian80(mp);

//根据平面坐标系中心点外扩10米生成缓冲区
string str = @"{‘type’: ‘Point’,‘coordinates’: [" + mpNew.X + “, " + mpNew.Y + “]}”;
Geometry pt = Ogr.CreateGeometryFromJson(str);
Geometry circlePolygon = pt.Buffer(10, 30);//10 指缓冲区（圆）的距离（半径），30指点之间的密集度默认30这个没有具体单位，是这个gdal自定义的一个数值范围
string outJson = circlePolygon.ExportToJson(null);
SHAPE = outJson.Replace(”“type”: “Polygon”,", “”).Replace(“coordinates”, “rings”);
SHAPE = SHAPE.Replace(" “, “”);
IDbDataParameter[] idbDataParameter = new IDbDataParameter[1];
idbDataParameter[0] = new OracleParameter(”:SHAPE", OracleDbType.Clob);
idbDataParameter[0].Direction = ParameterDirection.InputOutput;
idbDataParameter[0].Value = SHAPE;
updataSql = string.Format(“update tb t set t.xzb = ‘{0}’, t.yzb= ‘{1}’,t.shape = :SHAPE where t.xzqdm ={3} and t.wybh = ‘{4}’”, mpNew.X, mpNew.Y, SHAPE, dr[“xzqdm”], dr[“wybh”]);
dbHelper.ExecuteNonQuery(updataSql, CommandType.Text, idbDataParameter);
}
}

public class EarthParam
{
///
/// 长轴 eg:6378137
///
public double A { get; set; }
///
/// 扁心率 eg:1 / 298.257223563
///
public double F { get; set; }

    public static EarthParam WGS84 { get; private set; }
    public static EarthParam XIAN80 { get; private set; }
    public static EarthParam BJ54 { get; private set; }
    static EarthParam()
    {
        WGS84 = new EarthParam() { A = 6378137, F = 1 / 298.257223563 };
        XIAN80 = new EarthParam() { A = 6378140, F = 1.0 / 298.3 };
        BJ54 = new EarthParam() { A = 6378245.0, F = 1 / 298.257 };
    }
}

///
/// 84经纬度坐标转西安80平面坐标
///
///
///
///
public static MapPoint Wgs84ToXian80(MapPoint sourcePoint, Param4 param = null, int ProjNo = 0)
{
if (sourcePoint == null)
{
throw new ArgumentNullException(“sourcePoint”);
}
MapPoint targetPoint = GaussProjCal(sourcePoint, EarthParam.WGS84, ProjNo);
if (param != null && param.PX != 0 && param.SX != 0)
{
return Tras4Param(targetPoint, param);
}
else
{
return targetPoint;
}
}
///
/// 高斯投影正算
///
///
///
///
public static MapPoint GaussProjCal(MapPoint sourcePoint, EarthParam param, int ProjNo = 0)
{
double x, y;
GaussProjCal(sourcePoint.X, sourcePoint.Y, param.A, param.F, out x, out y, ProjNo);
return new MapPoint(x, y);
}
//高斯投影正、反算
//3度带宽
//高斯投影由经纬度(Unit:DD)反算大地坐标(含带号，Unit:Metres)
private static void GaussProjCal(double longitude, double latitude, double a, double f, out double X, out double Y, int ProjNo = 0)
{
//int ProjNo = 0;
int ZoneWide; 带宽
double longitude1, latitude1, longitude0, latitude0, X0, Y0, xval, yval;
double e2, ee, NN, T, C, A, M, iPI;
iPI = Math.PI / 180; 3.1415926535898/180.0;
ZoneWide = 3; 3度带宽
//a = 6378245.0; f = 1.0 / 298.3; //54年北京坐标系参数
//a=6378140.0; f=1/298.257; //80年西安坐标系参数
//a=6378137m;f=1/298.257223563;//WGS-84坐标系
if (ProjNo == 0)
{
ProjNo = (int)Math.Round(longitude / ZoneWide);
}
longitude0 = ProjNo * ZoneWide;
longitude0 = longitude0 * iPI;
latitude0 = 0;
longitude1 = longitude * iPI; //经度转换为弧度
latitude1 = latitude * iPI; //纬度转换为弧度
e2 = 2 * f - f * f;
ee = e2 * (1.0 - e2);
NN = a / Math.Sqrt(1.0 - e2 * Math.Sin(latitude1) * Math.Sin(latitude1));
T = Math.Tan(latitude1) * Math.Tan(latitude1);
C = ee * Math.Cos(latitude1) * Math.Cos(latitude1);
A = (longitude1 - longitude0) * Math.Cos(latitude1);
M = a * ((1 - e2 / 4 - 3 * e2 * e2 / 64 - 5 * e2 * e2 * e2 / 256) * latitude1 - (3 * e2 / 8 + 3 * e2 * e2 / 32 + 45 * e2 * e2 * e2 / 1024) * Math.Sin(2 * latitude1)
+ (15 * e2 * e2 / 256 + 45 * e2 * e2 * e2 / 1024) * Math.Sin(4 * latitude1) - (35 * e2 * e2 * e2 / 3072) * Math.Sin(6 * latitude1));
xval = NN * (A + (1 - T + C) * A * A * A / 6 + (5 - 18 * T + T * T + 72 * C - 58 * ee) * A * A * A * A * A / 120);
yval = M + NN * Math.Tan(latitude1) * (A * A / 2 + (5 - T + 9 * C + 4 * C * C) * A * A * A * A / 24
+ (61 - 58 * T + T * T + 600 * C - 330 * ee) * A * A * A * A * A * A / 720);
X0 = 1000000L * (ProjNo) + 500000L;
Y0 = 0;
xval = xval + X0; yval = yval + Y0;
X = xval;
Y = yval;
}

///
/// 平面四参数转换方法
///
///
///
///
public static MapPoint Tras4Param(MapPoint sourcePoint, Param4 param4)
{
return Tras4Param(new MapPoint[] { sourcePoint }, param4)[0];
}

    public static MapPoint[] Tras4Param(MapPoint[] sourcePoints, Param4 param4)
    {
        if (sourcePoints == null || sourcePoints.Length == 0)
            throw new ArgumentNullException("sourcePoints");
        if (param4 == null)
            throw new ArgumentNullException("param4");
        double[][] parameter = new double[4][];
        parameter[0] = new double[] { param4.PX };
        parameter[1] = new double[] { param4.PY };
        parameter[2] = new double[] { param4.SX };
        parameter[3] = new double[] { param4.SY };
        int m = sourcePoints.Length;
        double[][] oldp = new double[m][];
        for (int i = 0; i < m; i++)
        {
            oldp[i] = new double[] { sourcePoints[i].X, sourcePoints[i].Y };
        }
        double[][] newp = Tras4Param(oldp, parameter, m);
        MapPoint[] targetPoints = new MapPoint[m];
        for (int i = 0; i < m; i++)
        {
            targetPoints[i] = new MapPoint(newp[i][0], newp[i][1]);
        }
        return targetPoints;
    }

    /// <summary>
    /// 四参数转换
    /// </summary>
    /// <param name="oldpoint"></param>
    /// <param name="parameter"></param>
    /// <param name="M"></param>
    /// <returns></returns>
    private static double[][] Tras4Param(double[][] oldpoint, double[][] parameter, int M)
    {
        int i;
        double[][] newpoint;                    //用于存储旧坐标经转换后的新坐标   
        newpoint = new double[M][];
        for (i = 0; i < M; i++)
        {
            newpoint[i] = new double[2];
        }
        for (i = 0; i < M; i++)
        {
            newpoint[i][0] = parameter[0][0] + oldpoint[i][0] * parameter[2][0] - oldpoint[i][1] * parameter[3][0];
            newpoint[i][1] = parameter[1][0] + oldpoint[i][1] * parameter[2][0] + oldpoint[i][0] * parameter[3][0];
        }
        return newpoint;
    }