保津峡谷鵜ノ川口から峡谷右岸基盤岩露頭へのアプローチ Approach to the basement rock outcrops on the right bank at the upper end of the Hozugawa Gorge (Unokawa-guchi)

はじめに

 保津峡谷口,または保津峡谷入り口についてのぼくの言及を,友人は地理学の谷口集落の連想からか,嵐山と誤解した。それを踏まえてこの投稿の英語タイトルを決定した。upper endが重要だ。日本語では地名があるので不要だろう。

 防災や土木用語としては請田サイトが使われている。時刻水位月表検索

観測所記号 観測所名 水系名 河川名
306041286606300 請田(ウケタ) 淀川 桂川

 しかし驚くほど情報がない。どこかに委託しているのか,サボっているのか。

 亀岡の水害を調べ始めると,請田,という地名は使われているが,実質がないようだ。

以上,2026年2月16日。

1. アプローチの選択

 請田神社から府道を使って保津峡谷に入ると,かつての山陰線(現トロッコ列車軌道)の上方には大きく切ったと思われる崖面が見える(図1にはその下部だけ)。この場が大山咋命(おほやまくひのみこと)の丹の海干拓のための主要な工事現場であったとぼくは当初考えていた。

図1 請田神社前から保津峡谷入り口を

 「古代タニハ『丹の海』とその排水プラグを地形発達史から復元 Part 1」関西大学文学論集75-4は本日(2026年3月12日)刊行されたが,この報告では,大山咋命が見た丹の海の湖面の海抜高度が98mであったことを明らかにした。図1の付近は海抜80mだからおよそ落差20mの滝があったとぼくは考えていた。その滝を数十年またはそれ以上の年数をかけて,開削をすることで,湖底は侵食されていったと想像したのである。開削工事については,上記論文のPart 2として発表する予定である。

以上,2026年3月2日。2026年3月12日修正。

 この観点ゆえに,図1に見える右岸の岩質を確認する必要があった。元の山陰線,現トロッコ列車軌道の上方にはコンクリートの大きな壁面が見え岩盤は隠されているので,河原で観察するしかない。

2. 右岸河原に辿り着くには

 冬の渇水で鵜ノ川を渡ることができた。図のルートでまずは鵜ノ川の右岸へ。そして,図2に見える堰を図3のように伝って人工島にたどり着いた。高齢でバランスを崩す可能性が濃厚で,堰の途中で,キャラバンからウェットブーツに履き替えた。

図2 鵜ノ川河口から保津峡を

図3 鵜ノ川右岸に降り立つルート

図4 堰き止めバリアを伝って人工島へ

 図5には,軌道面から一段低いテラスから川への梯子が見えるので,鵜ノ川の流れが早い場合にはテラス面を使うルートもありうると思っている。図7は保津川の本流路でこれが保津川下りのルートを確保している。図1は請田神社からの眺望であるが,この人工島で分割されて,保津川下りの航路が確保されているのがよくわかる。

図5 軌道面から一段低いテラスから川への梯子が

図6 右手の水域は鵜ノ川への水域区分

図7 保津川下り用の本流路

3. 岩石種は玄武岩か

 保津川右岸河床に立つと,左岸の黒色岩石の露頭が見える(図8〜10)。なぜか,左岸では節理が目立つ。

図8 左岸露頭繋ぎ写真 左

図9 左岸露頭繋ぎ写真 中

図10 左岸露頭繋ぎ写真 右

 アイフォーンのパノラマ機能で図11を撮影した。左手の河岸上方の朱色は請田神社の柵,中央手前にカーブして流れているのが保津川本流,右手のほぼ静水域は鵜ノ川河口からの水域。中央にはコンクリートで固めた人工島がみえる。

図11 人工島からのパノラマ写真

 図12は,差別侵食で残った枕状熔岩だと思う。図13,14の正面に見えているマッシブな壁面はコンクリート壁である。最下部にへの字の洞穴が見えるがここには人が積み上げた河床礫がコンクリート壁が侵食された結果,露出している。気になって現場に行ってわかったことである。

図12 根つきの枕状熔岩

図13 図14とのステレオ写真の左

図14 図13とのステレオ写真の右

 図15,16は,コンリート壁のすぐ下流側のステレオ写真である。ベッディングに見える。この産状が下流に続くのを左岸から見ていたので,熔岩ではなく泥岩と考えざるを得ないと考えていた。近接して観察すると,右下がりの「層理面」に間違いなく,節理ではない。ぼくの貧弱な知識だと海洋底の玄武岩とは到底思えなかったのであるが,図12のような部分もあって,海洋底での噴出物の降下堆積物と考えざるを得なくなった。というか,玄武岩であることに安堵したのである。

図15 図14とのステレオ写真の左

図16 図15とのステレオ写真の左

 図17,18は,図16の壁面の一部であるが,チャートも含まれているようである。

図17 チャートも含まれるか

図18 図17のズームイン

 図19〜21では熔岩と火山礫・火山灰相が見える。

図19 図20とのステレオ 左

図20 図19とのステレオ 右

図21 熔岩レンズと火山灰

 図19〜23は同じ場所で,図15からすると,200mほど下流のやはり右岸である。

図22 熔岩と火山礫相

図23 熔岩と火山礫相

 請田神社前(図25)の河床には玄武岩が見えるが,少し上手では泥岩が見えているようだ。この上の府道の露頭では,泥岩にチャートが混入している。

図24 請田神社前のうちの上流側

図25 請田神社前

おわりに

 これから水位が上がると右岸に到達するのは腰まで水に浸かる必要が出てくるかもしれない。請田神社側つまり左岸は道路から簡単に草木の繁る斜面を降りると岩盤に辿り着くことができる。図25の水位観測点の下流側から降りるのが簡単のようだ。(調査日:Feb. 15, 2026)

以上,2026年3月10日。2026年3月12日修正。




Start GIMP instead of Adobe Photoshop

はじめに

 awesome Adobe applications to sharewares で示したように,Adobe Photoshopから離れて,GIMPを使用したい。

 GIMP使用当初はPhotoshopに比べてステップ数が多いと感じたのであるが慣れてきて,かつてGIMP使用を断念した時よりはよほど肩が軽くなった。Photoshopから拒絶されたので致し方ない。

1. ダウンロード

 https://www.gimp.org/downloads/ に,GIMP, Current Stable Version

The current stable release of GIMP is 3.0.2 (2025-03-23), のダウンロードサイトがある。Donationを以前したことがある。使えるようになったら改めてしたいと思う。

 ver. 2.1を削除してこのver. 3.0をインストールしたのであるが,ツールボックス表示のタイミングがうるさくなくなったという印象を持つが,戸惑うことがあり,Perplexityにお世話になった。

2. 最近は何をしていたのか

 Adobe Photoshopで最近何をしてきたのかを考える。今後,必要に応じて追加してゆく。

2.1 スクリーンショットや写真からWeb用コンテンツを作成

 現在の最頻度の作業は,macやWindowsでの作業をスクリーンショットで貯めて,それらを軽くしてWebページへ掲載してきた。フィールドなどで撮影した写真も同様の利用をする。iPhoneで撮影した写真は,mac写真アプリに並べてコメントを付し,jpgで出力してそれらをスクリーンショットとほぼ同様の作業をしてきた。

 スクリーンショットはフィールド写真と比べるとかなり軽量で多少処理過程に違いはある。総じてPhotoshopでの作業は次のようであった。

1  必要な部分を切り取る場合がある。

2 写真の場合はスライダーを使って彩度を高くする。

3 ランドスケープの場合は,解像度300dpi,縦1200dpi,横1600dpi。

  ポートレート の場合は,解像度300dpi,縦1600dpi,横1200dpi。

4 容量を軽くする。Web用だとスクリーンショットの場合最大100 kbyteほど,写真だと500kbytesほどにしてきた。

 さて,Perplexityでの質問文を考えてみる。

 スクリーンショット: GIMPでの作業過程を教えてください。スクリーンショットの一部を切り取って,その形がランドスケープの場合は,解像度300dpi,横1600dpiに揃えたい。そしてファイル容量100kbytesほどのjpgで出力したい。

 写真: GIMPでの作業過程を教えてください。iPhoneで撮影した写真の彩度を高めて,その形がランドスケープの場合は,解像度300dpi,横1600dpiに揃えたい。そしてファイル容量500kbytesほどのjpgで出力したい。

3. スクリーンショットをWebページ掲載画像ファイルjpgに変換する手順

 では,スクリーンショットの場合を聞いてみた。回答は簡潔でそのままでは実行できない。以下に。

作業手順

3.1 スクリーンショットの切り抜き

 1.1 GIMPで,メニューのファイル > 開く/インポート

図1 moto > Dropbox > スクリーンショット

 1.2 設定したパスでのファイルを選択,ダブルタップ

図2 ファイルが開いたところ

 ファイルを直接,作業画面にドラッグ&ドロップするのが簡単ではあるが(一般的にはこれが迅速),図1のようにパス指定すると,その後は,パスを覚えてくれている。図1の右手にはファインダーでのスクリーンショットフォルダ内が見えており,対応している。日を跨ぐと,更新日時欄の表示には時刻が表示されない。名前欄の幅が小さいので,スクリーンショットの長いファイル名が表示されることはない。

 目的のファイルをダブルクリックすると選択できる。ファイルを選択したあとに右下隅の「開く」ボタンをタップしてもよい。
 次に,Keep the image’s color profileパネルが表示される。この下端には「維持」と「Convert」ボタンがあり,「維持」ボタンをタップすることになる。ただ,このパネルは,「次回から確認しない」に☑️を入れて良いだろう。ぼくの能力では,このパネルを生かすことはできない。

 図2の作業枠は,画像サイズがこれを超えると,枠域は自動で拡大される。Photoshopには無い融通性だ。

3.2 矩形選択ツールで必要な範囲を選び,切り抜く

図3 切り抜いたところ

 他のツールを使っていない場合は,自動的に矩形選択ツールが使え,選択範囲でクリックすれば,切り抜くことができる。直前のファイル作業で他のツールを使っている場合は,次のような手続きが必要である。

2.1 メニューのツール > 選択ツール > 矩形選択ツール
 2.2 画像 > 選択範囲で切り抜き

3.3 画像の拡大・縮小

図4 画像 > 画像の拡大・縮小

 メニューの画像 > 画像の拡大縮小

  図4で解像度(水平と垂直は連動)を300dpiに幅(高さも連動)を1600pxに設定。確定後、「OK」ボタンはないので,「拡大・縮小」ボタンをタップ。

図5 画像の拡大・縮小を調整して「拡大・縮小」ボタンをタップ

3.4 エキスポート

 4.1 メニューのファイル > 名前を付けてエキスポート

図6 メニューから「ファイル」→「名前を付けてエクスポート」を選択

 図6に示したように,「名前を付けて保存」ではなく,「名前を付けてエクスポート」である。 

 4.2 保存先をダウンロードフォルダにしてJPEGでエキスポート過程に

図7 保存先をダウンロードフォルダにしてJPEGでエキスポート

 JPEG形式で保存しないと容量などの点からウェブコンテンツにならない。図7の下段にはファイル形式の設定メニューが見えるが,「🔽ファイル形式の選択」欄はもともと閉じている。これを開いたのが図7である。このウィンドウのほぼ中央ぐらいに,JPEGが見えるので,それを選択するのである。図7ではupload用フォルダーが選択されているが,GIMPは記憶してくれないので都度このパス指定が必要になり時間を取られるので,左端のフォルダー表示で,Downloadsは簡単に選択できるのでここにエキスポートしたい。

 4.3 ファイル容量を100kB程度に

図8 品質を100kB程度に

 図7の画面の次に図8の画面が現れる。左上隅には,品質調整のスライダーがある。スライダー上の左方部をタップすると品質が低下する。この図8では,147.4kBが見えるが,後背の画像をみつつ,67.2kBにして,エキスポートした。

 エキスポートした画像を見たが全く問題なかった。

 なお,一つの画像に関わる作業が終われば,

メニューのファイル > ビューを閉じる

を選ぶ。保存するかどうか聞いてくるが,すでにエキスポートしているので,保存しない,を選択する。保存にするとGIMP固有のファイルとなり,引き続き作業が可能になる。

線を描く

 図6には,注意喚起のために,赤線を入れている。この入力法を次に。

メニューのツール > ツールボックス

とするとツールボックスが出る。その鉛筆を選んで色も指定する。

鉛筆の太さの変更は次のようにする。

メニューのウィンドウ > ドッキング可能なダイアログ > ツールオプション

を開いて,サイズ,を調整する。

以上,2025年4月5日。

テキストボックスの選択解除

 なぜか,複数のテキストを作成したあと,その一つのテキストボックスが解除されない。メーンメニューの選択 >「 選択の解除」,は灰色になってしまっている。

 Perplexityに質問してその回答の一つを実行して成功した。

————————————————

  1. フローティング選択範囲の確認と解消 • フローティングレイヤーがある場合は、以下の手順で処理します:

    1. レイヤーダイアログを開きます(「ウィンドウ」→「ドッキング可能なダイアログ」→「レイヤー」)。
    2. フローティングレイヤーが表示されている場合は、右クリックして「新しいレイヤーに変換」を選択します。
    3. これによりフローティングレイヤーが通常のレイヤーに変換され、「選択を解除」が使用可能になります

————————————————

 この回答の通りには実行できない。

2.1の [ レイヤーダイアログを開きます(「ウィンドウ」→「ドッキング可能なダイアログ」→「レイヤー」)。] を実行した。図9の図で三つのテキストのうち,「/wp………」のテキストボックスが選択状態になっている。右上に出現しているツールボックスでも,「/wp………」が選ばれている。それで,最下部の「スクリーン」を選択した。

図9 図で三つのテキストのうち「/wp………」のテキストボックスが選択状態

 図10の図中の「/wp………」のテキストボックスの選択が外れている。右端のツールボックスでは,見えにくいが「スクリーン」が選ばれている。幸なるかな。

図10 図中の「/wp………」のテキストボックスの選択が外れた

以上,2025年4月6日。




続国土地理院の基準点成果の利用 a sequel to using control point survey results made by the Geospatial Information Authority of Japan or GSI

はじめに

 国土地理院の基準点成果の利用 は,2022年秋に作成したものである。これを元に,文学論集Vol. 74, No. 4に投稿すべく記述してみて,2年間のブランクもあって,一つの章を追加する必要性を感じ(2年前に残したことであった),それをここにまとめることにした。制限ページ数に近く,ここで書き込んだもののダイジェストを文学論集Vol. 74, No. 4の第Ⅹ章にしたいと思った。

0. 「Ⅹ LiDARルート測量とライダーROI測量を繋ぐ」(コピペ)

2024年秋現在,筆者の利用の観点からすると,最も優れた3D LiDAR ScannerアプリはScaniverse[i] (v.4.0.2) である。Polycam[ii]と,これまで述べてきたMetascanは,マニュアルを読む限り,パワフルなアプリという印象を受けてきたが,1セッションのLiDAR測量継続時間はいずれも10分間に限られる。Polycam Pro [for iOS v.3.5.18 (225)]の説明書には10分間の限界は書かれていないが,10分を超えると,「前のエリアに戻って再開する」と「最初からやり直す」の二択を迫られる。戻っても再開されない。Scaniverseには,この制限に関連する情報は記されていないが,幸い経験的にこのような限界はないように思う。

さて,ここでは本章タイトルを実行する。許されたページ数は残りわずかなので,これまでのような順を追った説明は省く。これまで,基準点4カ所を通るLiDARルート測量結果を平面直角座標系に載せる手法を示してきたが,ROIのライダー測量結果を平面直角座標系に載せる手法は示していない。しかし,賢明な読者はこれから述べる内容をすでに見抜いておられるのではないかと想像している。

本報告の目的は,ROIのLiDAR測量結果を平面直角座標系に載せる手法を示すことにある。測量ルートの選定は当然,ROIを通過するかその近接地にすることになる。ルート測量の途中でROIのLiDAR測量を実施することも可能かも知れない。ただ,ROIでは丹念にLiDAR測量を実施するので通常,ルート測量とは別にしたいと考えるだろう。

