はじめに

かつてGIS(地理情報システム)の補助ツールを開発したとき、指定した秒単位で直線(線分)を分割してプロットする必要が生じたのですが、
探し方が悪かったのか、そのロジックを検索しても見つからないことがありました。
仕方がないので文系アラフォーながら渋々と初歩の数学から勉強し直して、そして解法に辿り着いたので、
それをここに紹介いたします。

考え方

順を追って説明していきますが、要するにまずベクトルの傾きを求めて、次に大きさの比率を求めたものを乗ずると言うだけです。

プロット間隔(秒)をミリ秒に換算する

通常、地図系に於ける秒はミリ秒で演算されるため。

線分に注目する

下記の枠内の線分に注目すると一つのベクトルデータ。すなわち線分となる。

つまり図例のような折れ線は線分の集まりであるのでループ処理で線分を取り出していくこととなる。

演算の単純化

プログラムにおいてロジックの単純化がもたらすメリットは計り知れない。
ということで演算を単純化する手法を考えた。

2013/07/12 追記：

単純化しない場合、余弦定理を利用した解決が後段で必要になります。

その方法とは全てのベクトルのベクトル方向を無視し、始点(0,0)からの単一ベクトルとみなすことである。
具体的にはベクトル方向をフラグ化し、ベクトル量を絶対値にする。

ベクトル方向の単純化

ベクトルの向きは緯度と経度、それぞれ以下のように2つのフラグで管理できる。

ベクトル量の絶対値化

ベクトル量については線分の始点と終点の差分を絶対値で取得すればよい。

なお、ベクトルの大きさにはそもそも符号が付かないので微妙におかしい表現だが、気にしない方向で。

直角三角形を求める

三角比を使うため直角三角形を導出する。
ここまでの処理で始点を0とする1次関数座標が形成されたので、三角形を投影する。

三角形の斜辺を求める

投影した三角形について三平方定理により斜辺(図例では辺b)の長さを求める。

傾きを求める

投影した三角形について三辺の長さが求まったので∠CABの角度(Cosθ)を求める。

Cosθ=c/b

傾きをラジアン値にする

傾きは求められたが、前項で求まったのは度数法である。
しかしながら通常の言語では度数法で演算が出来ないため、Cosθを逆余弦関数(.NETではMath.Acos)を使いラジアン値に変換する。

こうして求まった角度がベクトルの傾きとなる。

プロット間隔の比率を求める

まず始点(0,0)から指定されたプロット間隔rを半径とした円と、求まった傾きを持つ直線を仮想的に考え、その交点Pを求める方法を考える。

Step.1

まず点Pから垂線と水平線を考え、それぞれのX軸とY軸との交点について考える。

Step.2

従前に求めた角度と三角関数(.NETではMath.CosとMath.Sin)を用いることで座標の比率が求められる。

求めた比率をオフセット値にする

求めた比率に指定間隔を乗じるとオフセット値となる。

ベクトル方向ごとに加算、ないし減算を行う

2つのフラグを判定しオフセット値を加算(減算)していく。

ベクトル量が終点を越えたらそこで処理完了となる。

課題

演算処理で求めた角度(=ベクトルの傾き)については演算誤差が発生するため、長い描画オブジェクト(km単位)だと終点での誤差が目立ちます。
これを修正するためには、ベクトルを内分して指定間隔で傾きを求め直す補正を行うことで誤差を縮小できると思います。

