2012-01-17 16 views
4

KMLで円を描くGPS座標を中心とした地球上の円上の点を計算する方法は?

地球儀上の点の座標(小数点以下の形式など)を取得し、その点を中心とする円に近似する多角形の座標を生成するにはどうすればよいですか?

20以上のデータポイントを持つポリゴンは円のように見えます。データポイントが多いほど、サークルの見栄えがよくなります。

私はKMLを生成するプログラムを作成していますが、ポリゴン頂点の座標を計算する方法はわかりません。データ入力の

例:

緯度、経度、円の半径(フィート)、NumberOfDataPoints

26.128477

、-80.105149、500、20

+0

あなたのプログラムにはどのような言語を使用していますか? – arx

答えて

6

これがある場合、私は知りません最も単純な解決法であり、世界は球であると仮定します。

を定義:

Rは球(すなわちアース)の半径です。

rは円の半径(同じ単位)です。

tは、t = r/Rラジアンのように、球の中心で長さrの大きな円の円弧で囲まれた角度です。

ここで、球が半径1で原点の中心にあると仮定します。

Cは、円の中心を表す単位ベクトルです。

北極を一周する円を想像し、円の平面が地球の中心から北極までの線と交差する点を考えてみましょう。明らかに、この点は北極のどこかにあるでしょう。

KはC(すなわち、サークルの面がCと交差する場所)、「下」は、対応する点であるので、K = COS(T)C

S(すなわち3次元空間で測定された円の半径はありません球上の)は、S = SIN(T)

が今は必要中心K、半径Sとを有する3次元空間内の円上の点K.

This answer

を通って垂直な平面内にある(無視回転物)は、平面の基底ベクトル(すなわち、法線KまたはCに直交するベクトル)を見つける方法を説明します。クロスプロダクトを使用して秒を探します。

コールこれらの基底ベクトルUとVの

// Pseudo-code to calculate 20 points on the circle 
for (a = 0; a != 360; a += 18) 
{ 
    // A point on the circle and the unit sphere 
    P = K + s * (U * sin(a) + V * cos(a)) 
} 

球座標に各ポイントを変換し、あなたが行われています。

退屈なので、私はこれをC#でコーディングしました。結果は妥当である:彼らは円の中にあり、球の上にある。ほとんどのコードでは、ベクトルを表すstructが実装されています。実際の計算は非常に簡単です。

using System; 

namespace gpsCircle 
{ 
    struct Gps 
    { 
     // In degrees 
     public readonly double Latitude; 
     public readonly double Longtitude; 

     public Gps(double latitude, double longtitude) 
     { 
      Latitude = latitude; 
      Longtitude = longtitude; 
     } 

     public override string ToString() 
     { 
      return string.Format("({0},{1})", Latitude, Longtitude); 
     } 

     public Vector ToUnitVector() 
     { 
      double lat = Latitude/180 * Math.PI; 
      double lng = Longtitude/180 * Math.PI; 

      // Z is North 
      // X points at the Greenwich meridian 
      return new Vector(Math.Cos(lng) * Math.Cos(lat), Math.Sin(lng) * Math.Cos(lat), Math.Sin(lat)); 
     } 
    } 

    struct Vector 
    { 
     public readonly double X; 
     public readonly double Y; 
     public readonly double Z; 

     public Vector(double x, double y, double z) 
     { 
      X = x; 
      Y = y; 
      Z = z; 
     } 

     public double MagnitudeSquared() 
     { 
      return X * X + Y * Y + Z * Z; 
     } 

     public double Magnitude() 
     { 
      return Math.Sqrt(MagnitudeSquared()); 
     } 

     public Vector ToUnit() 
     { 
      double m = Magnitude(); 

      return new Vector(X/m, Y/m, Z/m); 
     } 

     public Gps ToGps() 
     { 
      Vector unit = ToUnit(); 
      // Rounding errors 
      double z = unit.Z; 
      if (z > 1) 
       z = 1; 

      double lat = Math.Asin(z); 

      double lng = Math.Atan2(unit.Y, unit.X); 

      return new Gps(lat * 180/Math.PI, lng * 180/Math.PI); 
     } 

     public static Vector operator*(double m, Vector v) 
     { 
      return new Vector(m * v.X, m * v.Y, m * v.Z); 
     } 

     public static Vector operator-(Vector a, Vector b) 
     { 
      return new Vector(a.X - b.X, a.Y - b.Y, a.Z - b.Z); 
     } 

     public static Vector operator+(Vector a, Vector b) 
     { 
      return new Vector(a.X + b.X, a.Y + b.Y, a.Z + b.Z); 
     } 

     public override string ToString() 
     { 
      return string.Format("({0},{1},{2})", X, Y, Z); 
     } 

     public double Dot(Vector that) 
     { 
      return X * that.X + Y * that.Y + Z * that.Z; 
     } 

     public Vector Cross(Vector that) 
     { 
      return new Vector(Y * that.Z - Z * that.Y, Z * that.X - X * that.Z, X * that.Y - Y * that.X); 
     } 

     // Pick a random orthogonal vector 
     public Vector Orthogonal() 
     { 
      double minNormal = Math.Abs(X); 
      int minIndex = 0; 
      if (Math.Abs(Y) < minNormal) 
      { 
       minNormal = Math.Abs(Y); 
       minIndex = 1; 
      } 
      if (Math.Abs(Z) < minNormal) 
      { 
       minNormal = Math.Abs(Z); 
       minIndex = 2; 
      } 

      Vector B; 
      switch (minIndex) 
      { 
       case 0: 
        B = new Vector(1, 0, 0); 
        break; 
       case 1: 
        B = new Vector(0, 1, 0); 
        break; 
       default: 
        B = new Vector(0, 0, 1); 
        break; 
      } 

      return (B - minNormal * this).ToUnit(); 
     } 
    } 

    class Program 
    { 
     static void Main(string[] args) 
     { 
      // Phnom Penh 
      Gps centre = new Gps(11.55, 104.916667); 

      // In metres 
      double worldRadius = 6371000; 
      // In metres 
      double circleRadius = 1000; 

      // Points representing circle of radius circleRadius round centre. 
      Gps[] points = new Gps[20]; 

      CirclePoints(points, centre, worldRadius, circleRadius); 
     } 

     static void CirclePoints(Gps[] points, Gps centre, double R, double r) 
     { 
      int count = points.Length; 

      Vector C = centre.ToUnitVector(); 
      double t = r/R; 
      Vector K = Math.Cos(t) * C; 
      double s = Math.Sin(t); 

      Vector U = K.Orthogonal(); 
      Vector V = K.Cross(U); 
      // Improve orthogonality 
      U = K.Cross(V); 

      for (int point = 0; point != count; ++point) 
      { 
       double a = 2 * Math.PI * point/count; 
       Vector P = K + s * (Math.Sin(a) * U + Math.Cos(a) * V); 
       points[point] = P.ToGps(); 
      } 
     } 
    } 
} 
0

私は、それをしないPythonで小さなオープンソースパッケージをPolycirclesを書かれています。地理空間計算にはgeographiclibが使用されます。

enter image description here

関連する問題