そこでルート測量の際に,ROI付近に一時的な二次的基準点を4カ所設置すれば良いと考えられる。ルート測量して,引き続きまたは翌日にでもROIのLiDAR測量を実施することもできる。なんらかの理由で繰り返しLiDAR測量したい場合にもこの4カ所の二次的基準点を繰り返し使うことができるのである。


[i] Scaniverseは,Niantic(ナイアンティック)https://nianticlabs.com/products?hl=ja による開発中のアプリ。NianticはGoogle社内のスタートアップ企業でGoogle EarthもGoogle mapも開発している。現在Googleから独立したらしい。かつては有料だったが無料化されて,先進的なGaussian Splatting機能が付加されてもなお,無料のままだ。

[ii] https://poly.cam Polycam Proを2024年10月に月契約($26.99)してみた。

1. Oct. 26, 2024実験

 曇りの日で雨が降りそうな天気であった。自分の影が映らなくて良いと感じた。Splattingに期待もしていた。Polycam (Pro)が,LiDAR測量の10分間限定を超えるのではないか,という期待があった。先に結果を言えば,PolycamもScaniverseもSplattingはルート測量には適さないということである。LiDAR測量こそ,ルート測量に適していることが判明した。諦めがついた。図1を参照してほしい。

図A 実験ルート

 後述するが,DOIに当たるのが,打越池南縁の歩道である。DOIを公共の基準点とつなぐ目的でルートを設定している。
 南西部の基準点3-61から出発し,打越池南縁のDOIに階段②から入って③の阪急バス停「青松苑」から出て,⑤の交差点を通過して,④の基準点3-42に到達して交差点⑤に戻り,⑥の基準点3-40に北進する。そして,⑦の基準点3-41,②を通過して,①の基準点3-61に戻っている。

文学論集図51に登録済み。

 階段②からバス停③の間がDOIという設定であるが,そのルートのGoogle map (2024年) 画像を図A下部に示している。画像右下隅に見えるスケールの長さは10mである。③と④の間の300mほどの距離の間に,レーベル L-1〜-4を設置しているのである。

 図B,Cは打越池南縁プロムナードのステレオ写真である。ベンチの向こう先端の左手のレンガとアスファルト境界付近の地面にLabel-1を置いていた。左右両方の金網塀には葛が絡まっている。

図B 打越池南縁プロムナード ステレオ左

図C 打越池南縁プロムナード ステレオ右

Label 1〜4を図D〜Gに示す。Label 1は階段に近く,Label 4はバス停に近い。

図D Label 1

図E Label 2

図F Label 3

図G Label 4

以上,2024年10月24日記

2 Polycam社そのものの社会的信頼度は低い

 Splatting機能をもつPolycamとScaniverseではあるが,ぼくの昨日のルート測量には使えないことがわかった。Scaniverseは無料で,Polycamは有料版でないと撮影はできるがその後の計算処理や出力ができない。2024年10月14日付でアップルから3300円が請求されていた。アカデミックで申請して半額の筈であるがフルに請求されている。1カ月当たりの料金である。別途ぼくの銀行からもDebit Cardから出金しているので(Debit Cardのみ受付と要求された),二重の請求である。銀行とVISAからの二重請求である。Debit Cardは即時送金なので,アップルにクレームをつけるしかない。11月11日に銀行引き落としなので,書類を揃えてアップルケアに連絡することになるだろう。
 今,書類は全部揃えた。明日,アップルケアに電話しよう。ぼくのVisaカードの11月11日の支払い明細(3300円),アップルからの領収証(3300円),銀行の取引明細(2103円)。3300円/$26.99= 122円/$,2103円/$13.5= 156円/$,である。

 図H, Iは,購入の際に見た価格やアカデミックディスカウントの情報。図9のhttps://poly.camから入ってアカデミックの契約をして,Debit Cardのみの受付と要求されたので,ぼくの銀行から送金したことになっている。そうしないと,進めなかったのである。

図H Proは$26.99/月

図I acacemicは $(26.99/2) /月

以上,2024年10月27日記。

2.a アップルケアに連絡

 さきほど,アップルケアに問い合わせた。相談員との話の前に,サブスクリプション,の問題として宣言し,iPhoneだけで読めるメッセージが届く。そこにあるリンクから申請もできるが,相談員とコンタクトを持てるとアナウンスがあり,相談員が電話に出るのを待った。事情を説明し,了解いただいた。二日間の審査があって,ぼくのクレームの是非が審査されるという。
 どう審査されるのかを問いただした。PolycamのウェブサイトでのぼくのDebit Cardでの支払いをアップルは把握しようとしないことが判明。ぼくの手持ちの書類を受け取る形も持っていない。アップル関係者が,PolycamのWebサイトでぼくのように,アカデミックディスカウントのリンクに入ってDebit Cardでの支払いをしてみないとわからないと,相談員に伝えた。相談員はかなり優秀でぼくの意見は理解していたようである。Polycamにぼくがクレームをした方が良いという提案もあった。

2.b Polycamにメール

 poly.camと入力すると,https://poly.cam/library に入ることになる。これがトップページと考えて良いだろう。この最上段のリンクの一つに,contact salesリンクがある。個人情報を入れてコメント欄に入力することになる。そのコメントは次のよう。

Please return ¥3300 = $26.99 per month for Polycam Pro to me.
My apple account: xxxxx@xxxxx.ac.xxxx
————————————————

  1. I bought “Polycam for Education” using my Debit Card you need on Oct. 14, 2024.
    “¥2103=$13.5 on Oct. 14, 2024” is shown in the Description of a trade of my city bank.

I had received your message: Polycam for Education
Thank you for applying for Polycam for Education!
Here’s your promo code for 50% off Polycam Pro: EDU0523AXJ
Here’s how to redeem your promo code:
Create or log in to your account on our website at https://poly.cam

  1. I recognized the following supposedly impermissible thing yesterday.

a. I found a mail from Apple Inc. to moto@kansai-u.ac.jp: a receipt of ¥3300 = $26.99 per month for Polycam _ LiDAR 3D scanner (order no: MML2T4M8MN)

b. I checked the billing details of my credit VISA card and found ¥3300 = $26.99 from Apple Com Bill.

————————————————
I made a phone call to the Apple Care, and advised to call you. Thank you. Oct. 28, 2024.

 発送のちすぐに対応するとのことだが,すでに2時間あまりが経過したが返事はない。現地時間は活動時間帯ではないかもしれない。

 このコメント発信ののち,当方のacademyメールをチェックしたら,アップルからの3300円の領収書が来ていた。びっくりして,アップルケアに電話した。どうも,ロボットに,サブスクリプション,と宣言した時点で,アップルから10月14日の領収証が再送されたらしい。

2.b アップルから返金があった

Polycamからの返事は全くなかったが,アップルからOct. 14, -¥3300の連絡があった。このアップルに登録しているクレジットカードの請求明細は変更されていない。次回の明細にOct. 29 -¥3300となるのであろう。

図I-J 1 アップルから払い戻し通知書届く

図I-J 2 返金リクエストの結果メール

図I-J 3 Polycamの欄に返金済みとある

以上,2024年10月30日。

2.c 一ヶ月後にまたデビットカードから泥棒された

 Polycamのサイトにログインしても,脱会のリンクがない。返事も来ない。それで,ユーザー登録そのものを削除した。その際,アクセスできない趣旨の警告があった。iPhone 12 Proからアプリを削除し,iPhone 12 Proのシステム設定のアプリのリストも削除した。前述の1回だけ,Polycamにアップロードしただけであった。

 アクセスしようにも,次の図I-J 4, 5のように既存ユーザーがアクセスするリンクがない。コンタクトしても返事はないのでもう,やる気がしないし。ハイリスクの印象だ。

図I-J 4

図I-J 5

 にもかかわらず,Nov. 7, 2024にキャンセルしないとまた請求するよというメールが届いた。宣伝メールやこの種のメールが届くが怖くて無視。返事はない。そして,Nov. 14に

図I-J 6

図I-J 7

 図I-J 6は,Polycamからのメールだ。その一分後にりそなから,引き落としの連絡だ。

 このPolycamから離脱すべく,りそなのデビットカードに連絡した。デビットカードを停止にしても勝手に引き出される可能性があるという。結局16桁のVISA番号を廃棄するしか方法がない。銀行に出向いて,デビットカードから退会して,キャッシュカードだけにするしかない。銀行に出向いての手続きの予約した。デビットカードは引き落とされたら戻らない。保証がない。かなり危険なツールである。クレジットカードはそういう意味で安全が担保されている。

 ただ,このウェブサイトのコンテンツ作成にともなって銀行口座からのクレジットカード引き落としを最近調べた。少なくとも数年にわたって二ヶ月に1回6000円あまりが抜かれていて,覚えがなく,二ヶ月前かにこの銀行に行って,この引き落としを差し止めた。この会社からはなんの連絡もない。トヨタのTSキュービックカードのようだ。トヨタ販売会社に問い合わせ,車利用のなんらかの諸費用の発生の所在を聞いたが,思い当たらないという。そういう意味でトヨタも危険だ。トヨタは気づいていないのかも知れない。誰かが騙っている。銀行の調べでは引き出している事務所は移転していて,その移転先はわからないという。それでも銀行は気にしない。

 このりそな銀行のデビットカードの相談の後,アップルケアに電話して,アップルストアの問題点を詳細に伝えた。まあ上に上げると言っていたが,二日前の日経ニュースで,アップルからの不当な手数料回避のためにソフト会社が抜け道を模索しているらしい。Polycamそのものが特に悪というわけでなく,システム作りが甘いだけでなないかと思う。返事もないし,Webページに書かれている環境を実現できていないのだろうとは思う。とにかく,POLYCAMはやめた方がいいのは確かだ。Apple Storeはこの面では保証外だ。手数料を下げて,ソフト販売会社の手綱をしっかりと握ってほしいものである。

以上,挿入,2024年11月15日。

3 Polycam – LiDAR 3D スキャナーによるスプラッティング

 スプラッティング splattingに対応する日本語はない。この技術が我々消費者に提供されたのは昨年秋である。NeRFという技術に基づいている。たとえば,NeRFとは?従来の技術との違いとユースケース、これからの課題について 2023/10/28によれば,「NeRF(Neural Radiance Fields) は、カルフォルニア大学の研究員が2020年3月に発表された新しい技術で、ニューラルネットワークを使用して、複数の写真から3Dシーンを生成する技術」とされている。「NVIDIAは,2022年1月に高速化を実現した「Instant NeRF」を発表し、数時間かかっていた処理が数秒で実現させることに成功しました」とある。

 Polycam Proで図Aのルートを21分間でsplatting スキャンした。動画を選んだ。スキャン完了後,そのファイルは1.97GBとある。「uploadと処理」か「後でアップロードする」という選択が出るので,高速のWi-Fi環境がある自宅で実行すべく,さしあたり,「後でアップロードする」を選択した。
 帰宅後,「uploadと処理」を選択すると,図Jが現れる。写真枚数は最低20〜最大2000とあり,今回のは717枚を自動撮影したらしい。動画を選んだのであるが,写真も撮られる。実際,スキャン中,シャッター音がひっきりなしにあった。基準点や交差点での信号待ちの際にはパシャパシャとシャッター音がしていた。図JのProcessing modeでまずは,Gaussian splatボタンを選択したら,図Kが現れた。図Jの詳細メニューは無くなっている。「アップロードと処理(1.97GB)」のボタンをタップした。
 図Lはprocessingの途中である。18:24にスタートして,20:03でアップロード中25%となっている。図Mはprocessingの終了直後である。4時間20分ほどprocessingにかかったということである。サーバーへのアップロードが終わったということであり,次にサーバーでsplattingの解析が始まっているのである。図Nはライブラリーに戻った際の表示で,対象ファイルは,「処理中」とある。しかし,図Oのように,23:06には「キャプチャー処理に失敗しました」というメッセージが現れた。

 つまり,ぼくの今回のルートスキャンのような作業にはPolycam ProのSplattingには対応していないということである。狭い範囲のDOIでの詳細なスキャンを実施すれば結果が出る可能性はもちろんあるだろう。

図J 

図K

図L

 

図M

図N

図O

4. Scaniverseによるスプラッティング

 図Aのルートを24分間かけてsplattingした。スキャン中(図P),more slowlyと忠告が出続けた。最後まで中断されることはなかった。more slowlyの忠告は図Pのスキャン中出ていたがスクリーンショットには映っていない。画像はゴミのようなものが続いた。以前タニハでこれを使った際にはかなりゆっくりとスキャンして一応,出力できたが,スキャン中は雲のようで映像がはっきり見えなかった。Polycamの映像がクリアなのと大きな違いである。無料だからしかたがないとは思う。
 図Qはスキャンが完了した際のスクリーンショットである。Process Laterを選んだ。Polycamと違って,Wi-Fiでアップロードせず,iPhone内で完結する。ここがメリットとも言える。自宅に戻ってProcessingを実行したが4回ほどトライしたがerrorから抜け出せない。ルート測量には適さないということだ。

図P

図Q

図R

以上,2024年10月28日。

5. Scaniverse によるLiDAR測量のおよそ

 結局,ルート測量はLiDAR測量に限定されるということらしい。写真測量はかなりの情報が必要なのに比べて,LiDAR測量は,直接ベクトル長や方向角を捉えるために,空間把握が「高速」なのである。図1ではルート途中に打越池南縁のROIが含まれるのであるが,ルート測量の途中なので,このROI部分を詳細にスキャンするという発想もあったが,それは誤りである。ルートスキャン測量と詳細スキャン測量は別に考えた方がよい。前者で公共の基準点4点とROI内のLabels 1 to 4をスキャンして,公共の基準点4点からROI内のLabels 1 to 4の座標点を求めて,その上で,Labels 1 to 4を含むROI内の詳細スキャンを実施するのが適当ということになる。

 10月26日にLiDARによってルート測量をScaniverseで実行した結果をアイフォーンのライブラリー表示でまずはみたいと思う。このルート測量結果は図1のものである。図Sでは,①つまり3-61基準点を出発しこの①に戻っている。同じ3-61基準点の3D距離が26.38mにもなっている。ズレているのである。図Tではこの3-61基準点が二つのヘビの鎌首のように見える。巣閉距離も垂直距離もずれている。図Uでは,図の下方から上方に進み,上端付近では,プロムナードが終わってバス停があるアスファルトの歩道に出ている様子が見える。バス停そばのLabel 4付近は立ち止まってスキャン時間が多少長くなっているので,プロムナードがはっきり見えるが,そのあと早足で歩道に出ているのである。画像が切れていても,場所情報は繋がっている。これはLiDAR測量の貴重な長所とも言えるのである。
 図Vはライブラリー表示の画面であるが,下部の生垣などの画像はScaniverseのSplattering測量の結果が見えるが,画像はLiDAR測量と違ってはっきりと再現されている。詳細スキャン測量はSplatteringが有効とわかるだろう。

図S

図T

図U

図V

6. LiDARルート測量結果をグリッド平面座標に

 これからの作業は,macでこのWebページを作成,計算処理も実施しつつ,Windows PCでCloudCompare操作をしてその作業過程のスクリーンショットを撮影し,それをmacでのWebページ作成に利用するので,macとWindows PCの間をサーバー接続環境にしている。

6.a 基準点を使ってアラインメント

1 iPhoneアプリScaniverseのライブラリーに格納されている筆者作成のライダー測量結果を点群系のPLYファイルでmacにアップロードし,共有設定したWindows PCにも保存する。ファイル名:IshimaruOct26-2024byScaniverseLiDAR.ply

2 Windows PCのCloudCompareから,このPLYファイルを開く。

3 macで表6のルート5地点①,④,⑥,⑦,①の (Y,X,Z)の基準点3D座標値から,csvファイルを作成して,両共有フォルダーにコピーしておく。次の表6にまとめている。

4 CloudCompareのコンソールのPLYクラウド上で,ルート5地点のpoint list pickingを実施する。point list picking:Tools > Point list picking,またはトップのツールバーの左から5番目のアイコンをタップ。