参考コード(C#)

平易な表現でコーディングしてますので解読に問題はないと思います。

private List CreateDividedPolyLine(List TargetLine, int iPitch)
{
    // 引数 TargetLine:処理したい折れ線のリスト
    // 引数 iPitch    :与える間隔(ミリ秒)

    Point newPts = new Point();                  //作成した頂点座標
    List newPtsList = new List();  //作成した頂点座標リスト


    //【対象オブジェクトを頂点単位でループして順次頂点データを生成する】
    Drawing.Point[] pts = TargetLine.ToArray();
    int ptsMax = pts.Length - 1;
    //対象オブジェクトの頂点-1の数だけループ
    for (int i = 0; i <= ptsMax - 1; i++) {
        System.Drawing.Point crntP = pts[i];        //現在の頂点
        System.Drawing.Point nextP = pts[i + 1];    //次の頂点
        bool vrtFlag = false;                        //垂直線フラグ
        bool hrzFlag = false;                        //水平線フラグ

        //(1).最初の頂点は無条件に格納する
        newPts.X = crntP.X;
        newPts.Y = crntP.Y;
        newPtsList.Add(newPts);

        //(2).水平線か垂直線か判定する
        if (crntP.Y == nextP.Y) {
            hrzFlag = true;
        } else if (crntP.X == nextP.X) {
            vrtFlag = true;
        }

        //(3).オフセット値を求める
        int xOffset = 0;
        int yOffset = 0;
        if (hrzFlag == true) {
            // 水平線はXのみ操作
            xOffset = iPitch;
        } else if (vrtFlag == true) {
            // 垂直線はYのみ操作
            yOffset = iPitch;
        } else {
            // 斜線は座標から仮想的に三角形を投影してポイント生成処理をする
            // 
            //     　C         ∠CABに対して
            // 　b ／|         Sin=a/b
            // 　／　| a       Cos=c/b
            // ／　　|         Tan=a/c
            //A￣￣￣ B
            //   c
            //
            //(処理1)高さaと底辺cを求める
            int a = Math.Abs(crntP.Y - nextP.Y);
            int c = Math.Abs(crntP.X - nextP.X);
            //(処理2)三平方定理で斜辺を求める
            double b = Math.Sqrt((Math.Pow(a, 2)) + (Math.Pow(c, 2)));
            //(処理3)三辺が求まったので∠CABの角度(=Cosθ=c/b)を求める
            double theta = c / b;
            //(処理4)求まった角度Cosθをラジアン値に変換する
            double rad = Math.Acos(theta);
            //(処理5) 指定されたプロット間隔を半径とした仮想的な円を考え、それと仮想三角形との
            //    交点からXY値(の比率)を取得してプロット間隔と乗ずる
            xOffset = Convert.ToInt32(Math.Round(iPitch * Math.Cos(rad), 0));
            yOffset = Convert.ToInt32(Math.Round(iPitch * Math.Sin(rad), 0));
        }

        //(4).ベクトル方向を判定
        bool lonFlag = false;        //経度の増加フラグ(正方向=True , 負方向=False)
        bool latFlag = false;        //緯度の増加フラグ(正方向=True , 負方向=False)
        if (crntP.X <= nextP.X) {
            lonFlag = true;
        }
        if (crntP.Y <= nextP.Y) {
            latFlag = true;
        }

        //(5).ベクトル方向に応じたオフセットを付与
        int x = crntP.X;
        int y = crntP.Y;
        do {
            if (lonFlag) {
                x += xOffset;
            } else {
                x -= xOffset;
            }
            if (latFlag) {
                y += yOffset;
            } else {
                y -= yOffset;
            }

            // 一つの線分を超えたらループ脱出
            if (lonFlag == true && latFlag == true) {
                if (x >= nextP.X && y >= nextP.Y)
                    break;
            } else if (lonFlag == true && latFlag == false) {
                if (x >= nextP.X && y <= nextP.Y)
                    break;
            } else if (lonFlag == true && latFlag == false) {
                if (x <= nextP.X && y >= nextP.Y)
                    break;
            } else {
                if (x <= nextP.X && y <= nextP.Y)
                    break;
            }

            // 生成したポイントを追加
            newPts = new Point();
            newPts.X = x;
            newPts.Y = y;
            newPtsList.Add(newPts);
        } while (true);
        newPts = new Point();
        newPts.X = nextP.X;
        newPts.Y = nextP.Y;
    }

    //最後の頂点は無条件に格納
    newPtsList.Add(newPts);

    //値を返す
    return newPtsList;

}