図Wは実行後の結果で,Points #0 to #4の座標値はcsvファイルで出力している。

図W ScrnS 791

5 上記のpoint list pickingの結果のcsvファイルを元に次の表7を作成した。

 表6と表7に関わってハンドホール蓋上面と基準点真鍮鋲の関係を解説する。LiDAR測量では基準点の真鍮鋲はハンドホール内にあって測量プロセスの観点から直接には測量できない。そこで基準点の高度として,真鍮鋲ではなく,ハンドホールの蓋上面の高度を使用する必要がある。LiDAR測量結果の点群数はたとえば図Wでは429,321点に及ぶので,たった4点の基準点の高度を変更した方が適切なのである。そのため,表6の「”reference entities” csvファイル出力用」に真鍮鋲とハンドホール蓋の高度関係を組み込むことになる。

6 さて,align実行の下準備は完了した。LiDARルート測量結果クラウド(IsimaruOct26-2024byScaniverseLiDAR)と表6の最後の3列からなるcsvファイル Ishimaru5points_revised.csvをcontrolキーで同時に選択して,tool > Registration > Align (point pairs picking),またはトップのツールバーの左から9番目のアイコンをタップする。
 図XではSelect aligned entitiesのパネルでIsimaruOct26-2024byScaniverseLiDAR – Cloudを選んでいる。そして,OK。

図X ScrnS 795

 次に,aligned entities (クラウド上で基準点対応点をあらかじめ求めた表7の5点)のまずは第一点A0を入力する。コンソールにはそのA0点が表示される。

図Y ScrnS 797

 次にはreference entities(基準点に蓋上面からの真鍮鋲までの深さを足したZ’値を反映した表6の5点)のA0に対応するR0をまずは入力する。同様にA4, R4まで入力したのが図Zである。alignボタンをタップした結果が図Zである。

図Z ScrnS 801

 図Zは図Aの3-61地点付近を拡大表示している。R4とR0は同じ位置に重なっている。A0はこのルート測量の出発点で,A4は終着点であり,いずれもR4(R0)基準点であった。できれば,A0とA4を一致させたかったが,難しいと考え,aligned entitiesのA4とreference entitiesのR4を表右手の❌をクリックして削除した。その上で,alignボタンをタップした。その結果が図AAである。

図AA ScrnS 802

 ☑️した結果が図BBである。RMSが見える。A点とR点がこの画面でほぼ近いことが見える。

図BB ScrnS 803

 図CCはこのアラインメントを終了した後の表示である。この図は図AAとほぼ一致している。この程度の処理で筆者の目的からすると十分である。緑色のPoint #付きのアイコンはクラウド上のピッキング記録であり,そのそばにある黄色の四角はレファランスの基準点である,ペアで最も近接しているのは⑥3-40であり,最も離れているのは,Point #4と3-61である。3-61で最初に利用した①のPoint #0は,#1, 3とほぼ同様の距離である。

図CC ScrnS 805

以上,2024年11月2日。

6.b Labels 1〜4の平面直角座標系座標値を求める

 前節では,基準点3-61, 3-42, 3-40, 3-41の4地点とそれらを通過するLiDARルート測量を実施し,LiDARルート測量結果の点群point cloudをグリッド座標系に載せた。その点群上のLabels 1〜4をCloudCompareのpoint list pickingツールを利用して,グリッド座標値,つまりは,平面直角座標系6の座標値を求めることになる。LiDARルート測量の際には前もって,打越池南縁のプロムナードにあたるDOIに,Labels 1 to 4の4点のラベルを設置している。それゆえ,この点群には,4ラベルが記録されているのである。

 ルート測量の際に早足で歩いたために,点群密度はDOIのスキャン作業結果と比べてかなり低くなっている。4点個々のラベルとその付近では立ち止まってスキャンしていたが,ラベルが小さいこともあって,点群では設置場所を特定できなかった。それゆえ,メッシュ系のFBXを観察してラベルの位置を特定し,点群と対比して点群上でも特定することができた。

 FBXのコンソールでの表示は,ラベルが明瞭に見えるところまでズームインし,スクリーンショットを撮り,その後,ズームアウトして,さらにスクリーンショットを撮影した。FBXの像とPLY点群の像のアウトラインはかなり異なっているので,ズームアウト像は必須である。
 この経験からの反省点を次に挙げる。ラベル紙は7cm弱四方でクラウドで見るには小さすぎた。ラベル用としては,A4用紙の短辺23cm四方が適当かもしれない。その上でラベルとその周辺のスキャンに時間をかけた方が良いようだ。このような反省点はあるが,DOIでのラベルの平面直角座標系6の座標値は図DDのように得ることができた。

図DD ScrnS 814

 図DDの右上のpicking listの#番号は5〜8になっている。前節での作業を継続させたためである。

 次の表8にはこの作業での結果をまとめている。平面直角座標系6での作業では,座標は(Y, X, Z)であるが,CloudCompareでのグリッド座標系では,図DDの表にあるように,(X, Y, Z) になっている。
 前節では,4基準点を使って(スキャンに含めて),LiDARルート測量結果である点群をグリッド座標に載せたのであるが,DOIでのLiDARスキャン結果はこの4ラベルを使って(スキャンに含めて),点群をグリッド座標に載せる必要がある。

7. 点群からメッシュにむけて

 ここでは,国土地理院の基準点成果を利用して,調査フィールドDOIのLiDAR測量結果をグリッド座標系に載せる手法を示した。LiDAR測量の結果の景観再現性は,外観的には点群よりもメッシュの方が良い場合があるが,CloudCompareではver. 2.6以降,メッシュFBXの場合もこの種の手法が可能となっているが,異なる座標系での作業を試したことがなかった。マニュアルにもその実例はない。

 点群の平面図(Set top view)は,測量を実施した場の平面図と一致しているが,メッシュでは正面図(Set front view)の裏返し 3D image flipping が点群の平面図と対応している。そのため,pickingして表示される座標値は,(X, Z, -Y)になっており,基準点座標との関係を求める場合,注意が必要だと思うが,挑戦して,さいわい成功した。

以上,2024年11月2日。

7.1 メッシュと点群の座標軸の違い

 点群は感覚的によく合っているが,メッシュは開いたところがSet top viewにあたっているが,正面図だ。

図EE ScrnS 816

 Set top viewにすると,図FFのように,平面図になる。

図FF ScrnS 817

 裏返したのが図GG。右下の座標軸を見ると,水平右方向がX,水平下方向がZ,そして垂直方向がYだ。

図GG ScrnS 818

 点群とメッシュの比較をするために,図Aの各基準点対応のメッシュ座標値をPoint List Pickingしてみよう。Point #1 to #4だ。図Aの,①,④,⑥,⑦,①だ。図819〜823を見ると,ハンドホールがはっきりと見える。点群と全然ちがう。

図819 ①

図820 ④

図821 ⑥

図822 ⑦

図823 ①’

7.2 メッシュ上でのPoint List Picking

 Point List Pickingの結果はcsvで出力した。その結果を表9に示すが,前述のように,点群 (Y, X, Z’)は,メッシュ (Y, Z’, -X) になっている。鬱陶しいだろうが,表6も追加したものを次に。

7.3 アラインメント align (point pairs picking)

 メッシュentitiesと基本基準点 entitiesとの対応関係をとることになるが,点群 (Y, X, Z’)と違って,メッシュ (X, Z, -Y)とする必要性がある。表9はalign entitiesであるが,新たに表10を用意した。これでアラインメントの準備はできた。これが正しいかどうかはやってみないとわからない。

 図824でアラインメントの準備完了。コンソールに,A3, R3が見えるが,align entities,reference entitiesに,A0〜A3, R0〜R3が入力されたことを示している。

図824

 alignの結果を図825に。大成功だ。

図825

図826

 RMS: 8.86だ。点群と同様だな。クラウドでの基準点位置極めは,メッシュの方が精度が高かったが,RMSはあまり変わらない。

 図827ではalignの表は無くなった。

図827

 さて,ROIに4カ所の二次的基準点を設けるべく,Lavels 1〜4に,picking実施。点群ではラベルが見えなかったが,はっきりとその中心も見える。さいわいなるかなあ。

図828 Label 1, Piciking #5

図829 Label 2, Piciking #6

図830 Label 3, Piciking #7

図831 Label 4, Piciking #8

 図833では,このmeshの座標軸の配列を理解するのにいいかと思う。コンソールの右下隅に,R(X)G(Y)B(Z)軸が見える。

図834 Set view to “back” isometric

 表11では,平面直角座標系6の座標軸とメッシュの座標軸の対応関係を示している。

 以上,アラインメント機能を使って,グリッド座標系(平面直角座標系6)とメッシュの座標軸との関係を明確にすることができた。メッシュを使ってこそ,基本基準点をROIに二次的基準点を作成することができるし,環境把握も妥当なものとなることがわかった。

以上,2024年11月5日。




3D Gaussian Splatting

はじめに

 いやー,驚いた。ChatGPT ver. 4の契約をするかどうか,ネットサーフィンしていたら,このgaussian splattingテクノロジー情報をみつけた。ぼくのiPhoneのScaniverse 3D map scan アプリにもすでにgaussian splattingが実装されている。このgaussian splattingの技術が実用化したのは昨年(2023年)夏とま新しい。原理を理解するのは今のぼくには難しい。これまでのフォトグラメトリなら簡単な幾何学に過ぎないが別世界の技術だ。解析速度もより速く,解像度もより高い。LiDAR測量に拘ってきたが,通常の写真や動画の方がより有効な気配だ。

メモ: SfM (Structure from Motion)またはフォトグラメトリ(Photogrammetry)とは、被写体をさまざまなアングルから撮影し、そのデジタル画像を解析、統合して立体的な3DCGモデルを作成する手法。

 ぼくが最も信頼する3D scan soft, Scaniverseにまずは眼が向いたが,最強だが有料のpolycamが劇的に使い易くなったようだ。他にLuma AIもある。全部,iPhone11以降に対応している。以下,ぼくの知りたい情報をピックアップしつつ,学びたいと思う。そして,実際に試した結果もここに報告したい。

以上,2024/04/12。

1. Niantic Lab

Niantic (ナイアンティック)って知らなかった。Google社内のスタートアップ企業で,Googleから独立したらしい。そして,3D ScanningソフトScaniverseを吸収したようだ。かつては有料だったが無料化されて,先進的なGaussian Splatting機能が付加されてもなお,無料のままだ。ぼくのiPhone 12 Proのアプリは,Apr. 1, 2024にもアップデートされていたNianticはGoogle EarthもGoogle mapも開発したGoogle内のティームだったようだ。

2. Scaniverse

Scaniverseが3D Gaussian Splattingをサポート 〜 スマホからフォトリアルな3Dシーンを迅速に作成 〜 2024.3.20 を参照して欲しい。LiDARではなく,ムーヴィーか,シーンを重ねて撮影した写真で,Gaussian Splattingが可能だ。次のアドバイスがある。

メッシュ生成とスプラット生成をどのように使い分けるべきですか?

スプラットが適している場合:
照明と反射を備えたフォトリアリスティックな結果を求める場合
背景を含む全体のシーンをキャプチャしたい場合

メッシュが適している場合:
スキャン結果を他の3Dソフトウェアやゲームエンジンで使用したい場合
背景のない独立した3Dモデルを求める場合
正確な測定が必要な場合
明確に定義された被写体がないシーンをモデル化したい場合(例えば、建物内のすべての部屋を含むスキャンなど)

 写真の3D合成は,他のアプリでは,有料で,企業のサーバーにアクセスする必要があるが,このScaniverseでは,iPhone内で実行することができる。しかも無料だ。まだ,完成形ではないと,アプリ開発者が考えて,ユーザーの意見などを吸収して行こうとしている。

Scaniverseでスプラットをサポートしている端末を教えてください。

Gaussian Splattingは現在iOSでのみ利用可能です(iPhone 11以降、比較的新しいiPad[*])。Android版の開発中ですので、お待ちください!

[*] iPad: 第9世代および第10世代、iPad Pro 11インチ: 第3世代および第4世代、iPad Pro 13インチ: 第5世代および第6世代、iPad Mini: 第6世代、iPad Air: 第4世代および第5世代。

3. Scaniverse with Gaussian Splattingの使い方

3.1 Gaussian Splatting

Scaniverseが3D Gaussian Splattingをサポート 〜 スマホからフォトリアルな3Dシーンを迅速に作成 〜 2024.3.20 の説明を引き続き,次にコピペ。

3D Gaussian Splattingとは何ですか?

Gaussian Splatting(ガウス・スプラッティング)は、3D空間の再構成における革命的なラスタライゼーション技術であり、既存のNeRF(Neural Radiance Fields)およびフォトグラメトリ技術の課題を組み合わせて克服します。一連の画像を使用して、Structure from Motion (SfM) が点群を推定するために使用されます。各点はガウスに変換され、位置、方向、スケール、色、透明度などの一連のパラメーターを含みます。簡単に言えば、スプラットは3Dシーンの照明条件や材料の特性をまったく新しい方法で再現することができ、これにより、光沢のあるまたは反射する表面をより良くレンダリングすることが可能になります。

この技術は、計算コストを低減し、リアルタイムレンダリングとフォトリアルな品質の可視化を可能にします。そして、Scaniverseでは初めて、サーバー処理なしでデバイス上で直接利用できるようになりました。これは、あなたのスキャンが以前よりもさらに鮮明で生命感あふれるものになることを意味し、インターネット接続なしで全てをあなたのスマートフォンから行うことができます。

Scaniverseで高品質のスプラットを作成するには?

Scaniverseアプリをダウンロードし、赤い「New Scan」ボタンをクリックして「Splat」を選択するだけです。その後、画面上の指示に従って、あらゆる角度からデバイスの360度ビデオを撮影してください。良好で均一な照明条件の下、均等な動き、最小限の背景ノイズ、できるだけ多くの角度から作成されたスプラットが最も良い結果を生むと考えています。

新しいScaniverseでは様々な対象物に対してスプラットを作ることができます。ただし、現在、日常的なオブジェクトには適していますが、より大きなエリアのスキャンについてはアプリのパフォーマンスを改善中であることをご了承ください。

ビデオキャプチャが完了すると、キャプチャしたデータを点群として表示します。ビデオに追加するためにスキャンを再開するか、または「Process」を押して少し待つと結果が表示されます!

最新のアップデートでは、さらに「Enhance」を押すことで、高品質な結果を得るために追加の1分間隔でモデルをトレーニングすることもできるようになりました。

生成したスプラットで何ができますか?

デバイス上で直接スプラットを閲覧し、個人のライブラリに保存したり、他のScaniverseユーザーと共有したり、友人や家族と共有するためのビデオリールを作成することができます!また、PLY形式でモデルをエクスポートすることも可能です。

Gaussian Splattingの3Dソフトウェアサポートはまだ限られていますが、Unityには無料のアドオンが、Unrealでは商用プラグインが利用可能です。

メッシュ作成には提供されているトリミング機能は、スプラット向けにも近く提供される予定です。

過去にスキャンした結果からスプラットを生成できますか?

はい、Rawデータが保存されていれば可能です。そのためには、対象のスキャンを開き、右上隅にある「…」アイコンを押してください。「Reprocess Scan」を選択し、その後で「Splat」を選択します。

私がスキャンしたデータはクラウドに送られますか?

いいえ、すべてのデータはデバイス上でローカルに処理されます。スプラットを共有することを選択しない限り、データはデバイスを離れることはありません。

Gaussian Splattingには協力なGPUが必要だと思っていましたが、スマートフォンだけで実現できるのですか?

可能です。現代のスマートフォンは、ポケットサイズのスーパーコンピュータです。いくつかの注意深い最適化を行うことで、簡単にスプラットのトレーニングを処理することができます。

以上,1:42,Apr. 10, 2024記。

3.2 Gaussian Splattingのために

 今後記述予定。

4. Polycamがより魅力的に

 新たなPolycam情報を見るまでは,Polycamには無料版が無いがリーズナブルの価格で,LiDAR測量については,Scaniverseよりも,より広域であっても測量誤差が抑えられる,という認識であった。ところが,現在では無料版が提供されており,ここで関心の高いGaussian Splatting機能が追加されている,ということである。

 今後,特にこのGaussian Splatting機能に関して,ScaniverseとPolycamを実際に使って試したいと思っている。

4.1 アプリPolycamとは

 iPhone用アプリは,アップルストアで。Polycam 3D Scanner, LiDAR, 360 4+ Easy 3D modeling & floor plans
Polycam Inc. #73 in Photo & Video 4.7 • 31.3K Ratings Free Offers In-App Purchases,とある。

 これをダウンロードして開くと,まずは試してみる機会があって,LiDARではなく,写真を選んで,ぼくの愛(あい)椅子を撮影してみた。自動的に連続シャッターが始まって,計126枚が撮影され,Gaussian splatを選択して,アップロードした。かなり時間がかかる(1時間ほどか)。これを実行すると有料になってしまうのではとちょっと不安だった。

 図1,2のように,無料版では3D objectは取り込めない。この場合,ファイルサイズは31.6MBとなる。

図1 Polycamのサーバー内で3D表示された愛椅子

図2 取り込むにはPro契約が必要

 図1の背もたれのトップは薄い茶色のガスがかかっているように見える。夜間の二カ所からの天井灯による照明なので,このような現象が生じているからかも知れない。ここでは点群PLYで出力したかったのであるが,図1のような再現性は望めないであろう。いますぐProにアップグレードする気になれない。

 図3,図4は,サーバー内の3D オブジェクトを観察したものである。図3の脚部では写真撮影が失敗した結果であろう,オブジェクト情報が粗い。図4の脚部は明瞭で問題が無いが,背面上部には図1で述べた薄い茶色のガス状のものが見える。狭い空間で撮影していることや照明の真下で撮影しているために,撮影者から乱反射があった可能性がある。図4の脚部は明瞭な像を結んでおり,図3でと同様,iPhoneを逆さに撮影したのであるが,前面がテーブルとの間と異なり,比較的窮屈でない場での撮影が功を奏したものであろう。この撮影は,椅子の回りを二回ほど廻ってのものであった。

図3 椅子の前面

図4 椅子の背面

 試行錯誤しつつ,撮影を繰り返すことで,よりよい3Dオブジェクトが得られるであろうが,この程度のターゲットでも解析に1時間かかるところがなかなか,厳しい。

図5 価格区分

 図5では,ぼくに関係するのは,FreeとProだ。Proの場合は,Polycam for Educationが適用される可能性があるが,まだ調べていない。

1 LiDAR利用ではFreeとProに違いがない。Freeの,Unlimited LiDAR,は,このPolycamでは極めて重要だ。というのいは,これまで使用してきた3D scannerアプリでは使用時間の制限がある。Scaniverseも10分間以降だったか,旧に解像度が落ちるが,PolycamではFreeが無い時代だが,解像度が落ちなかった。Freeの時代になっても,これは変わらないだろう。

2 FreeとProの違いは,静止画や動画に現れている。Freeでは静止画150枚限定だ。Proでは,写真測量で静止画2000枚まで,gaussian splattingで1000枚まで。他のアプリでは有料版でシャッター撮影は自分で実行してゆくタイプで,静止画200枚までだった。静止画2000枚は驚くべきである。ただ,自動でどんどん撮影される?ので,2000枚と言っても,droneなどとは違って,歩いての場合,撮影範囲は広くはないだろう。

3 Proでの,3Dモデルの出力ファイル形式は,obj, fbx, stl, pdfなどの12+種類がある。アプリの説明には,メッシュとしては,OBJ, GLTF, DAE, FBX, USDZ, STL,点群ファイル形式としては,DXF, PLY, XYZ, LAS, PTSが列挙されている。ぼくがよく使うのは,点群ファイルとしてPLY, メッシュファイルとしてOBJ, FBXである。
 フリーヴァージョンで用意されている出力形式は,GLTFファイル形式だけだ。スキャン時にはProである必要性が無いことが多いだろうが,iPhoneに保存された3Dオブジェクトの出力の際には,Pro契約が必要になる。

4 Proの価格は,$18/月,$100/年,だから,月3000円か,年15000円。ぼくは使用頻度は極めて低いので,使い始める時に月単位で契約すれば良い,と思う。

4.2 使用実例に学ぶ

4.2.a 洞窟測量の実例コンテンツ

 洞窟測量におけるPolycam(3Dスキャンアプリ)操作マニュアル 関連の書誌は次のよう。
小竹祥太,林田敦,柏木健司,2023. iPhone Pro によるLiDAR 測量を活用した洞窟調査: 七久保の道穴(群馬県下仁田町青倉川上流)を例として.群馬県立自然史博物館研究報告, No. 27, pp. 171-179.
Kotake, S., Hayashida, A. and Kashiwagi, K., 2023. Cave survey using LiDAR measurements with the iPhone Pro: an example of Nanakubo-no-michiana Cave(limestone cave)in the upper reach of Aokura River, Gunma Prefecture of central Japan. Bull. Gunma Mus. Nat.Hist., No. 27, pp. 171-179.

<本論文の付ファイル01,02,03は群馬県立自然史博物館Webpageの下記アドレスで公開>
付ファイル01 洞窟測量におけるPolycam(3D スキャンアプリ)操作マニュアル
http://www.gmnh.pref.gunma.jp/wp-content/uploads/bulletin27_17_01.pdf
付ファイル02 Blender 3.0.1 で使うアドオンCave Mapper マニュアル;導入編
http://www.gmnh.pref.gunma.jp/wp-content/uploads/bulletin27_17_02.pdf
付ファイル03 七久保の道穴の3Dウォークスルー動画
http://www.gmnh.pref.gunma.jp/wp-content/uploads/bulletin27_17_03.mp4
 以上の3urlは,暗号化されていないので,Google関連ソフトなどではアクセスできない。Firefoxでは問題無くアクセスできる。

付ファイル02 Blender 3.0.1 で使うアドオンCave Mapper マニュアル;導入編,に関連しては,https://github.com/CaveMapper/CaveMapper/wiki/02_CaveMapper_Manual
に公開されている。これを見ると,
Blenderでの,モデルの位置合わせ,位置合わせしたモデルの一体化,データの軽量化,断面の作成,からなる。ぼくはCloudCompareでこの種の実行はしてきたので不要ではある。

 なお,ネット上で,修正前の群馬県自然史博物館のコンテンツのurlをただコピペした場合,次の画面が出てきて,ブラウザーを閉じることができなくなった。macのパワーキーを二回押し込んで,何とか離脱。大きな声と動画でかなりビビった。初めて経験した。アップルからのリアクションとは思えない下品な雰囲気が漂っていた。群馬県立自然史博物館のWebサイトはやられているようだ。

4.2.b これから参考になる点は

1 iPhone Pro によるLiDAR 測量を活用した洞窟調査: 七久保の道穴(群馬県下仁田町青倉川上流)を例として——————
 これを見ると,iPhone ProのLiDAR測量は便利やなあ,という感想,ということで,ぼくの体験や感想も同じだ。
p. 172右段の方法,では,「iPhoneに搭載されたLiDARは,搭載カメラが撮像可能な明るさ環境で,スキャン対象の色情報の取得が可能である.洞窟記載の基礎データとなる洞内測図を作成するため,洞内ではiPhoneにヘッドライトを固定してLiDAR測量を実施した」とあった。「色情報の取得」というか,明るくないとLiDAR測量はできない。
p. 174左段では,「LiDAR測量に要した時間は事前の下見を含み約15分であった.七久保の道穴周辺の地形断面の概略を把握するため,沢床から町道,洞口,および洞口上方の尾根まで,LiDAR測量を分割して実施した.全体を先ずは見渡し小崩壊地形を避けてルートを設定し,沢床から斜面を傾斜方向に登りながら測量した」,とあって,洞穴外の地形もLiDAR測量されている。ぼくもこの種の作業を実施したが,iPhoneのLiDAR測量能力は凄い。
p. 177の議論では,「iPhoneを利用したLiDAR測量は,iPhoneそのものが小型かつ片手で操作が可能であり,狭い空間を含む洞窟での測量において高い利便性を有する.一方,LiDAR測量に伴うiPhoneのバッテリーの消耗は大きく,複数の洞窟を対象とした筆者らのこれまでの試行では,3-4時間程度でバッテリーが低レベルの状態となった.考えうる現時点でのバッテリー対策としては,予備のモバイルバッテリーの準備,ないし車が近くにある場合はシガーソケットからの充電などが考えられる.また,複数名でLiDAR測量を行う場合,スキャンする空間を計画的に割り当てることで,一定規模の洞窟空間をカバーするLiDAR測量が可能である.計測範囲の5 mを越える平面的に広い空間については,洞床と洞壁を含みくまなく歩いてスキャンすることで,データ取得は可能である.一方,天井の高い空間については,自撮り棒や伸縮性ポールの先端にiPhoneを取り付け,長さと方向を調整しながらLiDAR測量を行うなどが想定できる」とある。まあバッテリーは安いので複数持参すれば問題はない。

 補助的な測量や公共基準点からの位置づけなどがされているのではと思ったが,無かった。

2 付ファイル01 洞窟測量におけるPolycam(3D スキャンアプリ)操作マニュアル——————
 この報告にはページ番号が振られていない。
 4枚目の,4.Polycam の応用的な使用方法,で,次の記述は参考になった。
「4.1 スキャンの一時停止および再開
Polycam のスキャンには一時停止機能がある.スキャン中を示す赤色の録画ボタンの右横の一時停止ボタン(図1 の④)をタップする.すると,画面下に「capture paused」が表示される.再開には「○に▷記号の入るアイコン」をタップし,終了時には赤い録画ボタンをタップする.一時停止機能は,狭洞など円滑な動きが難しい空間での測量に効果的で,一時停止後にスキャンし易い体制に体を入れ替え,スキャンを再開する.
Polycam には,スキャン終了後にスキャンを再開できる機能がある.再開する場所ないし対象物付近で過去のスキャンデータを開き,画面下部の「追加」アイコンをタップする.「追加」アイコンが見当たらない場合は,画面下を左右にタップしてアイコンを移動させる.過去にスキャンした3D データとスキャン中の3D 形状に,相互にオーバーラップする形状が認識されると,「relocalization successful! You can now proceed to scanning new areas」と表示され,白い録画ボタンを押すとスキャンが再開される(図1 の⑩).洞窟の分岐に加え,スキャンし忘れた箇所の補完に対応できる.なお,「追加」機能はデータサイズの増加につながり,ファイルサイズが大きくなることでアプリのエラーや強制終了が生じることがある.また,過去の3D データにオーバーラップする範囲が認識されない場合,「rescan the place you originally scanned to extend, targeting visually rich areas」として「iPhone を動かして開始」という表示が画面上に示される(図1 の⑪).我々の数度の経験では,同じ場所のスキャン再開は問題なくできたので,iPhone を様々な方向に動かしてもこの表示が消えない(再開に移行しない)場合は,スキャン対象である場所ないし対象物を間違っているか,同じ場所であっても改変が生じた可能性を考慮する必要がある.
4.2 スキャンの分割
一度にスキャンする空間の大きさやそれに伴うデータ量に依っては,取得データのサイズが大きくなり過ぎることでPolycam の処理にエラーが発生し,またはPolycam が強制終了することがある.どの程度の広さや長さを目安に,スキャンする空間を分割するのが適切ないし最適であるかは,調査目的に係
るスキャンの精度に紐付けされる取得データ量に左右され,数字で示すことは難しい.なお,一定以上の空間を対象に分割してスキャンする際,後の3D モデル結合処理の際の利便性を考慮して,目安として2-3 m 程度の重複を必ず確保しておく必要がある.また,特徴的な形状が含まれる空間で重複を取ることが望ましい.メモ:この赤字化した部分だが,ぼくの場合,紙で作ったマーカー(四分割対角黒白塗)を貼っている。
4.3 データの再処理
Polycam では,処理済みの3D データを繰り返し処理できる機能が備わっている.処理済みの画面下に並ぶアイコンのうち,左端の「処理」をタップ(図1 の⑦)すると,処理画面に移動する(図1 の⑤).ここで右下の「カスタム」をタップすると,深度範囲,ボクセルサイズ,簡略化,自動クロップ,LOOP CLOSURE について,処理条件を調整できる(図1 の⑫).処理条件を調整後に「処理」をタップすると,画面中央のダイアログボックス内に「Reprocess capture? This will overwrite the existing capture and erase your edits.」と表示される.以前に処理した3D 画像データは上書きされるので注意が必要である.処理を続ける際は「はい」をタップする.すると,新しい条件での処理が開始される

 赤字化した内容のうち,特に,「Polycam のスキャンには一時停止機能がある.スキャン中を示す赤色の録画ボタンの右横の一時停止ボタン(図1 の④)をタップする.すると,画面下に「capture paused」が表示される.再開には「○に▷記号の入るアイコン」をタップし,終了時には赤い録画ボタンをタップする.一時停止機能は,狭洞など円滑な動きが難しい空間での測量に効果的で,一時停止後にスキャンし易い体制に体を入れ替え,スキャンを再開する」,は気が付かなかった。

4.3 Scaniverseとの比較

 すでに両方を使った当方の感想をぼくは3年ほど前にまとめているが現在,大きく変わった。他のサイトでもScaniverseとPolycamの比較があるが,やはり古い。現在のところ,現在の両アプリを比較したサイトは存在しないようだ?この節については,今後,使用した結果を示したいと思う。

 

5. Luma AI

 Luma AIは,iPhone用のフリー3Dスキャナーである。NERFで3Dモデルを作る「Luma AI」無料での使い方や商用利用の可否を紹介! には,使用法が簡潔に書かれているが,小さなオブジェクトに限定されている印象だ。「Luma AIはNeRF(Neural Radiance Fields)という技術を使って3Dモデルを作成しています。NeRFとは簡単にいうと、複数枚の写真を用いて機械学習から三次元の状態に復元する技術のことです。Luma AIはNeRFを用いたスマートフォンからの3Dキャプションアプリの筆頭といえます。」,とあり,写真測量法やGaussian Splattingとも異なる。

 他のサイトも見たが,ぼくにとっては魅力を感じなかった。

以上,2024/04/21。

 




Specs desired for CloudCompare

はじめに

 狭い範囲のdrone測量3Dmapをスムーズに観察できない。mouse製のGeforce(NVIDIA社が設計開発しているGraphics Processing Unit (GPU) のブランド名)を銘打っているゲームマシーンであるが。メモリは64GB,もう十年ほど前に購入したものだ。人生が残り少ないので,新たな投資はどうかと思うが,ちょっとネットでみた。

1. 参考になったウェブページ

1 フォトグラメトリ・点群処理向けマシン選定のポイント では,代表的なソフトウェアと要求スペック,の欄には,

 CloudCompare 基本構成 高度な構成
メモリ 32GB以上 64GB以上
CPU 4コア以上 8コア以上
GPU メインストリーム品 高性能ビデオカード(Quadro等)

まあ,ほぼこれを満たしている。でもダメダメなのだなあ。

2 mouse社の,グラフィックボード(グラボ) グラフィックス グラフィックカード ビデオカード GPU これは次回,購入の際に参考になるだろう。

3 Handling large areas, Viewport Performance,

途中。嫌になった。中断。

以上,Apr. 6, 2024記。




光波測量によるネットワークをdroneマップに投影する方法 how to match distance-measured points from different setting bases with a drone-measured map

1. 課題

 フィールドワークでの測量に問題が生じた。2020年測量のdrone3Dマップと,2018年の光波測距儀を使って得られたネットワークが対応しない。droneマップは光波測距儀で得られたgcpとの対応関係が取られている。そのgcpを使ってdroneマップをコントロールしているので問題はない。この場所には自衛隊のレーダーサイトがあって,コンパスが正常に作動していない危惧があり実際に西偏角度が想定されるものと異なっている。

 光波測距儀の設置点は14カ所で個々の設置点で幾つか測点があって,前視と後視で,位置関係を一応把握したのであるが,手抜きしてGPS計測値を光波測距儀の地球座標値として使用した。これが今抱える問題の元凶だ。そこで,2018年の測量時の測点の現場写真からdrone3Dマップと対応できる測点を利用して,ls_1〜_14の光波測距儀の設置点を求めたいと思う。

2. ls_2の設置点を求める

 ls_2の放射星形区で3Dマップ上での2カ所を現場写真から特定できた。その一つは,図1のls2_20であり,もう一つは図2のls2_11である。後者は光波ネットワークとも一致している。それを図3に示す。

 ls2_20の方は,図2に示すようにかなりズレている。図1の現場写真ではa-rotの下に見えるが,図2では3mほど内陸にずれ込んでいる。

 求めたいのは図1の,赤字で「かわぐち」としたls_2の設置点の3Dマップ上の座標値である。

図1 ls_2放射星形区の位置特定
図2 CloudCompare画面上のls2_20, ls2_11の位置関係 左手が海側
図3 この図の中央に,ls2_11が

3. ls_2の3Dマップ上の座標値を求める

3.1 計算法

 当初,3D空間での三角形平面を想定して計算すべきと考えた。しかしながら,熟考するとその必要性が無いことが明らかになった。光波測距儀の構造は台座を水平にして,光波測距儀が装着されている。コンパスで磁北をゼロ度として,台座上を時計回りに回転して,光波測距儀をターゲットに向けて,直線距離を求める。その測距儀の水平角と垂直角を記録する。これだけだ。今,問題なのは3Dマップでのls_2設置点を決定することであり,2D空間に過ぎない。

 連立方程式は簡単にできる。文字表現の簡略化したい。現在の設置点をO’点として,3Dマップに載せうる設置点をO点(Xo, Yo)とする。ls_2放射星形区の2点のうちの1点について,現在の位置をA’点として,3Dマップに載せうる設置点をA点(Xa, Ya)とする。他の1点については,現在の位置をB’点として,3Dマップに載せうる設置点をB点(Xb, Yb)とする。

 光波測量の結果から,次式が成り立つ。赤字表現しているのが未知数である。他の数値はすでに得られた計測値または3Dマップから読み取りできる。
|OA|=|OA’|=√((Xa-Xo)^2+(Ya-Yo)^2) ⋯⋯⋯⋯⋯式(1)
|OB|=|OB’|=√((Xb-Xo)^2+(Yb-Yo)^2) ⋯⋯⋯⋯⋯式(2)
 上記の計算式が成り立つ理由は次のようである。光波測量の結果に基づいて位置ベクトルを3Dマップにプロットするのであるが,ls_2の設置点が決まらないだけであって,その測定結果のうち,水平距離,垂直距離,さらには視準線の間の挟角は変わらない。

 式(1), (2)の連立方程式から,(Xo, Yo)の解が2個現れる。A点を中心として半径|OA|の円と,B点を中心として半径|OB|の円との交点は2点形成されるからである。3Dマップを見ている作者からするといずれが正解かは自明ではある。

 式(1),(2) は何れも円の方程式である。両円は交わって,その二つの交点を通る直線が成立する。その観点で整理したいと思う。式(2)と(3)はそれぞれ,次のように整理できる。
(Xo– Xa)^2+(Yo– Ya)^2 = |OA|^2 ⋯⋯⋯⋯⋯式(3)
(Xo– Xb)^2+(Yo– Yb)^2 = |OB|^2 ⋯⋯⋯⋯⋯式(4)

 二つの交点を通る直線の式は,式(3) – 式(4),で次のように,得ることができる。
2(Xa-Xb)Xo + 2(Ya-Yb)Yo +〔(|OA|^2 – |OB|^2) – (Xa^2 -Xb^2) – (Ya^2 – Yb^2)〕= 0⋯⋯式(5)

 両円の交点の計算は,式(3)と式(4)の何れかの円と式(5)の直線との交点を求めれば良い。ここでは,式(3)の円と式(5)の直線の交点を求めることにする。
 式(3)の,中心座標(Xa, Ya),半径 |OA|の円,と,式(5)の直線の交点をM, Nとして,交点を表現したいのであるが,かなり複雑になるので,式(5)を簡潔化して次に表現する。

e=2(Xa-Xb), f=2(Ya-Yb), g=〔(|OA|^2 – |OB|^2) – (Xa^2 -Xb^2) – (Ya^2 – Yb^2)〕として式(5)は次のようにする。
eXo + fYo + g = 0⋯⋯式(6)

以上,ここまで,Jan. 25, ’24記。

 さて,Microsoft Excelの計算式の書式で,交点MNの座標値を次に示す。

Xm=(e*Df*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Xa⋯⋯式(7)
Ym=(f*D+e*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Ya⋯⋯式(8)

Xn=(e*D+f*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Xa⋯⋯式(9)
Yn=(f*De*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Ya⋯⋯式(10)

 なお,D=abs(e*Xa+f*Ya+g)⋯⋯式(11) とする。

 式(7)と式(9),そして,式(8)と式(10),の違いは第二項と第三項の間の符合に見える。

この計算プロセスを踏まえたMicrosoft Excelファイルは次節に用意する。

3.2 測量値などの準備

 ls_2放射星形区で,2カ所 ls2_11と_20の3Dマップ上の位置決めができている。ls2_11の座標値は変わらず,ls2_20だけを現場写真(図1)と対照して決定する必要がある。図4で特定している。光波測距儀による測量結果では,altitudeは1.171mだから,図4の右上に示した表では二行目が近い。まあ,これよりも海側下方であろうが,ほぼ対応しており,図4右上表の,Index 107222891を採用したい。
X: -96035.546875, Y: -625538.500000, Z: 1.243085

図4 ls2_20の想定位置でPicking listsを作成

 以上から,式(3), (4)の既知数値は次のようになる。なお,注意すべきなのは世界測地系の (x, y)=(northing, easting)と3Dマップの(X,Y)=(easting, northing)は入れ替わっているということである。学校教育で使用するデカルト座標系に替わっている。繰り返すが,ls2_11の座標値は3Dマップと一致している。両円の半径は光波測距儀の読値HD horizontal distanceに対応している。

ls2_11を中心とする円関係: (Xa, Ya)=(-96021.025, -625506.207), |OA|=42.830⋯⋯式(12)
ls2_20を中心とする円関係: (Xb, Yb)=(-96035.547, -625538.500), |OB|=65.508⋯⋯式(13)

式(12), 式(13)は式ではないが,便宜上,式番号を与えた。「3.1 計算法」の計算プロセスをMicrosoft Excelに埋め込もうと思う。 

3.3 Microsoft Excel計算シート

 計算シートを作成し,実行した。しかしながら,ls2_20の3Dマップ上の認定を失敗したようで,求める結果を得ることはできなかった。そのプロセスは次の節にしめす。Microsoft Excelの計算シートでは,まずは次のテキストボックスを配置した。

————————————————

解 説

 3Dマップ上で特定された2地点から光波測距儀の設置点を求める手法をここに示す。

現在,設置点としてきたものをO’点として,3Dマップに適切に載せうる設置点をO点(Xo, Yo)とする。この例では,ls_2放射星形区の2点のうちの1点について,現在の位置をA’点として,3Dマップに載せうる設置点をA点(Xa, Ya)とする。他の1点については,現在の位置をB’点として,3Dマップに載せうる設置点をB点(Xb, Yb)とする。

 結局のところ,A点(Xa, Ya)を中心とする半径|OA|の円と,B点(Xb, Yb)を中心とする半径|OB|の円の交点を求める問題になる。連立方程式は次のようになる。

(Xo- Xa)^2+(Yo- Ya)^2 = |OA|^2 ⋯⋯⋯⋯⋯式(1)
(Xo- Xb)^2+(Yo- Yb)^2 = |OB|^2 ⋯⋯⋯⋯⋯式(2)

 交点MNの座標値を次に示す。

e=2(Xa-Xb), f=2(Ya-Yb), g=〔(|OA|^2 – |OB|^2) – (Xa^2 -Xb^2) – (Ya^2 – Yb^2)〕⋯⋯式(3) として,

Xm=(e*Df*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Xa⋯⋯式(4)
Ym=(f*D+e*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Ya⋯⋯式(5)

Xn=(e*D+f*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Xa⋯⋯式(6)
Yn=(f*De*sqrt((e^2+f^2)*|OA|^2-D^2))/(e^2+f^2)+Ya⋯⋯式(7)

なお,D=abs(e*Xa+f*Ya+g)⋯⋯式(8) とする。————————————————

 図5に,この計算部分を掲載する。

図5 Microsoft Excelの計算シートの計算部分

 図5のダウンロードファイルは用意しているが,ls_2放射星形区内の2点の選定について,試行錯誤の状態のため,成功の後のファイルを掲載したい。

3.4 計算結果の検証1

 図5の計算結果は,Ⅳ〜Ⅶに見えるが,位置関係から見ると,(Xm, Ym)を採用できる。内積計算でも評価できるが省略する。当方のMicrosoft Excelファイル Okinoerabujima川口木庭_2018Jun14 参考ファイルJan2024.xlsxの沖永良部島測量ls No.2では,ls_2の光波設置点はGeographicaで得られたGPS表示を採用しているが,その数値を,上記(Xm, Ym)に替えた。このスプレッドシートは,northing X, easting Yに設定しているので,図5の座標定義と逆になっているので注意が必要である。

  ls_2の設置点を上記(Xm, Ym)に替えたのあるが,適切な結果が得られなかった。光波の視準線の磁北設定も変更する必要があった。もちろん,この視準角のゼロ位置を変更する必要性は当初から考えていたが,結局,3Dマップの真北設定でのXY座標値の土俵に入っているのだから,光波の測定値も磁北の偏角を考慮する部分は削除する必要があったのである。

 「図6 光波計測値の表に図5の計算結果を入力して磁気偏角をゼロに」を参照頂きたい。この計算過程を実施する前に,もともとのls_2の計算シート(沖永良部島測量 ls No. 2)をコピーして,この右隣に配置した。そして,最上行の「ls-site northing X」, 「ls-siti easting Y」カラムに,上記(Ym, Xm)を入力し,さらに,磁気偏角をゼロにするために,HRの角度を,補正済み真北からの角度,にそのまま,コピペした。その変更を赤字で示している。Loc. 11はls2_11に,Loc. 20はls2_20に該当している。もともとのls_2の計算シート(沖永良部島測量 ls No. 2)の値と比べると,ls2_11は,(x North, y East) = (0.19, 0.215),ls2_20は,(x North, y East) = (0.46, 1.901),となった。

図6 光波計測値の表に図5の計算結果を入力して磁気偏角をゼロに

 この結果から,次の表に示したcsvファイルを作成し,CloudCompareに読み込んだ。

site_name	X easting	Y northing	altitude	Red	Green	Blue
ls2_5	-96002.656 	-625419.660 	1.872 	0	0	255
ls2_6	-96018.265 	-625483.024 	0.466 	0	0	255
ls2_7	-96022.656 	-625483.303 	0.153 	0	0	255
ls2_8	-96024.513 	-625483.271 	-0.493 	0	0	255
ls2_9	-96032.904 	-625480.370 	-0.113 	0	0	255
ls2_10	-96015.073 	-625486.809 	2.183 	0	0	255
ls2_11	-96020.810 	-625506.010 	2.159 	0	0	255
ls2_12	-96040.427 	-625525.617 	-0.510 	0	0	255
ls2_13	-96039.064 	-625525.209 	0.437 	0	0	255
ls2_14	-96032.476 	-625524.511 	0.432 	0	0	255
ls2_15	-96030.971 	-625523.683 	0.917 	0	0	255
ls2_16	-96029.524 	-625529.911 	2.325 	0	0	255
ls2_17	-96039.251 	-625539.972 	-0.491 	0	0	255
ls2_18	-96037.813 	-625540.109 	0.665 	0	0	255
ls2_19	-96034.701 	-625539.333 	0.733 	0	0	255
ls2_20	-96033.646 	-625538.043 	1.171 	0	0	255
ls2_21	-96033.215 	-625538.873 	2.464 	0	0	255
ls2_22	-96032.121 	-625554.452 	0.614 	0	0	255
ls2_23	-96026.456 	-625577.528 	0.404 	0	0	255

 図7を見ると,ls2_20などの位置はより左手(つまり海側)に移動しているが,例えばls2_18は,ls2_17付近に位置する必要がある。図8は新プロットだけを表示したものであるが,垂直関係に問題はないようであるが,図7のls2_19は本来,ls2_18が載っている緑藻の面の陸端に配置すべきものであり,光波測点が3Dマップとうまくフィットしていないことが見えるのである。

図7 新旧光波測点の比較 旧よりも新の方がより現場写真に近い

図8 新測点だけ

 結果からすると,ls2_20の3Dマップでの特定に誤りがあったようである。現場写真を見直して,ls2_19の方が特定しやすいように思われるので,これまでls2_20で実行したことを,ls2_19で試みることにした。

以上,14:36,Jan. 28, 2024記。

3.5 計算結果の検証2

 ls2_19の図10での想定位置位置を,図9の現場写真から決定した。図10の表の2の位置が適当と思われる。図10の鮮明度は低く,2階調にしている。このデータは次のようである。
X easting =-96036.617188, Y northing = -625539.000000, Z = 0.711838

図9 ls2_19の現場写真

図10 図9から想定できるls2_19をlist-picking

 ls2_11をCloudCompareで観察する限り,大きなズレはみにくいが,なにせ,drone撮影の影になっており,特定するのはかなり難しい。ls_19を採用してもうまく行かない場合は,ls2_11の代わりをみつける必要性がある。大変だなあ。

3.6 計算結果の検証3

 さて,これまで触らなかったls2_11に代わって,ls2_6(Z=0.466 m a.s.l.)の現場写真を観察して,3Dマップ上の座標値を求めた。図12右上のpicked points list表には2点が見えるが,下方の方を採用した。
 なお,図12のls2_6の海抜高度は高いように見えるが,Cloudを透過して表示されているからである。list picking作業で,水平位置が海側に移動することになる。

X easting = -96021.109375, Y northing = -625482.125000, Z = 0.445537である。

図11 ls2_6の現場写真

図12 CloudCompare上

 このページは一昨日以来,触らなかった。ls2_6付近の表示の差を図13,14に示している。図2では北縁の測点はls2_11までしか見えていない。図13, 14の場所は図2よりも北方にある。図11のプリズム三脚の位置がls2_6にあたるが,図13では3Dマップ上のls2_6ポイントは,図11よりも更に内陸側に位置してる。図14を見ると,ls2_6は,いい按配に見える。それゆえ,両円の中心として,半径|OA|を作るA点をls2_11からls2_6に替えることで,大幅に改善されたことがわかる。

図13 ls2_19で補正した場合

図14 ls2_19だけでなく,ls2_6も変更した場合

 このように,A点=ls2_6,B点=ls2_19,とすることで,図14のように,かなり改善を見た。では,ls2_19付近はどうなったのか。図7と比較すると,図15ではより改善され,光波の測点はより海側に移動した。

図15 ls2_20や_19はより海側に移った

 A点をls2_11からls2_6に替えることで,ls2_6付近も,ls2_19付近も改善されたのである。ぼくの性格から,ls2_20とls2_6の組み合わせはどうかとやってみたが,劣化した。今のところ,ls2_6, ls2_19の組み合わせがいいようだ。

図16 改めてls2_6付近を

 図15で気になるのは海側にずれ込み過ぎていること,図16で気になるのは,ls2_6. _7, _8が左方にずれ込んでいることである。

 ls2_19付近の3Dマップと現場写真の対照はたぶん,これ以上の改善は難しい。ls2_6付近はどうか。これは悩んだ上での決定であったが,droneで撮影した空中写真の広域の方で見えるのではないか,と考えている。明日は,空中写真を判読して3Dマップと対照して,ls2_6の位置決めをしたいと思っている。

以上,23:18,Jan. 29, 2024記。

3.7 計算結果の検証4

 現場写真から,ls2_6をほぼ特定した。この海側(写真の上方)に波のスプラッシュが見えるが,これは図16のls2_7(Sw-spP)に対応するだろう。ls2_6は,Lw-spPだ。宇和美崎狭域3Dマップで,この場所を特定する必要がある。宇和美崎広域写真で特定したls2_6の位置だ。spratidal platformの陸端だ。

図17 宇和美崎広域の 358part(左), 359 (右)の実体視 ls2_6の位置は〇

 図18には,宇和美崎狭域3Dマップ上で新たに特定したls2_6は,Point #1にあたる。

X= -96021.5703125,Y= -625484.5625, Z= 0.445600688457

図18 宇和美崎狭域3Dマップ上にls2_6を特定したPint #1

 やっと得たls2_6が新たなA点となる。B点(ls2_19)はすでに確定している。

 両円の交点をもとめ,それをls_2の設置点としてls2_5〜_23の(X, Y)を再計算してcsvファイル ls_2star_ls2_new6_19.csvを用意して,CloudCompareに読み込んだ。図19では,ls2_6付近のもので,右寄りのものが新たな結果で,フィールド写真とよく合っている。

図19 ls2_6付近 右側が改善された測点

 図20では,ls2_19付近のプロットで,過去で一番改善されたのであるが,いまだ,多少,海側にずれ込んでいる。

図20 ls2_19付近 右側が改善された測点

 これで,ls_2放射星域については終了しても良いが,ls2_19の位置決めを空中写真判読でも確認したいと思う。

以上,Jan. 30, 2024記。

3.8 計算結果の検証5

 ls2_19を写真(図9)を中心に,宇和美崎写真広域355, 356の実体視を加えて,検証した結果,あらたな位置を特定できた。図21には,新たにls2_19が特定されている(図中のPoint #0)。これまでの最新のプロットと対照すると,ls2_19などがより内陸側に移動できそうである。X= -96036.671875, Y= -625539.3125, Z= 0.754368662834となる。光波でのZ値は0.733m a.s.l.である。

図21 新たにls2_19を特定した(図中のPoint #0) これとこれまでの最新のプロットと対照している

1 ls2_new6とともに,_new19が更新され,両円の交点を新たに得ることができる。
2 得られた交点をls_2光波設置点としてls2_5〜_23の(X, Y, Z)を求める。
3 Microsoft Excelシートからcsvファイルを作成。
4 CloudCompareで読み込んで宇和美崎狭域3Dマップ上にプロットする。

 この結果を次に。図22には,平行に配置されている同名のものの右下が新たなls2_19によるもので,本来,ls2_19であるべき場所のPoint #0には到達できていない。何故だろう。

図22

  図23の結果は,図21の結果と極めて類似している。平行にズレている。さて?

図23 平行に配置されている同名のものの右下が新たなls2_19による結果

以上,Jan. 31, 2024記。

4. 計算法を確認したい

4.1 問題の所在

 これまで述べてきたように,図22のPoint#0は,現場写真,空中写真を見つめて,3Dマップにls_19の位置と決定した場所に,CloudCompare上でPicking points list機能を使って,プロットした地点である。ベクトルOBの関連では,Oがls_2設置点で,BがこのPoint#0にあたる。両円の交点がOだから,|OB|=66.509mにあたる。光波測量の計算表に示しているように,|OB|は,HD horizontal distanceに当たっている。そういうように設定した訳だから,交点とls2_19の距離は,|OB|=66.509mの筈であり,宇和美崎狭域3Dマップ上のPoint#0と光波計算表のls2_19は一致する筈である。なのに,なぜ一致しないのだろうか。

 改めて,宇和美崎狭域3Dマップ上に描いた,ls_2光波設置点(図24のls_2setting)と,ls2_5〜23の分布を見たい。

図24 宇和美崎狭域3Dマップ上に描いた,ls_2光波設置点(図24のls_2setting)と,ls2_5〜23の分布

以上,Jan. 31, 2024記。

 両円の交点を求めるこれまでのプロセスを,図25のように図化した。何故に今採用している計算プロセスを経て得られたls2_19が現場写真から特定したPoint #1と一致しないのか,それを特定するためである。

図25 両円の交点を求める図を作成して

 これからのプロセスを確認すべく,流れを次に示す。

確認1 図25での両円の半径が正しく光波測距儀のHDと一致しているか。
 図25を評価する際に二つの面から調べる必要がある。点A, Bの座標値は,現場写真と空中写真から特定して3Dマップにpicking pointsの座標値にあたっている。これをここで,再掲したい。そしてMicrosoft Excelで計算して|OA|,|OB|が光波のHDに一致することを確認する。
 両円の交点を求めてそれをls_2設置点の座標値として,光波の計算プロセスで得られたls2_6, _19の座標値を使った場合,|OA|,|OB|が光波のHDに一致するか。まあ,一致しないだろう。しかし,一応,この結果を確認したいと思う。

確認2 1の検証結果は,光波の計算プロセスで歪みが生じたということになるだろう。では原因は何か。

確認3 まだやってないからわからないけど,まずは3Dマップが真北を捉えていない。gcpを使って適当に配置している。そして光波で得たネットワークを構成する放射星形区も真北がわからない。両方に原因がある。しかし,3Dマップに光波放射星形区を載せようということだから,3Dマップの擬似真北に光波の放射星形区も合わせないといけない。じゃあ,どうするか。
 図25のベクトルOAの時計回りの角度は,点ls_2settingの座標値を(Xo, Yo)として,現場写真と空中写真から特定された点ls2_19photoの座標値を(Xb, Yb)として,直角三角形の,角 [(Xb, Yo)(この頂点が直角)-(Xo, Yo)-(Xb, Yb)]を余弦定理から求めれば良い。その角度と光波測量時のHR(時計回りの回転角)を対照すれば,西偏用の角度補正のセルを変更すれば,ls_2設置点からの測点の(X, Y)座標値を得ることができるはずだ。

 楽しみだ。今日はタニハで疲れた。明日にした方がいいな。

以上,0:20,Feb. 2, 2024記。

 次にMicrosoft Excelで検証した結果を示す。
————————————————
1. OAとOBを半径とする円の中心と1交点
X Y Z
A相当:ls2_6photo -96021.570 -625484.563 0.446
B相当:ls2_19photo -96036.672 -625539.313 0.754
光波設置点:ls_2setting -96052.835 -625474.797 -0.235

2. |OA|と|OB|の距離と光波HDとの関係 計算結果 光波HD 評価
sqrt((Xa-Xo)^2+(Ya-Yo)^2) 32.754 32.754 当然,一致
sqrt((Xb-Xo)^2+(Yb-Yo)^2) 66.509 66.509 当然,一致
評価: 両円から交点を求めるプロセスに誤りはない。

3. 光波計算過程を得て得られたls2_new6, ls2_new19
X Y Z
ls2_new6 -96021.186 -625483.235 0.466
ls2_new19 -96037.622 -625539.543 0.733
光波設置点:ls_2setting -96052.835 -625474.797 -0.235

4. |OA’|と|OB’|の距離と光波HDとの関係 計算結果 光波HD 評価
sqrt((Xa-Xo)^2+(Ya-Yo)^2) 32.754 32.754 そうか,やっぱり一致
sqrt((Xb-Xo)^2+(Yb-Yo)^2) 66.509 66.509 そうか,やっぱり一致
————————————————

 両円の交点の計算過程でも,これを使った光波放射星形計算過程でも,|OA|, |OB|いずれも変質していないことがわかった。すぐ上の確認1〜確認3での思考プロセスに問題がないことがわかった。

4.2 3DマップのXY座標系に光波放射星形区のX’Y’座標系を合わせる

 用語が定まらないが,このタイトルで考えたいと思う。

 現場写真とdroneで撮影した空中写真を読み込んで,ls2_6とls2_19を中心とする円の交点を使ってきたのだが,点群の分解能と現場写真の質の点から,ls2_19の信頼性が高く,これを優先して計算過程を進めることになる。

 図25の3Dマップの点群を非表示にして,Yo軸からls2_19軸までの時計回りの角度を求めたい。黄色で表示している△BOIに注目して,余弦定理で計算した。

図26 3Dマップ上でのls2_19軸の時計回りの角度を求める

————————————————
5. 余弦定理で直角三角形BOIの∠Oを求める
辺b=辺OI=|Xb-Xo| 16.163
辺i=辺BO=|OB| 66.509
辺o=辺IB=|Yo-Yb| 64.515
cosO: (b^2+i^2-o^2)/2*b*i 0.243
arccos(cosO) 1.325 radian
arccos(cosO)+p()/2 (Y軸からの時計回り) 2.896 : 3Dマップ上での時計回りのOB角 radian
ls2_19のHRは166度46分40秒を度表示で 166.778
radianに変換 2.911 :光波計算表でのHR radian
光波計算表でのHR補正値radian -0.015 radians
radiansからdegreesに変換 -0.841 degrees
degree + minute/60 + second/3600 -0°50′27.6″ :この値を磁気偏角にして実行する。
————————————————
 上記の計算過程で,光波計算表で,実行した。
————————————————
6. 座標回転を考慮した光波計算過程で得られた地点Bと地点Aの(X, Y)座標の結果
現場写真3Dmap 光波座標回転後 差
OB: ls2_new19 Xb -96036.672 -96036.674 -0.002
OB: ls2_new19 Yb -625539.313 -625539.313 -0.000
結果評価: OB軸で座標回転したから当然一致する
OA: ls2_new6 Xa -96021.570 -96021.066 0.505
OA: ls2_new6 Ya -625484.563 -625482.769 1.793
結果評価: 次の項で述べる。
————————————————

この結果を使って,ls2_5〜_23の点情報のテキストファイルを作成して,CloudCompareに読み込んだ。その結果を次に示す。図27はls2_19付近,図28はls2_6付近である。解像度が低くて見にくいが,いずれも海岸地形の断面図作成用に設置した測点に完全に一致している。

 このウェブページ作成の目的を完璧に達成できているのである。

図27 ls2_19付近

図28 ls2_6付近’

 ls2_19はベクトルOBを重視して,ls2_6の方も見事に達成できたのである。

 で,やっぱり,気になるのは,ls2_6つまり,OAの方である。計算してみた。
————————————————
7. OA軸についても座標回転の偏差を求めるかどうか
 一つの放射星形区で,現場写真から3Dマップとの関係での2地点を特定して,ベクトルOA, OB を中心とする半径HDa, HDbの円の交点を求めてそれを光波測距儀の設置点とするのであるが,地点A, Bを 特定して初めてこの作業が可能である。ここで示している例2.ではOBから光波測量結果をまとめた のであるが, OAについてもこの作業をする意味はない。地点A,Bでは,Aの方が信頼性が高く,Aの結果から得た偏角を採用して光波測量値を処理するしかない。
 光波のOA, OBの視軸がなす角度は,HRb’ – HRa’から得ることができる。radianでは,2.896-1.817=1.079である。一応,ここでのHRb – HRa を求めてみよう。直角三角形 BOAの∠BOAを次に求めてみる。

地点A, Bの座標値を再掲する。
X Y
A相当:ls2_6photo -96021.570 -625484.563
B相当:ls2_19photo -96036.672 -625539.313

辺b=辺OA 32.754
辺a=|OB| 66.509
辺o=辺AB 56.795
十進角(º)
cosO: (b^2+a^2-o^2)/2*b*a 0.521
arccos(cosO) = (HRb-Hra) 1.023 58.590
挟角の差 (HRb-Hra) – (HRb’ -Hra’) radian -0.056 -3.233
なお,(HRb’ – Hra’)は, 1.079 61.822
OA軸のY軸からの時計回り角度 radian 1.873 107.339
ls2_6の時計回り角度 radian 1.817 104.106

光波計算表でのHR補正値radian 0.056 3.233
:OAを採用すると補正角
————————————————

 以上で,一応納得の行く結果となった。
以上,Feb. 3, 2024記。

5. ls_2放射星形区の結果を使って次はls_1区

 これまでの努力で,ls_2は完成した。これに隣接する,ls_1とls_3にアタックすることになるが,まずは_1を。ただし,今のところ作成できた3Dマップは宇和美崎狭域のもので,ls_1の範囲の全域がカバーされていない。広域の3Dマップは未完成だ。
 放射星形区の間では,必ずしも2点での前視後視を実施していないが,ls_1とls_2の間では実行しており,Locs. 5, 6がそれに該当する。3Dマップとの対応関係はls_2で成功したので,ls_1の光波設置点は,ls2_5, 6の(X,Y)値を使って実行できるはずである。

5.1 新たなls_1の座標値を求める

 前掲の「3.ls_2の3Dマップ上の座標値を求める」,のプロセスをまずは実行することになる。(X, Y)値はls_2で最終的に得られたLocs 5, 6の値をls_1のそれとして代用することになる。HDについては,ls_1からのHDを使用する 。
easting   northing HD (ls_1)
ls2_5(ls1_5) 点A -96006.388 -625419.18372.459 (過ってls_2のもの),72.465 (正しくls_1のそれ)
ls2_6 (ls1_6) 点B -96021.066 -625482.76932.754(過ってls_2のもの),33.157(正しくls_1のそれ)

「3Dマップ上の2点から光波測距儀設置点を求める.xlsx」には問題があった。新しく2円の交点を内積で求める(2円交点EnToEn_revised.xlsx)を作成して両円の交点は2点得られた。Google Earthにプロットした旧来のls1_5, _6,GPSのls_1の配置からみて,Xm=-96053.17945, Ym= -625474.51596をls_1設置点と決定することができた。

5.2 新たなls_1の光軸回転角を求める

 「Okinoerabujima川口木庭_2018Jun14参考ファイルJan2024.xlsx」に,上記で得られた新たなls_1設置点の座標値を,入力した。すでにこのページで議論してきたように,さらに,光波測量の放射星形区の測量結果の座標回転が必要になる。「4.4 3DマップのXY座標系に光波放射星形区のX’Y’座標系を合わせる」,で述べた作業を実施する。△BOAを点O(新ls_1)の回りに回転するが,OA軸,OB軸の何れを使用するかで多少異なるが,現場写真を見ると,未完の宇和美崎広域3Dマップがたとえ後に完成しても,点A (ls1_5)を3Dマップにプロットするのは難しい。そこで,OB軸の回転角度を求めることにした。結果として,△の頂点と辺の名称はほぼそのまま使用できる。使用Microsoft Excelファイルは,「3Dmap上で光波星形区をフィット両円交点と座標回転.xlsx」(前掲のファイル名をより一般化してみた)にあたる。その計算シートを図29に示す。

 第29行(赤線の1)にはには余弦定理を使って求めた十進度104.413が見える。第31行(赤線1’)では,それをDMS表示に変換している。そして,その値をls_1星形区のLoc. 6 (ls1_6) のHR(DMS)で差し引いた結果が,赤線2に示す-02:53:53.6である。DMSの引き算の結果は➕️である必要があり,それを,第21〜27行右手に見えるテキスト枠に示している。第34行はその演算行にあたる。

図29 ls_1両円交点と光波軸の回転

 この「偏差入力値」を使って,ls_1星形区のLocs. 1〜6の再計算を実施する。

 この後の計算の流れは,Microsoft Excelで度分秒表示角の差を求める ,「3. 光波測量結果に」に示している。その図5を,ここでは図30として再掲した。

図30 「図5-hh:mm:ss.0のDMSセル分割」

 これまで述べてきた手法で,ls_3以降の星形区の光波結果を再計算したいと思っている。

 次の図31はCloudCompareで表示したところである。ls1_1〜_6とls2_5〜_13である。

以上,Feb. 19, 2024記。




犀川開析谷の水災環境 the potentiality of flood damage in the dissected valley bottom lowland of the Sai-kawa River

1. 自然環境を無視した住宅開発

 日本では,とくに,明治維新以降,西洋の高い土木技術の導入と人口圧や都市化の要請のなか,本来住めない場所でも,住宅開発が積極的に進められてきた。その一つが,犀川開析谷である。金沢市のWebサイトを参照したが,コンテンツの質は他の行政区同様であった。ハザードマップの質は日本のどこでも低い。

2. ハザードマップから

 さて,金沢市のハザードマップを参照して,過去の浸水実績を見たのが次の図1である。図1右手はご自宅の場所をGoogle mapから特定して図1のハザードマップでの位置を同定している。浸水域がモザイク状になっているのは,報告書を参照していないが,排水溝のトラブルによるものであろう。つまり,破堤によるものではない。犀川の水路幅が小さくて上流の域の排水を運ぶような機能が期待されていない。犀川は,川の名前が付いているが,自然の川ではなく,単なる水路になっている。では,安心していいのだろうか。

図1 浸水実績区域の分布図

 図2は,ご自宅の想定水害水位である。0.5m以上〜3.0m未満,河岸浸食となっている。さて,図1とどう違うのか。図1は比較的最近の実績である。リスクは図2をみた方がいいだろう。破堤の可能性は当然ある。その破堤をした際のリスクを図2は示している。すぐ家のそばの河岸が破壊されるということではなく,上流,場合によっては,下流であっても,破堤した場合,0.5以上,3.0m未満の浸水があるということだ。

 それゆえ,集中豪雨などがあって,非難注意報が出た場合,逃げた方が良いことになる。津波について金沢市のハザードマップは,かなり海岸付近に限定されているような図しかないが,10mぐらいの津波が来たら,この水路を遡って,超波することはありうる。ただ,10mの津波が生じるかどうか。5mぐらいの津波だったら,ここまで上がってくるのは難しいかも知れない。

 なお,逃げる場所だけど,行政(業者)の指定先は,全く間違っている。信用してはいけない。今回の能登半島の津波の際の指定避難所に到達できず,津波に呑まれている。東日本大震災の場合も避難場所は,海岸線のそばの海岸砂丘上であった。東北大の津波研究者の推薦場所だった。責任取らずにのうのうとTVにも今なお出ている。ハザードマップや行政の避難場所を信じてはいけない。

図2 浸水想定図

3. 破堤の場合のお宅に係わって想定されるルート

図3 可能性の高い破堤

図3にぼくが描いてみた破堤ルートである。集中豪雨があった場合,青色の線の上流側の部分から溢水して破堤する。その中心軸は,お宅を通っている。もちろん広範に洪水は生じるが,破堤場所として可能性が高い地点を考えて見た。

以上,Jan. 16, 2024記。途中です。




CloudCompareでCloud layersを handling cloud layers in CloudCompare

1. このページ作成の意味

 もともとの関心はGISソフトでの3Dクラウドとgcpの重ね表示の可能性を探ることであった。ここに至る経緯は次のページに示している。

 このページの作成過程で,gisソフトでは実現しないことがわかった。MeshLabに期待したが,これも重い3Dファイルを受け付けない点で問題があった。結局,CloudCompareに頼らざるを得なかった。CloudCompareではMeshLabでどうにもならないものでも何とか作業ができるのである。

 で,問題はレイヤー構造で使用出来るかであるが,Danielの情報で,使えることがわかったのである。ただ,図1のように,Plugins: Cloud layers,は灰色になっている。

図1 (図7 使用しているCloudCompareのPlugins: Cloud layers は灰色)再掲

 Help > About Plugins,で開いた画面が図2である。Cloud layersを選択するが,これ以上進めない。Referencesとして,サードパーティのWiggins Tech website が見えるが,ここをクリックするとこの注記が消えてしまう。

図2 Cloud layersは表示されているが

2. Plugins: Cloud layersのインストールに向けて

図3 Webサイト 3Ds-scan.de にプラグインが

 図3を参照して欲しい。なにか不可解なのだが,CloudCompareのサイトではなくて,3ds-scan.de で公開されている。ぼくがインストールしているのは,図1の v 2.12.4 (Kyiv)なので,2.12.0はより古いバージョンだ。

 キーウで現在活動する企業(サードパーティ)のようだ。

 このサードパーティのリンクをCloudCompareを使用しているWindowsで立ち上げた。ドイツ語が表示される。言語を換える選択肢があり,日本語に替えた。その画面が図3のものである。このサイトは2016年11月から始まっている。図3の「Cloud Layers」プラグインをクリックするとmp4のムーヴィーが始まるがCCもなく,音声も出ない。

 図3の表示のすぐ下には,ダウンロードリンクがある。CloudCompareのミラーサイトの認証を受けているようだが,ヴァージョンが2.12.0のままだ。更新されていない。CloudCompareはフリーだけどサポートが要るというメッセージが色つきになっている。

 次の情報がある。ただし,当方の Editメニューでは,’Edit > Cloud > Create single point cloud‘しか有効でない。Edit > Scalar fields についても下に示された選択肢はすべて灰色になっているのである。

Here are the new menu entries:

  • Edit > Cloud > Create single point cloud

    • to create a cloud with a single point (set by the user)

  • Edit > Cloud > Paste from clipboard‘ (shortcut: CTRL+P)

    • to create a cloud from ASCII/text data stored in the clipboard

  • Edit > Scalar fields > Split clouds (integer values)

    • will split the cloud in multiple sub-clouds based on the (integer) values of its active scalar field.
    • to be used with a classification SF typically.

  • Edit > Scalar fields > Add classification SF

    • shortcut to the ‘Edit > Scalar fields > Add constant SF‘ method.
    • creates a ‘Classification’ SF with a constant (integer) value.

 このWebサイトに入ると,ぼくはすごく怪しい印象を持った。どこにもプラグインのダウンロードのリンクはなく,ただ,コンタクトせよという。しかも種々のページを見ていると,bluetoothに繋げようとする。

 Danielさんが,フォーラムで,
Maybe you can use the new ‘qCloudLayers’ plugin?
Another way, maybe more simpler, is to use the ‘Edit > Scalar fields > Split cloud (integer values)’ method. This will create as many (sub) clouds as classes present in your cloud.

という回答をしている理由もそこにありそうな気がしている。

3. 最新バージョンをインストールしてみる

 信頼できるGitHubギットハブで確認した。Version 2.12.0には,確かに,上記の魅力的なプラグインが新たに装備されている。では,何故,最新のステーブル版では灰色になってしまったのか。これがまずは最初のポイントだ。ver. 2.12.0のコードネームKyivを獲得したのもキーウの貢献ゆえだろう。

 V. 2.12.1〜.4まで,特に当該プラグインがinactiveになったのかは,記されていない。

 ではと,最新のアルファ版をダウンロードしてみるのがいいのだろう。GitHubよりももちろん,本家からがいいだろう。本家をみると,より進んだβ版になっている。
CloudCompare 3D point cloud and mesh processing software Open Source Project

 インストーラーは265MB以上ありダウンロード継続の表示が見えないので非常に取り扱い難い。圧縮版は64MBで軽く,7zで解凍して,ドラッグアンドドロップで,Program Filesに投げ込んだだけだ。ショートカットを作って,デスクトップに配置した。macの方は不思議なのだがダウンロードが始まったのは見えたが,時間をかけて見ると,消失している。2回繰り返したが同様の現象が見られた。ファインダーで検索しても見当たらない。図4にはGitHubのダウンロードページであるが,このオレンジでマークしたtar.gzファイルをダウンロードした方が早かったであろう。これは7zで解凍する。Windows64用である。

図4 GitHubからダウンロード

図5 CC v.2.13βにみえるPluginsリスト

図6 Cloud layersを選択すると,Cloud has no scaler fieldsのエラーが

図7 Edit > Scaler fieldsのsubmenu

 図5〜7は,新たにダウンロードしたCloudCompare v. 2.13βの表情である。CloudCompareを選択して,v. 2.12.4と比べてみた。確実に改善されている。灰色だったものがactiveになっている。図5のように,開いたPluginsの下方にCloud layersが見える。これを選択した結果が,図6で,Cloud has no scaler fieldsのエラーが現れる。図7には,Editのメニューが見えている。このScaler fieldsを選択した際のリストが見えているが,このうち,Add constant SFなどが関係しているようである。このSFは,Scaler fieldsの略号と思われ,この定義をした上で,Cloud layersを見ることができるのだろう。

図8 Cloud layersの一例

 図8を見れば,Cloud layerが何ものか,わかるであろう。これは実際にやってみないとわからない。海岸地形は緑藻や紅藻のなどの分布が地形的な位置づけを反映しているので,おおよその区分を実施する上で参考にはなるだろう。

 とはいえ,ぼくはこのようなCloud layerを求めていた訳では無いので,残念っという思いの方が強い。

 このCloudCompareのCloud layersはこの程度のものなのか。とにかく,Danielの意見に戻らないといけない。

4. gcpクラウドを作成すべく既存のツールに期待したのだが

 次の,Edit > Cloudのうち,single point cloudしか,アクティブになっていない。もしクリップボードにデータがあれば,アクティブになる可能性がある。

  • Edit > Cloud > Create single point cloud

    • to create a cloud with a single point (set by the user)

  • Edit > Cloud > Paste from clipboard‘ (shortcut: CTRL+P)

    • to create a cloud from ASCII/text data stored in the clipboard

2章の末尾のDanielさんのアドバイスを次に再掲する。

Maybe you can use the new ‘qCloudLayers’ plugin?
Another way, maybe more simpler, is to use the ‘Edit > Scalar fields > Split cloud (integer values)‘ method. This will create as many (sub) clouds as classes present in your cloud.

 メニューを見ると,Edit > Scalar fields > Split cloud (integer values),はグレイになっている。これも何らかの前処理があってのことである。さて。

以上,19:32,2024/01/07記。

 ちょっと,考えて見て,gcpのクラウドをエクセルやエディターで作成して,CloudCompareに取り込んで,メーンのdroneから得られたcloudの上のレイヤーとして表示するだけのことじゃないか,と考えついた。以前から頭にあったことが出てきたという感じである。

 で,その構造であるが,easting, northing, altitude,なのだけど,できれば,gcpの名称も欲しい。その観点で次に。

5. gcp cloudの書式は?

 どうしてこんな簡単なことをすぐにしなかったのか。CloudCompareでは幾つものメッシュやクラウドを取り込んで,DB Treeで,上下に移動してきたので,CloudCompareを使い始めてから,実はクラウドのレイヤーを操作してきたのである。

5.1 まずはEdit > Cloud > Create single point cloud

 結局,数点のクラウドをどのように作るかである。Edit > Cloud > Create single point cloud,って,何だろうか,ってわけで実際に使ってみた。x(easting: 地球座標系とは異なることに注意), y(northing), z値(altitude)を入力した(図9)。この値は,苦肉の策の,座標回転後のものである。

 この本題とは異なるが,2020年の調査結果をまとめたフォルダー毎,みつからない。別氏から送られてきたこれまでのdroneデータを全部削除する過程で,このフォルダーも殺したらしい。macのトラブルがあって,それに関するページを次のようにDec. 14に作成している。その翌日にTimeMachineでバックアップしたのだが,ここに殺したフォルダが残っていて,コピペした。幸い,座標回転までのエクセルデータがあった。このあたりから,別氏の再計算を待つ形で,キンドル本発行に注力したらしい。助かった。インフルエンザで多少熱があて鼻水だらだらながらの削除作業だったので,失策をしたらしい。怖いなあ−。

図9 Edit > Cloud > Create single point cloudの入力

 図10は,座標値の簡略化の誘いの画面であるが,これに乗ってはならない。個々の座標値入力の際には必ず出てくる。drone測量の結果であるPLYファイルをCloudCompareに取り込む際にもNoを選択しているので,地球座標系の数値(GIS用にXYは例によって逆転)のままである。同じ座標系でないとgcpは宇和美崎オブジェクトに載らない。

図10 座標値の簡略化の誘いに乗らないように

 図10左ペーン最上部のDB Treeには,a_rot, f_rot,が見える。cloud☁︎(図中は白色だが)アイコンがついている。階層は🗂️アイコンがついたdroneファイルと同じだ。これでは,全部選んだとしても,CloudCompareのbinファイルとしては保存できない。

 ここで取り上げているCreate single point cloudコマンドで実行すると,DB Treeには自動で#001, #002,などとナンバーリングされる。この名称をタップすると,名称変更が可能なように入力できるようになって,a_rot,f_rot,とした訳である。メーンペーンでは白い小さな丸が現れて位置を確認できるのである。左ペーン中央のPropertiesのWindowでは,Name: a_rot, Visable: ☑️,Show Name:◻︎,が見える。Show Nameに☑️すると,メーンペーンに「a_rot」が現れる。

 これはいい!,と他の点も作成してゆくのである。5点まで登録して,「待てよ」となった。とにかく,一度,保存してみようと。図11左ペーンに見えるように,差し当たり作成した5点をまとめて保存することにした。File > save,だ。メーンペーンに,Save Ascii file,のパネルが現れた。しめしめ。パネルの意味はなんとなくわかる。

図11 現れたパネル: Save Ascii file

 図12 Save Ascii fileパネルに適宜入力したけど。coordinate precisionを8から12に(xy座標の期待できる精度から)。Scalarは変更せず6で(a_rotなどの名称ではないかと考えては英数字6文字で十分だと)。テキストとしては扱いやすいcsvで(それでseparatorをcommaに),白いアイコンも欲しいのでColorsもそのまま。

図12 Save Ascii fileパネルに適宜入力したけど

 さて,もともとのPLYがあったフォルダーには,5点が別々のcsvファイルとして保存されていた。名称は,gcp_000000.csv〜_000004.csv,である。まあ,いわばブラウザーでまとめて選んだだけで一つのファイルになるという想定が問題ではあった。

図13 ブラウザーでの5点ファイル

 なんと,ぼくはこのブラウザーで,gcp_000000.csvをgcp_a_rot.csvなどと変えちゃったりしたのである。これらのファイルを開くと自動でエクセルが開くが,名称やアイコンの色情報もなく,ただただ,x,y,zの値のみだ。

図14 ブラウザーで5点のファイルを変更

 CloudCompareで,図15のように,この5点をまとめて,読み込んだ。

図15 この5点をまとめて読み込み

 図16のように,パネル「Open Ascii file」が個々の点について現れるのだ。

図16 Open Ascii file パネル

 図17のように,5点がそれぞれ,フォルダーを持っていて,その下位に配置されている。一点ごとに。これは不合理だ。ポイント位置が見えない。そこで,Show nameに☑️すると,ポイントの名称が所定の場所に現れるのでおよそは位置がわかる。これではどうしようもない。PropertiesのBox dimensionsでは,x,y,zともゼロになっている。ここの「10」程度の数値が必要だろう。

図17 大変なのだ

 以上の経験から次の方針: ポイント位置を保持するようなポイントの定義が必要だ。点をまとめて一つのフォルダーにまとめる必要がある。csvファイルに,ポイント名,座標情報,アイコン定義,が必要だ。

5.2 必要な属性を備えたgcpファイルを作成するには

 Create a new layer of just ground points (already classified) での質問回答は,Danielが開発者故のものであった。フランスの質問者は,まさにぼくのこのテーマなのである。ベーシックな質問であった。

 ASCII export variables – what are they? この情報に先に会っていたら,このページのここまでの試行錯誤は不要であったろう。図12のSave Ascii Fileパネルのパラメータに対するDanielの回答が次に。

Re: ASCII export variables – what are they?

Post by daniel » Sat Nov 05, 2022 3:03 pm

– coordinates precision: number of digits for coordinate values (X, Y, Z)
– scalar precision: number of digits for scalar field values (distance, intensity, etc.)
– separator: separator character between the above values (use a comma to create a CSV file for instance)
– order: order of the values (columns)
– header option: to add a description header (starts with ‘//’ and then you’ll find the name of each column)
– number of points: simply the total number of points on the first line (that’s expected if you want to create a valid ‘PTS’ file for some other tools)
– colors:
* use float values (between 0 and 1) instead of integer (between 0 and 255) for R, G and B components
* save the alpha (transparency) component as well

しかし,これではまだわからない。

 種々調べたけど,わからない。結局,自分でやるしかない。Danielさんに一度聞いてみたいものではあるが。xyz座標値は必須だが,マークbodyにxyz,それに,その色だ。そんなのはDanielさんの上のリストにはない。とにかく,エクセルでcsvファイルを作って,CloudCompareに読み込ませて,そこで指定してゆくのだろうなあ。

 図18のように,マックのエクセルでcsvファイルを作ってウィンドウズの所定のフォルダーに入れた。それをCloudCompareで開いたら,図19のように,念願の表Open Ascii Fileが出た。これが欲しかった。やってみないとわからない。図19の葡萄茶色の行はぼくが作成した列タイトルだ。適当に,Placemark, x-easting, y-northing, z-altitude, body-x, -y, -z, R, G, Bとしたのが,このテキストから適当に解釈して,それぞれ, この葡萄茶色の行のすぐ上の灰色の行のように,Label, cood.X, Y, Z, SF Scalar 3個,Red (0-255), Green (0-255), B (0-255), と配列している。

図18 macのエクセルで作ったcsvファイルをウィンドウズに配置

図19 Open Ascii Fileが出たあ

図20 タイトル行を削除

図21 図20の右方への続き

 自ら設定した列タイトルに合わない場合,図22のように,変更することができる。

図22 用意された見出しを出してみた

Nx, Ny, Nzは,Per-point normal vectorsでここでは関係ない。SF-Scalarは,Per-point scalar values,であるが,ぼくはbox-xなどとしてしまったものでこれに対応している。ぼくが設定したRGBは,A unique color for the whole entityに対応しているようだ。Quaternion WXYZは複素数のようだ。とにかく,わからない。

 とにかく,図22の右下の,apply allのボタンを押した。図23の右上には正方形が表示されている。box-x,-y,-zとして全部10としたものである。図23で画像の拡大縮小に関わらずこの正方形のサイズは変わらない。アホだな。図24のように全7点が点アイコンもこの画面では変色しているが,指定の黄色だ。

図23 右上の正方形はScalar field

図24 全7点が問題なく表示されている

 結局のところ,Scalar fieldの役割はわからないままだが,この3列を廃棄すればいいことがわかった。

6. drone cloudとgcpとの対応関係を見る

 上記のように, csvファイルの作成要素はわかったので,マックで編集して,図25のようにCloudCompareで読み込んだが,その結果が図26である。この画面では黄色というよりも明るい緑色になっている。Google Earthでもよく体験することだ。

図25 csvを読み込んだところ

図26 gcpが読み込めた

 これ以降については,研究に流れ込むので,プライベートページに移動する。

以上,22:11,2024/01/08記。




QgisとGrassGISの何れが3Dメッシュと点群ファイルを生かすことができるか which is better Qgis or GrassGIS for the display and analization of 3D mesh and cloud data ?

1. 期待していること

 GrassGISをかつては使ってきたがご無沙汰している。Qgisは1回だけ使って精度面で問題を多く感じてアプリも削除した。Qgisは大衆化が著しく,悪貨が良貨を駆逐する,勢いである。日本の行政はGrassGISを飛び越して,専らQGIS派である。行政のコンテンツがQgisを意識したものになっている。

 さて,3Dメッシュまたは点群ファイルへの対応はどうだろうか。もちろん,GrassGISではこのアプリ公開以来,早くから米軍やNASAのLiDARデータの処理が実行されてきたのであるが,ぼくはこのコンテンツには触れなかった。

 いま,求めるのは,3Dモデルの3D表示能力はもちろんである。CloudCompareはその点では問題はない。しかし,CloudCompareでは,3Dモデルに,他のベクトルデータをオーバーレイすることができない。これはかなり不便なことである。それで,GISソフトに期待してのこのページ作成である。

 繰り返すと,3Dモデル表示が自在で,他のベクトルデータをオーバーレイできるか。この観点で,両アプリを調べたいと思う。

 

2. GrassGISではどうか

2.1 ケース1

 まずは次の回答を紹介する。

Import Point-Cloud to Grass-Gis
Data Processing – Discussion and Q&A

MichaelL

Feb ’20

Got it. I know Lidartools (which LAStools is dependent on) was not ported for QGIS 3 so that what I was getting at. 10M points isn’t all that big so if you have at least 16GB of RAM you shouldn’t have any problems, but considering you are having problems with a 500 point subset then we know that combined with the scaling issue that it’s probably not your plug-ins. Have you tried downloading and importing an XYZ version?

You might try PDAL.
https://grass.osgeo.org/grass78/manuals/v.in.pdal.html 23

 QGIS 3では対応できない。GLASSではXYZ形式であれば扱えるのではということである。RGB情報を破棄してしまうのは余りに不便でGLASSでも使い物にならないということだ。

2.2 r.in.pdal – Creates a raster map from LAS LiDAR points using univariate statistics

 このタイトルのリンクは,GLASSのマニュアルページである。

The r.in.pdal module loads PDAL library supported point clouds (with emphasis on LiDAR LAS files) into a new raster map using binning. The user may choose from a variety of statistical methods which will be used for binning when creating the new raster.

と,あって,精度が落ちるという点で,この種の処理をぼくは求めていないので,脚下だ。なお,ここでのテーマはLiDAR測量の生データの処理であり,いま,agisoftで,plyファイルなどが手に入るので,このケース1は求めるものではない。

2.3 ケース2

 それでは,まずは,plyファイルをGrassGISが取り込めるかについて調べてみる。

v.in.ply – Creates a vector map from a PLY file.

DESCRIPTION

v.in.ply imports a vector map in PLY vector format. A PLY file always holds a number of vertices which are imported as points. PLY vertices can have a number of properties in addition to their coordinates. These properties are stored in an attribute table. For larger PLY files with many vertices (> 1000) it is highly recommended to not use DBF as database driver, but SQLite (default in GRASS GIS 7), PostgreSQL or MySQL, because the DBF driver is rather slow and can consume a lot of memory. The database driver can be set with db.connect.

NOTES

v.in.ply is designed for large point clouds with the possibility to have only coordinates, and no attribute table (for speed reasons).

と,あって,興味深いことを次に。

 GISデータの基本コンテンツと考えていたDBFが,GrassGIS7以降,デフォルトとして,使われていないことである。SQLite (default in GRASS GIS 7)。これには驚いた。ぼくは時代遅れになっていた。
 上記Descriptionを見て大変期待できると思ったが,上記Notesを見ると,only coordinates, and no attribute table (for speed reasons),とあって,ショックではあった。PLYファイルからRGB情報が削除されるのである。

 こう考えると,3D表示の可能性も見えないし,GISを使うメリットが感じられなくなったのである。

 v.ply.rectify – Imports PLY points, georeferences and exports them. も重要だ。

DESCRIPTION

v.ply.rectify imports a PLY point cloud, georeferences and exports it. The first three vertex properties must be the x, y, z coordinates with property names “x”, “y”, “z”, in this order.

A text file with Ground Control Points (GCPs) must exist in the same folder where the point cloud is located, and the textfile must have the same name like the point cloud, but ending on .txt instead of .ply.

The text file with GCPs must have the following format with one GCP per line: x y z east north height status with x, y, z as source coordinates and east, north, height as target coordinates. The status indicates whether to use a GCP (status is not zero) or not (status is zero). Entries must be separated by whitespace or tabs. Decimal delimiters must be points.

The georecitifictation method used is a 3D orthogonal rectification where angles are preserved. 3D objects are shifted, scaled and rotated, but no shear is introduced. Please read the output of the module, in particular the root mean square (RMS) errors.

v.ply.rectify optionally exports the georeferenced point cloud not only with real coordinates, but also with shifted coordinates (-s flag) for display in meshlab or similar software that can not deal with real coordinates. The exported PLY point clouds will be in the same folder like the input PLY point cloud.

 次の観点は重要である: x y z east north height status with x, y, z as source coordinates and east, north, height as target coordinates。通常のGISの座標系と異なり,x,y,z=easting, northing, heightになっている。別氏によればメタシェープでの処理も同様のようだ。

3. Qgisではどうか

 QGISで色付きの3Dメッシュを作成する は参考になったが,かなり,面倒な過程を経ている。わざわざ,他のカラー付きの画像を利用するのも理解できない。ここでは結局のところ,MeshLabが使われている。Qgisには3D表示能力がないようだ。この図1の右手のペーンには,レイヤーが配置されている。複数の3Dレイヤーを重ねて表示できるのである。不勉強でこの基礎的な機能についてさえ,知らなかった。

図1 上記SKラボ.netのページの「MeshLabで表示してみる」

 Blenderでも可能である。Blenderで点群(.ply)を色付きでレンダリングする Sep.23, 2023
出来上がった3Dオブジェクトを見ると,いかにも寂しい。

 で,SKラボ.netの,QGISで色付きの3Dメッシュを作成する,に出会った,と思った。しかし,元々RGB情報を持つplyを直接使わないのか,わからない。plyをインポートする際にRGB情報を廃棄してしまうからではあるが。結局,ぼくには無用のものであった。

4. MeshLabで

 結局,GLASSもQgisも使えないことがわかった。MeshLabに期待したい。MeshLabの操作法を学んで,このページに示したいと思う。

以上,19:15,2024/01/06記。

図2 show layer dialog boxを表示した

 図2のように,layer dialog boxを開いて,240103uawmisakibase.plyを
File > Import Mesh,で呼び出した。2分ほど。All files opened in 102,894 msecと表示されている。vertices: 198,925, 288だ。約2億個。

図3 応答なし

 Filters/Sampling/Clustered Vertex Subsampling, で,Cell Size1.0を例えば0.3などに変更してApplyしようとしたが,メーンメニューのFiltersをタップするだけで,応答無し,でフリーズしてしまう。重すぎるのだ。

 仕方無く,CloudCompareでFiles > Openを実行した。

5. CloudCompareでdown-samplingを

 試行錯誤してやっと見つけた。すごく簡単ではあった。Edit > subsample を実行すれば良い。開いたパネルで,min space btw points =,入力する数値であるが,2.7mにすると,出来上がったオブジェクトは7073点でスカスカだった。0.01mにすると全く変化なし。同じだ。で,0.100(10cm)にすると(図4は実行開始時, 図5は完成後),3,135,591点で,適当であった。元々2億点(198,925,288)からかなりスマートになった。

図4 min space = 10cmで実行

図5 3,135,591点からなるオブジェクト

 図6には保存した新たなPLYファイルが見えていて,5GBから238MBへとよほど軽くなっている。

図6 保存した新たなPLYファイルは238MBに

 では,MeshLabへ,という流れになるはずではあったが,偶然,Danielの凄いアドバイスを見つけてしまった。

Create a new layer of just ground points (already classified) である。

以上,0:15,2024/01/07記。

6. CloudCompareでgcpも表示することに

 上記の,Create a new layer of just ground points (already classified) ,の内容を次にコピペする。

Create a new layer of just ground points (already classified)
Joined: Tue Aug 09, 2022 10:34 am
Create a new layer of just ground points (already classified)
Post by okay » Tue Aug 09, 2022 10:36 am

Hi - I’m trying to make a separate layer of just the ground points in my point cloud. In other words, I want to be able to turn on and off a plot of just the ground points. I know about the CSF plug, but this isn’t necessary to use since the point cloud is classified. I can only figure out to how filter by value, but not by classification. Thanks.
Location: Grenoble, France
————————————————
Contact: Contact daniel
Re: Create a new layer of just ground points (already classified)
Post by daniel » Wed Aug 10, 2022 5:56 pm
Maybe you can use the new ‘qCloudLayers’ plugin?
Another way, maybe more simpler, is to use the ‘Edit > Scalar fields > Split cloud (integer values)’ method. This will create as many (sub) clouds as classes present in your cloud.

Daniel, CloudCompare admin
————————————————

 CloudCompareの開発者Danielさん,に,感謝。

 この情報を契機に,CloudCompareのPlugins Cloud layers,の存在を知った。図7のように,いま,使用中のCloudCompareの最上部には,CloudCompare v.2.12.4 (Kyiv) Stereo [64-bit],が表示されている。メーンメニューPluginsに見えるCloud layersが欲しいのだが,灰色だ。

図7 使用しているCloudCompareのPlugins: Cloud layers は灰色

 このページのテーマが遷ってきているので,別にページを立てることにする。後に,そのページのリンクを示す。

以上,13:28,2024/01/07記。

 




ドローン測量の限界の一例 an example of the limitation of drone measurements

はじめに

 共同研究者の一人beがドローンを飛ばして光波測距儀での測量についてはプリズムをボクが移動し,測点の記載を実施した。それはもう2020年春のことである。2018年と2019年にはkaさんと海岸地形の光波測距儀での測量などを実施した。両方のデータの整合性を取る必要があって,ドローン測量をここ1カ月以上に渡って,検討してきた。

 これはドローン測量のほぼファイナルな結果である。残念ながら期待した結果は得られていない。この結果に至るまでのやり取りなどや分析の詳細はプライベートページに掲載されている。まだ研究成果を発表するまでには時間が必要で,公開できない内容を多く含むためである。

図1 ドローン測量の限界の一例

ortho.tifの限界

 図1は,Agisoft Metashapeで作成されたオルソ写真と光波測距儀で求めた8点の位置のずれを示している。ビューワーはGoogle Earth (Pro)である。グーグルアースで表示されている画像はモザイク写真であり,写真画像で正確な位置を捉えることはできないので注意が必要である。

 画像とは関係なくグーグルアース画像での位置情報,繰り返すが画像とは関係のない,グーグルアース座標系での位置情報の精度は高い。グーグルアースプロの限界 の「GEProの位置精度」の章に示すように,地球面上の位置精度は1〜2mほどと言えるだろう。

 それゆえ,光波測距儀で得たgcp候補8点のうち採用した7点の位置は正しくグーグルアースにプロットされる。光波測距儀はノンプリズム対応のもので8点のうち2点はプリズムの設置が困難な場所であった。ドローン測量結果の位置決めはgcp7点を使用した。ドローン測量でえた空中写真で同地点で複数の写真にgcp位置を同定しており,ドローン測量結果はgcpとの対応が取られていて,グーグルアース上では,オルソ写真の画像上のgcp位置と一致する筈であった。ボクは当然一致すると考えていた。

 ところが,図1に示したように,gcp8点(図1の中央列のc点はドローン測量結果には使用されていないが球面位置情報は正しいので掲載した)のいずれも一致していない。図1では黄線分でズレを示しているが解像度が低くて見にくい。オルソ写真の位置と光波測距儀で決定した位置のズレを示している。参考データとしてGoogle画像で示されたスケールを使ってそのズレ値を例えば,e:1.3mなどと示している。Google画像での距離精度は低いがオルソ写真と光波測距儀での座標から得られた点のズレなので,精度はグーグルアースの精度ほど低くはないのであるが。

 で,図1の右端の画像はオルソ画像域全体のなかでの,ズレの傾向を示している。この画像の右下隅に方位ホウィールが見える。研究対象域は海岸線付近なので,gcpの点a,f, b,gに注目すると,1.85m±0.35mほど,北に変位している。

 これまで,種々の計算処理をしてきており,これ以上の対処は難しい。海岸地形の研究の観点からすると,このオルソ画像をいわばインデックスマップとして使わざるを得ないのである。このインデックスマップに,以前に測量したポイントの座標系を平行移動または回転などして行かざるを得ないのである。

今後の作業は思って

 このオルソ写真というか,ドローン測量結果が使えるかどうかは,別の観点から,評価する必要がある。この海岸で2018年または2019年に測量した結果,特に海抜高度が完全に一致しているかどうかである。球面位置はこのオルソに合わせるにしても,海抜高度に問題があれば,使えない。

 3DオブジェクトのPLYを使って,次に評価したいと思う。

以上,0:16,Dec. 24, 2023, to Christmas Eve。

グーグルアースから離脱か,Qgisか

 ズレ,1.0〜2.2mは,そういえば,グーグルアースの誤差範囲である。
参考:グーグルアースプロの限界
 図1ではgcpとオルソのズレは同じではないので,オルソ全体のズレだけでなく,gcpの個々のズレもある。グーグルアースの限界に基づくものと考えれば,納得もゆく。

 今日,be氏に電話して,Qgisにプロットしてもらうことにした。ぼくはGrassGISであるが,Qgisのユーザーは多く,簡便であり,Qgisにシフトした方がいいと思う。ゲームなどでも利用される3DオブジェクトのFBXやUSDZなどにも対応しているかも知れない。その確認もお願いした。

 ぼくはCloudCompareで,海岸線付近の高度情報を調べることになる。

以上,Dec. 24, 2023記。