説明

音響トモグラフィー計測システム及び音響トモグラフィー計測方法

【課題】どのような環境下においても、河川の流量等の物理量を高効率かつ高精度に長期計測する。
【解決手段】音響トモグラフィー計測システム100は、一方の河岸に設置された上流のトランスデューサ2Bと他方の河岸に設置された下流のトランスデューサ2Aとの間で指向性のない音波を送受信する。時間情報算出部21は、音波が送信された時刻を基準とする、各トランスデューサ2A、2Bにおいて受信される音波から復調されるデータと所定のデータとの相関波形に基づいて、音波の伝播時間に関する情報を算出する。平均情報算出部24は、音波の音線長Lと、平均伝播時間tm、伝播時間差Δtに基づいて、河川の横断面における断面平均音速cm及び断面平均流速umを算出する。物理量算出部25は、算出された断面平均音速cm及び断面平均流速umに基づいて、河川の横断面に関する物理量を算出する。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、河川の横断面における流量等の物理量を計測する音響トモグラフィー計測システム及び音響トモグラフィー計測方法に関する。
【背景技術】
【0002】
治水、河川の環境管理、水資源管理等の観点から、河川の流量を精度良く計測できる方法の登場が望まれている。近年では、河川の流量を計測する方法として、水位流量曲線(H−Q)法、ADCP(Acoustic Doppler Current Profiler)法、H−ADCP法、AVM法等の様々な方法が提案されている(例えば、非特許文献1参照)。
【0003】
この他、発明者等は、塩水が混じっているような感潮域における河川流量をモニタリングする方法を提案している(例えば、非特許文献2、3参照)。
【先行技術文献】
【非特許文献】
【0004】
【非特許文献1】二瓶、「超音波ドップラー流速分布計に基づく実河川流速・流量モニタリング」、ながれ26(2007)
【非特許文献2】川西ら、「次世代超音波流速計による感潮域の流量と水温・塩分の連続モニタリング」、水工学論文集、第53巻、2009年2月
【非特許文献3】川西ら、「次世代超音波流速計による感潮河川流量の長期連続モニタリング」、河川技術論文集、第15巻、2009年6月
【発明の概要】
【発明が解決しようとする課題】
【0005】
しかしながら、H−Q法は、感潮域では用いることができず、水位から流量を求めるための関係式のキャリブレーションが必要となる。
【0006】
また、ADCP法では、スポット観測に限定され、長期連続観測や、洪水時の観測は困難である。
【0007】
また、H−ADCP法では、塩水遡上時には音波の屈折、洪水時には浮遊土砂が観測の障害となる。
【0008】
さらに、AVM法では、更正係数が必要となり、洪水時には欠測となることがある。また、装置が大型となり高価になる。
【0009】
本発明は、どのような環境下においても、河川の流量等の物理量を高効率かつ高精度に長期計測することができる音響トモグラフィー計測システム及び音響トモグラフィー計測方法を提供することを目的とする。
【課題を解決するための手段】
【0010】
上記目的を達成するために、本発明の第1の観点に係る音響トモグラフィー計測システムは、
一方の河岸側に設置された上流の音波送受信器と他方の河岸側に設置された下流の音波送受信器との間で送受信される指向性のない音波によって河川の横断面に関する物理量を計測する音響トモグラフィー計測システムであって、
前記各音波送受信器から前記音波が送信された時刻を基準とする、前記各音波送受信器において受信される前記音波から復調されるデータと所定のデータとの相関波形に基づいて、前記音波送受信器間における前記音波の伝播時間に関する情報を算出する時間情報算出部と、
前記上流の音波送受信器と前記下流の音波送受信器との間で送受信される前記音波の経路の長さである音線長と、前記時間情報算出部によって算出された前記伝播時間に関する情報とに基づいて、前記河川の横断面における断面平均音速及び断面平均流速との少なくとも一方の情報を算出する平均情報算出部と、
前記平均情報算出部によって算出された情報に基づいて、前記河川の横断面に関する物理量を算出する物理量算出部と、
を備える。
【0011】
前記時間情報算出部は、
前記音波の伝播時間に関する情報として、前記音波送受信器間における音波の平均伝播時間と、前記音波送受信器間における音波の伝播時間差とを算出し、
前記平均情報算出部は、
前記音線長と、前記平均伝播時間と、前記伝播時間差とに基づいて、前記断面平均流速を算出し、
前記物理量算出部は、
前記断面平均流速に対して前記横断面の断面積を乗算することにより、前記横断面の流量を前記物理量として算出する、
こととしてもよい。
【0012】
前記時間情報算出部は、
前記音波の伝播時間に関する情報として、前記音波送受信器間における前記音波の平均伝播時間を算出し、
前記平均情報算出部は、
前記音線長から、前記平均伝播時間を除算することにより、前記断面平均音速を算出し、
前記物理量算出部は、
前記断面平均音速、塩分、水温及び水深の間の関係式を用いて、塩分及び水温のいずれかを前記物理量として算出する、
こととしてもよい。
【0013】
前記横断面の測量情報と水温と塩分の分布から求まる前記横断面の音速分布とに基づいて、音線解析を行うことにより、前記音線長を算出する音線解析部をさらに備え、
前記平均情報算出部は、
前記音線解析部で算出された前記音線長に基づいて、前記断面平均音速及び前記断面平均流速の少なくとも一方を算出する、
こととしてもよい。
【0014】
コンダクティビティ・テンプラチャ・デプス・プロファイラ(CTD)により計測される水温、塩分に基づいて基準時刻での前記断面平均音速を算出し、前記断面平均音速と、前記基準時刻での前記音波送受信器間における音波の平均伝播時間とに基づいて、前記音線長を算出する音線長算出部をさらに備え、
前記平均情報算出部は、
前記音線長算出部で算出された前記音線長に基づいて、前記断面平均音速及び前記断面平均流速の少なくとも一方を算出する、
こととしてもよい。
【0015】
前記時間情報算出部は、
前記上流の音波送受信器から前記音波が発せられてから、前記下流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりにおいて閾値に達するまでの第1の時間と、
前記下流の音波送受信器から前記音波が発せられてから、前記上流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりにおいて前記閾値に達するまでの第2の時間と、を求め、
前記第1の時間と前記第2の時間との平均を、前記音波送受信器間における音波の平均伝播時間として算出し、
前記第1の時間と前記第2の時間との差を、前記音波送受信器間における音波の伝播時間差として算出する、
こととしてもよい。
【0016】
前記横断面の測量情報と、水温と塩分の分布から求まる前記横断面の音速分布とに基づいて、音線解析を行うことにより、前記音線長を算出する音線解析部と、
前記音線解析部によって算出された前記音線長と、コンダクティビティ・テンプラチャ・デプス・プロファイラ(CTD)により計測される水温及び塩分から求まる基準時刻での前記断面平均音速とに基づいて、前記音波送受信器間における音波の平均伝播時間を算出し、算出された前記平均伝播時間と、前記下流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりと、前記上流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりとの関係に基づいて、前記閾値を決定する閾値算出部をさらに備え、
前記時間情報算出部は、
前記閾値算出部によって決定された前記閾値を用いて、前記第1、第2の時間を算出する、
こととしてもよい。
【0017】
前記物理量算出部によって算出された前記物理量の時系列データに含まれる異常値データをウェーブレット変換を行って検出し、
前記物理量の時系列データのうち、異常値データを除く他のデータに基づく統計モデルに対する前記異常値データの残差に基づいて前記異常値データを修正するスパイク修正部をさらに備える、
こととしてもよい。
【0018】
前記時間情報算出部は、GPS衛星からの時刻信号に基づいて、前記音波の伝播時間に関する情報を算出する、
こととしてもよい。
【0019】
前記所定のデータは、M系列法により生成された疑似乱数から構成されている、
こととしてもよい。
【0020】
また、本発明の第2の観点に係る音響トモグラフィー計測方法は、
一方の河岸側に設置された上流の音波送受信器と他方の河岸側に設置された下流の音波送受信器との間で送受信される指向性のない音波によって河川の横断面に関する物理量を計測する音響トモグラフィー計測方法であって、
前記各音波送受信器から前記音波が送信された時刻を基準とする、前記各音波送受信器において受信される前記音波から復調されるデータと所定のデータとの相関波形に基づいて、前記音波送受信器間における前記音波の伝播時間に関する情報を算出する時間情報算出工程と、
前記上流の音波送受信器と前記下流の音波送受信器との間で送受信される前記音波の経路の長さである音線長と、前記時間情報算出工程において算出された前記伝播時間に関する情報とに基づいて、前記河川の横断面における断面平均音速及び断面平均流速との少なくとも一方の情報を算出する平均情報算出工程と、
前記平均情報算出工程において算出された情報に基づいて、前記河川の横断面に関する物理量を算出する物理量算出工程と、
を含む。
【発明の効果】
【0021】
本発明によれば、河川の横断面を音波の導波路と仮定して横断面の物理量を計測する。すなわち、一方の河岸側に設置された上流の音波送受信器と他方の河岸側に設置された下流の音波送受信器との間で送受信される指向性のない音波の伝播時間に関する情報が求められ、音線長と、伝播時間に関する情報とに基づいて断面平均音速又は断面平均流速が求められ、それらの情報に基づいて河川の横断面に関する物理量が求められる。
【0022】
音波は、水面や河床で反射、屈折しながら横断面を伝播し、様々な経路を経て音波送受信器で受信される。このため、音波の伝播時間に関する情報は、その横断面における断面平均音速や断面平均流速と相関関係のある値となる。したがって、音線長が既知であれば、音波の伝播時間に関する情報から、断面平均音速や断面平均流速を精度良く求めることができる。これは、洪水時であっても、感潮域であっても同じである。
【0023】
音波の伝播時間の情報は、受信される音波から復調されるデータと所定のデータとの相関波形に基づいて求められる。このようにすれば、SN比を向上することができる。
【0024】
また、伝播時間に関する情報と断面平均音速又は断面平均流速との間の関係式や、断面平均音速又は断面平均流速と河川の横断面の物理量との間の関係式には、更正が必要な係数が含まれていないので、どのような状況下にあっても、その演算式を用いて河川の横断面の物理量を精度良く算出することができる。
【0025】
以上のように、本発明によれば、一対の音波送受信器だけで、洪水時や感潮域であってもキャリブレーションを行うことなく、河川の横断面の物理量を計測することができるので、どのような環境下においても、河川の横断面の物理量を高効率かつ高精度に長期計測することができる。
【図面の簡単な説明】
【0026】
【図1】図1(A)は、本発明の第1の実施形態に係る音響トモグラフィー計測システムにおけるトランスデューサの配置を模式的に示す斜視図であり、図1(B)は、トランスデューサの配置を説明するための河川の横断面図である。
【図2】本発明の第1の実施形態に係る音響トモグラフィー計測システムの構成を示すブロック図である。
【図3】図3(A)は、一方の時系列データ作成部で生成される相関波形の一例を示すグラフであり、図3(B)は、他方の時系列データ作成部で生成される相関波形の一例を示すグラフであり、図3(C)は、閾値CCと平均伝播時間tmと、伝播時間差Δtとを示す図である。
【図4】本発明の第2の実施形態に係る音響トモグラフィー計測システムの構成を示すブロック図である。
【図5】本発明の第3の実施形態に係る音響トモグラフィー計測システムの構成を示すブロック図である。
【図6】図6(A)及び図6(B)は、閾値の求め方の一例を説明するための図である。
【図7】本発明の第4の実施形態に係る音響トモグラフィー計測システムの構成を示すブロック図である。
【図8】図7のスパイク修正部のメイン処理のフローチャートである。
【図9】図8のスパイク検出処理のフローチャートである。
【図10】図8のスパイク置換処理のフローチャートである。
【図11】スパイクデータが修正される様子を説明するためのグラフである。
【発明を実施するための形態】
【0027】
この発明の実施の形態について、図面を参照して詳細に説明する。
【0028】
(第1の実施形態)
まず、本発明の第1の実施形態について説明する。
【0029】
図1(A)、図1(B)には、本実施形態における一対の音波送受信器としてのトランスデューサ2A、2Bの配置が示されている。トランスデューサ2Aは、一方の河岸3A側に設置され、トランスデューサ2Bは、他方の河岸3B側に設置されている。
【0030】
トランスデューサ2Aは、トランスデューサ2Bに対して下流に設置されている。川の流れの方向を示すベクトルと、トランスデューサ2A、2Bとを結ぶ線分とのなす角度をθとする。
【0031】
トランスデューサ2Aは、トランスデューサ2Bに向かって音波を送信し、トランスデューサ2Bは、トランスデューサ2Aに向かって音波を送信する。また、トランスデューサ2Aは、トランスデューサ2Bから送信された音波を受信し、トランスデューサ2Bは、トランスデューサ2Aから送信された音波を受信する。
【0032】
トランスデューサ2A、2Bから送信される音波は、指向性がなく、放射状に水中を伝播する。
【0033】
図1(B)における横断面は、図1(A)において、トランスデューサ2A、2Bとを通過する横断面となっている。図1(B)に示すように、トランスデューサ2A、2Bから送信された音波は、水平に進むものもあるが、水面や河床で反射、屈折しながら進むものもある。音波はこのようにして横断面を伝播し、様々な経路を通って、他方のトランスデューサ2A、2Bに到達する。
【0034】
ここで、音波の経路の長さを、以下では音線長Lと呼ぶ。
【0035】
なお、図1(B)では、音波の経路は直線で図示されているが、河口などの塩水が入り交じった塩分濃度の違いにより階層化された感潮域では、音波は水中で屈折するので、その経路は曲線を描くようになり、さらに複雑なものになる。
【0036】
図2には、本実施形態に係る音響トモグラフィー計測システム100の概略的な構成が示されている。図2に示すように、音響トモグラフィー計測システム100は、上述のトランスデューサ2A、2Bに加え、送受信回路4A、4Bと、河川情報取得部5をさらに備える。
【0037】
トランスデューサ2Aは、音波を送信する送波部2A1(トランスミッタ)と、音波を受信する受波部2A2(ハイドロフォン)とを備える。トランスデューサ2Bは、音波を送信する送波部2B1(トランスミッタ)と、音波を受信する受波部2B2(ハイドロフォン)とを備える。
【0038】
送受信回路4Aは、GPS(Global Positioning System)アンテナ10、時刻管理部11、送信制御部12、送信回路部13、受信回路部14及び受信制御部15を備える。
【0039】
GPSアンテナ10は、GPS衛星(不図示)からのGPS情報を受信する。このGPS情報には、時刻情報が含まれている。
【0040】
時刻管理部11は、GPSアンテナ10によって受信されたGPS情報に含まれる時刻情報に基づいて、基準クロック信号を生成して送信制御部12及び受信制御部15に出力している。
【0041】
送信制御部12は、自局を識別するための識別コードを含むデータを送信回路部13に出力する。この識別コードは、疑似乱数発生法の1つであるM系列法により生成されたものである。送信回路部13では、その識別コードをD/A変換し、さらに位相変調してトランスデューサ2Aの送波部2A1から、特定の送信時刻に、識別コードに対応する音波が送信される。
【0042】
一方、受信回路部14は、トランスデューサ2Aの受波部2A2で受波された信号をA/D変換して復調し、受信制御部15に出力している。
【0043】
受信制御部15は、受信回路部14から出力されたデジタルデータと、予め記憶されている識別コード(M系列)との相関度を算出する。より具体的には、受信回路部14から出力されたデジタルデータをベクトルxi(i=1〜n;nは2以上の自然数)とし、識別コードをベクトルyi(i=1〜m;mは2のM系列の次数乗;m≦n)として、次式に示すベクトルxiとベクトルyiとの相互相関係数を、相関度として算出する。
【数1】


ここで、相関度は0〜1までの値をとり、その値が高ければ両者の相関性が高いということになる。算出された相関度は、そのときの時刻情報とともに随時出力される。
【0044】
送受信回路4Bの構成は、送受信回路4Aと同じである。
【0045】
河川情報取得部5は、CPU及びメモリを備える。メモリに格納されたプログラムを実行することにより、以下の各構成要素の機能が実現される。
【0046】
時系列データ作成部20A、20B、時間情報算出部21、閾値保持部22、音線解析部23、平均情報算出部24及び物理量算出部25を備える。
【0047】
時系列データ作成部20Aは、送受信回路4Aから出力された時刻情報とその時刻の相関度とに基づいて、相関波形(相関度の時間変化を示す波形)を生成する。図3(A)には、時系列データ作成部20Aによって生成された相関波形の一例が示されている。
【0048】
図3(A)において、時刻0は、トランスデューサ2Bから音波が送信された時刻を示している。また、時刻0で送信された音波は、その後、トランスデューサ2Aで受信されるが、音波が受信された時刻近傍において、相関波形に鋭いピークが表れている。
【0049】
一方、時系列データ作成部20Bは、送受信回路4Bから出力された時刻情報とその時刻の相関度とに基づいて、相関波形を生成する。図3(B)には、時系列データ作成部20Bによって生成された相関波形の一例が示されている。また、時刻0でトランスデューサ2Aから送信された音波は、その後、トランスデューサ2Bで受信されるが、その音波が受信された時刻近傍において、相関波形に鋭いピークが表れている。
【0050】
図3(A)のグラフと図3(B)のグラフとを重ね合わせたのが図3(C)のグラフである。図3(C)に示すように、時系列データ作成部20Bで生成される相関波形のピークは、時系列データ作成部20Aで生成される相関波形のピークよりも遅れている。この遅れは、川の流速により発生するドップラー効果によるものである。
【0051】
例えば、河川の横断面中の音波の平均音速(断面平均音速)をcmとし、川の平均流速の音線に沿った方向の成分(断面平均流速)をumとすると、上流のトランスデューサ2Bから下流のトランスデューサ2Aへの音波の伝播時間tAと、下流のトランスデューサ2Aから上流のトランスデューサ2Bへの音波の伝播時間tBとは、それぞれ次式のようになる。
【数2】


ここで、Lは上述のように音線長である。
【0052】
この場合、断面平均音速cmと、断面平均流速umは次式のように求められる。
【数3】


ここで、tm、Δtは次式のようになる。
【数4】


以下では、tmを平均伝播時間といい、Δtを伝播時間差という。すなわち、tmは、トランスデューサ2A、2Bで音波が受信されたときのそれぞれの伝播時間の平均であり、Δtは、トランスデューサ2Aで音波が受信された時刻と、トランスデューサ2Bで音波が受信された時刻との時間差である。
【0053】
時間情報算出部21は、図3(C)に示すように、トランスデューサ2A、2Bから音波が送信された時刻0を基準として、時系列データ作成部20A、20Bから読み込んだ相関波形から、伝播時間に関する情報として、平均伝播時間tm、伝播時間差Δtを求める。
【0054】
図3(C)に示すように、平均伝播時間tm、平均伝播時間Δtを求めるためには、閾値CCが必要である。閾値CCは閾値保持部22に保持されている。時間情報算出部21は、閾値保持部22から閾値CCを読み込む。時間情報算出部21は、上流のトランスデューサ2Bから音波が発せられて(時刻0)から、下流のトランスデューサ2Aで受信される音波の相関波形のピークの立ち下がりにおいて閾値CCに達するまでの第1の時間としてのtAを求める。さらに、時間情報算出部21は、下流のトランスデューサ2Aから音波が発せられて(時刻0)から、上流のトランスデューサ2Bで受信される音波の相関波形のピークの立ち下がりにおいて閾値CCに達するまでの第2の時間としてのtBを求める。
【0055】
さらに、時間情報算出部21は、上記式(5)を用いて、平均伝播時間tmと伝播時間差Δtとを求める。すなわち、時間情報算出部21は、tAとtBとの平均を、トランスデューサ2A、2B間における音波の平均伝播時間tmとして算出する。また、時間情報算出部21は、tAとtBとの差を、トランスデューサ2A、2B間における音波の伝播時間差Δtとして算出する。
【0056】
平均伝播時間tmと伝播時間差Δtとを求めると、今度は、断面平均音速cmと、断面平均流速umを求めなければならない。断面平均音速cm、断面平均流速umを求めるためには、音線長Lを得る必要がある。音線解析部23は、予め測量されている横断面の測量情報(水深等の情報)と、水温分布f(r、z)及び塩分分布g(r、z)から求まる横断面内の音速分布とに基づいて、音線解析を行うことにより、音線長Lを算出する。ここで、rとzはそれぞれ横断面の水平と鉛直方向の位置座標である。なお、水温分布f(r、z)及び塩分分布g(r、z)は、後述するCTD6によってある間隔で横断面内の離散値として測定される水温及び塩分を補間することにより求められる。
【0057】
より具体的には、音線解析部23は、水温分布f(r、z)及び塩分分布g(r、z)に基づいて、例えば、以下の式(Medwinの式)を用いて、音速分布c(r、z)を算出する。
【数5】


ここで、Dは水深である。
【0058】
このようにして、求められた音速分布c(r、z)を用いて、音線の水平からの角度(入射補角φと、音波の鉛直座標z、音波の伝播時間tが、次式(7a)〜(7c)から求められる。
【数6】


ここで、式(7a)は、スネルの法則を表したものである。音線解析部23は、上記式(7a)、式(7b)の連立方程式を解くことにより、音波の軌跡上におけるdz/drの値を求める。続いて、音線解析部23は、求められたdz/drを、次の平面曲線長を求める式に代入することにより、音線長Lを求める。
【0059】
【数7】


ここで、Rは、トランスデューサ間の直線距離である。
【0060】
平均情報算出部24は、このようにして求められた音線長L、平均伝播時間tm及び伝播時間差Δtを、上記式(3)、式(4)にそれぞれ代入して、断面平均音速cmと、断面平均流速umを算出する。すなわち、平均情報算出部24は、上流のトランスデューサ2Bと下流のトランスデューサ2Aとの間で送受信される音線長Lと、時間情報算出部21によって算出された平均伝播時間tm及び伝播時間差Δtに基づいて、河川1の横断面における断面平均音速cm及び断面平均流速umとの少なくとも一方の情報を算出する。
【0061】
物理量算出部25は、平均情報算出部24によって算出された断面平均流速umに基づいて、次式を用いて横断面の断面積A(h)を乗算して、横断面の流量Qを物理量の1つとして算出する。ここで、hは水位である。
【数8】

【0062】
また、物理量算出部25は、断面平均音速cmに基づいて、例えば以下の式(Medwinの式)を用いて、河川の水温Tm又は塩分Smを算出する。
【数9】


ここで、Tmは水温(℃)、Smは塩分、Dは水深である。
【0063】
CTD6は、コンダクティビティ・テンプラチャ・デプス・プロファイラであり、電気伝導度、温度、水深を観測する装置である。CTD6では、電気伝導度、水温等から塩分が算出される。
【0064】
物理量算出部25は、CTD6から水温Tm又は塩分Smのいずれかを入力する。物理量算出部25は、入力した断面平均音速cmと、水温Tm又は塩分Smのいずれかを上記式(10)に入力し、残った水温Tm又は塩分Smを算出する。
【0065】
ところで、相関波形が、図3(A)、図3(B)に示すように、単峰性の波形(単一ピーク)になるのは、横断面を通過した音線の到達時間がほぼ同じ場合である。これは、放射状に発射された音波の伝播距離がほとんど同じで、かつ、音速も一定である、塩水遡上のない単断面河川を測定するときに見られる現象である。なお、河川の水平面内を迂回した経路でトランスデューサ2A、2Bに到達する音波は、減衰が大きいため、無視することができる。
【0066】
一方、洪水時の複断面河川の高水敷にトランスデューサ2A、2Bを置いた場合、高水敷上のみを通過して届く音波と、低水路に入り込みその底面で反射して届く音波とでは、到達時間が大幅に異なるため、図3(A)、図3(B)に示すようにはならず、2つのピークが出現する場合がある。この場合、断面平均流速umを求めるためには、1番目のピークではなく、2番目(後)のピーク波形を用いて、平均伝播時間tm、伝播時間差Δtを求めるのが望ましい。
【0067】
塩水遡上のある感潮河川では、音波の伝播状況が複雑となるため、上述のような音線解析を行って各音線の到達時間を推定することにより、平均伝播時間tm、伝播時間差Δtを求めるためのピークを選択するようにしてもよい。また、複数のピークのうち、一番大きなものを選択するようにしてもよい。
【0068】
次に、本実施形態に係る音響トモグラフィー計測システム100の動作について説明する。
【0069】
送受信回路4A、4Bの送信制御部12は、時刻管理部11から出力される時刻情報を参照し、特定の時刻0になると、送信回路部13を介してトランスデューサ2A、2Bの送波部2A1、2B1に音波を送信させる。
【0070】
しばらく時間が経ち、トランスデューサ2A、2Bの受波部2A2、2B2で音波を受信すると、受信回路部14で音波が復調され、受信制御部15でその相関度が算出され、時刻管理部11から出力される時刻情報とともに出力される。
【0071】
河川情報取得部5の時系列データ作成部20A、20Bは、それぞれの相関波形を作成し、時間情報算出部21に出力する。時間情報算出部21は、閾値保持部22によって保持される閾値CCを用いて、各相関波形のピークの立ち下がりが閾値CCに到達する時刻tA、tBを求め、時刻tA、tBに基づいて、平均伝播時間tm、伝播時間差Δtを算出する(時間情報算出工程)。
【0072】
平均情報算出部24は、平均伝播時間tm、伝播時間差Δtと、音線解析部23から出力された音線長Lとに基づいて、断面平均音速cm、断面平均流速umを算出する(平均情報算出工程)。
【0073】
物理量算出部25は、断面平均流速umに横断面の断面積を乗算して、流量Qを算出する。また、物理量算出部25は、断面平均音速cmと、CTD6から出力された水温Tm又は塩分Smのいずれかに基づいて、水温Tm又は塩分Smを算出する(物理量算出工程)。
【0074】
本実施形態によれば、河川1の横断面を音波の導波路と仮定して横断面の物理量を計測する。すなわち、一方の河岸3B側に設置された上流のトランスデューサ2Bと他方の河岸3A側に設置された下流のトランスデューサ2Aとの間で送受信される指向性のない音波の伝播時間に関する情報が求められ、音線長Lと、平均伝播時間tm、伝播時間差Δtとに基づいて断面平均音速cm又は断面平均流速umが求められ、それらの情報に基づいて河川1の横断面に関する物理量が求められる。
【0075】
音波は、水面や河床で反射、屈折しながら横断面を伝播し、様々な経路を経てトランスデューサ2A、2Bで受信される。このため、断面平均音速cm又は断面平均流速umは、その横断面における断面平均音速cmや断面平均流速umと相関関係のある値となる。したがって、音線長Lが既知であれば、断面平均音速cm又は断面平均流速umから、断面平均音速cmや断面平均流速umを精度良く求めることができる。これは、洪水時であっても、感潮域であっても同じである。
【0076】
平均伝播時間tm、伝播時間差Δtは、受信される音波から復調されるデータとそのM系列のデータとの相関波形に基づいて求められる。このようにすれば、SN比を格段に向上することができる。
【0077】
また、平均伝播時間tm、伝播時間差Δtとに基づいて断面平均音速cm又は断面平均流速umとの間の関係式や、断面平均音速cm又は断面平均流速umと河川の横断面の物理量(流量Q、水温Tm又は塩分Sm)との間の関係式には、更正が必要な係数が含まれていないので、どのような状況下にあっても、その演算式を用いて河川の横断面の物理量(流量Q、水温Tm又は塩分Sm)を精度良く算出することができる。
【0078】
以上のように、本発明によれば、一対のトランスデューサ2A、2Bだけで、洪水時や感潮域であってもキャリブレーションを行うことなく、河川1の横断面の物理量を計測することができるので、どのような環境下においても、河川の横断面の物理量を高効率かつ高精度に長期計測することができる。
【0079】
(第2の実施形態)
次に、本発明の第2の実施形態について説明する。
【0080】
図4には、本実施形態に係る音響トモグラフィー計測システム100の構成が示されている。図4に示すように、本実施形態に係る音響トモグラフィー計測システム100は、音線解析部23の代わりに、音線長算出部26を備える点が、上記第1の実施形態に係る音響トモグラフィー計測システム100の構成と異なる。
【0081】
音線長算出部26は、CTD6により計測された基準時刻での水温Tm0及び塩分Sm0に基づいて、基準時刻での断面平均音速cm0を算出し、時間情報算出部21で算出された基準時刻でのトランスデューサ2A、2B間の平均伝播時間tm0とに基づいて、音線長L0を(=cm0×tm0)算出する。
【0082】
平均情報算出部24は、音線長算出部26で算出された音線長L0に基づいて、時刻tmにおける断面平均音速cm及び断面平均流速umを、上記第1の実施形態と同様にして算出する。
【0083】
本実施形態に係る音響トモグラフィー計測システム100は、音線解析により音線長を求めるのが困難なときに好適である。例えば、洪水の観測を、水門や橋などからトランスデューサ2A、2Bを下ろして使うときなどに特に好適である。
【0084】
(第3の実施形態)
次に、本発明の第3の実施形態について説明する。
【0085】
図5には、本実施形態に係る音響トモグラフィー計測システム100の構成が示されている。図5に示すように、音響トモグラフィー計測システム100は、閾値算出部27を備えている点が、上記第1の実施形態に係る音響トモグラフィー計測システム100と異なる。
【0086】
閾値算出部27は、音線解析部23によって算出された音線長Lと、CTD6により計測された基準時刻での水温T0及び塩分S0から求まる基準時刻での断面平均音速cm0とに基づいて平均伝播時間tm0(=L/cm0)を算出する。さらに、閾値算出部27は、算出された平均伝播時間tm0を用いて、時系列データ作成部20A、20Bから入力した相関波形に基づいて、下流のトランスデューサ2Aで受信される音波の相関波形のピークの立ち下がりと、上流のトランスデューサ2Bで受信される音波の相関波形のピークの立ち下がりとの関係に基づいて、閾値CCを決定する(図3(C)参照)。
【0087】
より具体的には、図6(A)に示すように、平均伝播時間tm0を基準として、マイナス方向に時間dtだけずらしたときの下流のトランスデューサ2Aで受信される音波の相関波形のピークの立ち下がりの値と、プラス方向に時間dtだけずらしたときの上流のトランスデューサ2Bで受信される音波の相関波形のピークの立ち下がりの値とが同じ値となる時間dt(=Δt/2)を探索する。図6(B)には、両者の値が同じ値となったときの時間dtが示されている。閾値算出部27は、図6(B)に示すように、そのときの相関波形の値を閾値CCとして決定する。
【0088】
時間情報算出部21は、閾値算出部27によって決定された閾値CCを用いて、時刻tA、tBを算出する。
【0089】
本実施形態によれば、計測対象の河川1の横断面の水温T0や塩分S0に応じて閾値CCを決定することができるので、より高精度に、音波の平均伝播時間tmや、伝播時間差Δtを求めることができる。
【0090】
(第4の実施形態)
次に、本発明の第4の実施形態について説明する。
【0091】
図7には、本実施形態に係る音響トモグラフィー計測システム100の構成が示されている。図7に示すように、音響トモグラフィー計測システム100は、スパイク修正部28を備えている点が、上記第1の実施形態に係る音響トモグラフィー計測システム100と異なる。
【0092】
スパイク修正部28は、物理量算出部25によって算出された物理量(流量Q、水温Tm、塩分Sm)の時系列データに含まれる異常値データ(スパイクデータ)をウェーブレット変換を行って検出する。続いて、スパイク修正部28は、物理量の時系列データのうち、異常値データ(スパイクデータ)を除く他のデータに基づく統計モデルに対する異常値データ(スパイクデータ)の残差に基づいて異常値データを修正する。
【0093】
図8には、スパイク修正部28のメイン処理が示されている。スパイク修正部28は、このメイン処理を所定の間隔で行う。
【0094】
図8に示すように、スパイク修正部28は、物理量算出部25から出力された物理量の時系列データの読み出しを行う(ステップS1)。続いて、スパイク修正部28は、スパイク検出処理のサブルーチンを行う(ステップS2)。続いて、スパイク修正部28は、スパイク検出処理において、スパイクデータが検出されたか否かを判定する(ステップS3)。スパイクデータが検出されていた場合(ステップS3;Yes)、スパイク修正部28は、スパイクデータの検出結果を一時記憶する(ステップS4)。続いて、スパイク修正部28は、スパイク置換処理のサブルーチンを行う(ステップS5)。スパイク置換処理実行後、スパイク修正部28は、ステップS2に戻る。
【0095】
このようにして、ステップS3においてスパイクデータが検出されないと判定されるまで、ステップS2→S3→S4→S5が繰り返される。
【0096】
スパイクデータが検出されなかった場合(ステップS3;No)、スパイク修正部28は、スパイクデータが修正された(又はスパイクデータが検出されなかった)物理量の時系列データを記憶し(ステップS6)、処理を終了する。
【0097】
図9には、ステップS2のスパイク検出処理のサブルーチンが示されている。図9に示すように、スパイク修正部28は、物理量の時系列データを読み出す(ステップS11)。続いて、スパイク修正部28は、読み出した物理量の時系列データが長期間のデータであるか否かを判定する(ステップS12)。この判定が肯定されると(ステップS12;Yes)、スパイク修正部28は、物理量の時系列データから低周波成分を除去する(ステップS14)。
【0098】
続いて、スパイク修正部28は、ウェーブレット変換用のフィルタを選択し(ステップS15)、時系列データに対してウェーブレット変換を行う(ステップS16)。続いて、スパイク修正部28は、所定の演算を行って閾値を算出し(ステップS17)、その閾値を用いてフィルタリング処理を行う(ステップS18)。その後、スパイク修正部28は、ウェーブレット変換されたデータに対してウェーブレット逆変換を行う(ステップS19)。
【0099】
続いて、スパイク修正部28は、ウェーブレット逆変換により得られた物理量の時系列データと、元の物理量の時系列データとを比較することにより、スパイクデータの位置(スパイク位置)を検出し(ステップS20)、スパイクデータ列の長さ(スパイク長)を検出する(ステップS21)。
【0100】
続いて、スパイク修正部28は、スパイクデータ列前の正常値データ(正常値データ)のデータ長(正常値データ長)を検出する(ステップS22)。そして、スパイク修正部28は、これまでの検出結果を記憶する(ステップS23)。ステップS23実行後、スパイク修正部28は、スパイク検出処理を終了する。
【0101】
図10には、スパイク置換処理のサブルーチンが示されている。図10に示すように、スパイク修正部28は、スパイクデータ周辺の時系列データの読み出しを行う(ステップS31)。続いて、スパイク修正部28は、読み出した時系列データが長期間のデータであるか否かを判定する(ステップS32)。この判定が肯定されれば(ステップS32;Yes)、スパイク修正部28は、そのスパイクデータ列前の正常値データを選択する(ステップS33)。
【0102】
続いて、スパイク修正部28は、選択された正常値データを用いて自己回帰モデル(ARモデル)の6つのパラメータを推定する(ステップS34)。続いて、スパイク修正部28は、選択された正常値データを用いて条件付き最尤推定法を用いたモデル推定を行う(ステップS35)。続いて、スパイク修正部28は、選択された正常値データを用いてAIC(赤池情報量規準)による最適モデルを推定する(ステップS36)。
【0103】
続いて、スパイク修正部28は、ステップS34、S35、S36でそれぞれ求められた統計モデルに対する、時系列データの残差の相関性をチェックし、残差の違いが所定の許容範囲であるか否かをテストする(ステップS37)。さらに、スパイク修正部28は、モデルの定常性をテストする(ステップS38)。
【0104】
続いて、スパイク修正部28は、選択されたデータが、ステップS36、S37の2つのテストに合格したか否かを判定する(ステップS39)。
【0105】
両テストに合格していれば(ステップS39;Yes)、スパイク修正部28は、モデルエラー(残差)を算出し(ステップS44)、そのモデルエラーに基づいてスパイクデータに置換されるべき予測値を算出する(ステップS45)。そして、スパイク修正部28は、スパイクデータを予測値に置換する(ステップS46)。図11には、短期的に変動する物理量について、スパイクデータ(異常値データ)が、予測値に修正される様子が模式的に示されている。
【0106】
続いて、スパイク修正部28は、全ての物理量の時系列データの選択が終了したか否かを判定する(ステップS47)。この判定が否定されると(ステップS47;No)、スパイク修正部28は、ステップS31に戻る。
【0107】
スパイク修正部28は、再び物理量の時系列データの部分的な読み出しを行い(ステップS31)、読み出した時系列データが長期間のデータであるか否かを判定する(ステップS32)。この判定が否定されれば(ステップS32;No)、スパイク修正部28は、正常値データに対して、3次多項式補間を行ってスパイクデータに置換されるべき補間値を求め(ステップS42)、スパイクデータを補間値に置換する(ステップS43)。その後、スパイク修正部28は、全てのデータの選択が終了したか否かを判定する(ステップS47)。この判定が否定されると(ステップS47;No)、スパイク修正部28は、ステップS31に戻る。
【0108】
その後、ステップS31→S32→S33→S34→S35→S36→S37→S38と進み、ステップS39での判定が否定された場合、スパイク修正部28は、テスト失敗が初回であるか否かを判定する(ステップS40)。テスト失敗が初回であれば(ステップS40;Yes)、スパイク修正部28は、正常値データをさらに追加選択し(ステップS41)、ステップS34に戻る。
【0109】
以降、データが追加された状態で、ステップS34→S35→S36で3つのモデルが推定され、残差の相関性テストと、定常性テストが行われる(ステップS38)。そして、両テスト合格判定が否定され(ステップS39;No)、テスト失敗が初回でない場合(ステップS40;No)、スパイク修正部28は、正常値データに対して、3次多項式補間を行ってスパイクデータに置換されるべき補間値を求め(ステップS42)、スパイクデータを補間値に置換する(ステップS43)。その後、スパイク修正部28は、全てのデータの選択が終了したか否かを判定する(ステップS47)。
【0110】
このようにして、ステップS47において全データが選択されたと判定されるまで、処理が繰り返され、全データが選択されたと判定されると(ステップS47;Yes)、スパイク修正部28は、スパイクデータが置換された物理量の時系列データを処理結果として記憶し(ステップS48)、スパイク置換処理を終了する。
【0111】
このように、本実施形態によれば、求められた物理量の時系列データに含まれる異常値データを修正することができる。これにより、高精度な計測が可能になる。このようにすれば、例えば、流量Q、水温Tm、塩分Smの短期変動などを正確に検出することが可能となる。また、この方法は、取得されたデータが少ないため、異常値データをリジェクトすることが適切でない場合にも有効である。
【0112】
なお、本実施形態では、自己回帰(AR)モデル、条件付き最尤推定法によるモデル、AICによるモデル等を、物理量の時系列データの統計モデルとして適用したが、これに限られず、カルマンフィルタによるモデルを適用してもよい。
【0113】
本実施形態に係る音響トモグラフィー計測システム100は、例えば、塩水遡上のある感潮河川において、音波の伝播の双方向性が崩れ、双方のトランスデューサ2A、2Bからそれぞれ得られる相関波形の相似でなくなるために発生するスパイクデータの低減に特に有効である。
【0114】
また、上記各実施形態では、時刻情報をGPS衛星からのものとした。このようにすれば、極めて正確に音波の伝播時間等を計測することができる。
【0115】
なお、上記実施の形態において、実行されるプログラムは、フレキシブルディスク、CD−ROM(Compact Disk Read-Only Memory)、DVD(Digital Versatile Disk)、MO(Magneto-Optical Disk)等のコンピュータ読み取り可能な記録媒体に格納して配布し、そのプログラムをインストールすることにより、上述の処理を実行するシステムを構成することとしてもよい。
【0116】
また、プログラムをインターネット等の通信ネットワーク上の所定のサーバ装置が有するディスク装置等に格納しておき、例えば、搬送波に重畳させて、ダウンロード等するようにしてもよい。
【0117】
また、上述の機能を、OS(Operating System)が分担して実現する場合又はOSとアプリケーションとの協働により実現する場合等には、OS以外の部分のみを媒体に格納して配布してもよく、また、ダウンロード等してもよい。
【0118】
なお、本発明は、上記実施の形態及び図面によって限定されるものではない。本発明の要旨を変更しない範囲で実施の形態及び図面に変更を加えることができるのはもちろんである。
【産業上の利用可能性】
【0119】
本発明は、流量、水温、塩分等の河川の横断面に関する物理量の計測に最適であり、特に、河口近くに塩分が入り交じった感潮域や、洪水時等における河川の物理量の計測に最適である。
【符号の説明】
【0120】
1 河川
2A、2B トランスデューサ
2A1、2B1 送波部
2A2、2B2 受波部
3A、3B 河岸
4A、4B 送受信回路
5 河川情報取得部
6 CTD(コンダクティビティ・テンプラチャ・デプス・プロファイラ)
10 GPSアンテナ
11 時刻管理部
12 送信制御部
13 送信回路部
14 受信回路部
15 受信制御部
20A、20B 時系列データ作成部
21 時間情報算出部
22 閾値保持部
23 音線解析部
24 平均情報算出部
25 物理量算出部
26 音線長算出部
27 閾値算出部
28 スパイク修正部
100 音響トモグラフィー計測システム

【特許請求の範囲】
【請求項1】
一方の河岸側に設置された上流の音波送受信器と他方の河岸側に設置された下流の音波送受信器との間で送受信される指向性のない音波によって河川の横断面に関する物理量を計測する音響トモグラフィー計測システムであって、
前記各音波送受信器から前記音波が送信された時刻を基準とする、前記各音波送受信器において受信される前記音波から復調されるデータと所定のデータとの相関波形に基づいて、前記音波送受信器間における前記音波の伝播時間に関する情報を算出する時間情報算出部と、
前記上流の音波送受信器と前記下流の音波送受信器との間で送受信される前記音波の経路の長さである音線長と、前記時間情報算出部によって算出された前記伝播時間に関する情報とに基づいて、前記河川の横断面における断面平均音速及び断面平均流速との少なくとも一方の情報を算出する平均情報算出部と、
前記平均情報算出部によって算出された情報に基づいて、前記河川の横断面に関する物理量を算出する物理量算出部と、
を備える音響トモグラフィー計測システム。
【請求項2】
前記時間情報算出部は、
前記音波の伝播時間に関する情報として、前記音波送受信器間における音波の平均伝播時間と、前記音波送受信器間における音波の伝播時間差とを算出し、
前記平均情報算出部は、
前記音線長と、前記平均伝播時間と、前記伝播時間差とに基づいて、前記断面平均流速を算出し、
前記物理量算出部は、
前記断面平均流速に対して前記横断面の断面積を乗算することにより、前記横断面の流量を前記物理量として算出する、
ことを特徴とする請求項1に記載の音響トモグラフィー計測システム。
【請求項3】
前記時間情報算出部は、
前記音波の伝播時間に関する情報として、前記音波送受信器間における前記音波の平均伝播時間を算出し、
前記平均情報算出部は、
前記音線長から、前記平均伝播時間を除算することにより、前記断面平均音速を算出し、
前記物理量算出部は、
前記断面平均音速、塩分、水温及び水深の間の関係式を用いて、塩分及び水温のいずれかを前記物理量として算出する、
ことを特徴とする請求項1に記載の音響トモグラフィー計測システム。
【請求項4】
前記横断面の測量情報と水温と塩分の分布から求まる前記横断面の音速分布とに基づいて、音線解析を行うことにより、前記音線長を算出する音線解析部をさらに備え、
前記平均情報算出部は、
前記音線解析部で算出された前記音線長に基づいて、前記断面平均音速及び前記断面平均流速の少なくとも一方を算出する、
ことを特徴とする請求項1乃至3のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項5】
コンダクティビティ・テンプラチャ・デプス・プロファイラ(CTD)により計測される水温、塩分に基づいて基準時刻での前記断面平均音速を算出し、前記断面平均音速と、前記基準時刻での前記音波送受信器間における音波の平均伝播時間とに基づいて、前記音線長を算出する音線長算出部をさらに備え、
前記平均情報算出部は、
前記音線長算出部で算出された前記音線長に基づいて、前記断面平均音速及び前記断面平均流速の少なくとも一方を算出する、
ことを特徴とする請求項1乃至3のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項6】
前記時間情報算出部は、
前記上流の音波送受信器から前記音波が発せられてから、前記下流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりにおいて閾値に達するまでの第1の時間と、
前記下流の音波送受信器から前記音波が発せられてから、前記上流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりにおいて前記閾値に達するまでの第2の時間と、を求め、
前記第1の時間と前記第2の時間との平均を、前記音波送受信器間における音波の平均伝播時間として算出し、
前記第1の時間と前記第2の時間との差を、前記音波送受信器間における音波の伝播時間差として算出する、
ことを特徴とする請求項1乃至3のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項7】
前記横断面の測量情報と、水温と塩分の分布から求まる前記横断面の音速分布とに基づいて、音線解析を行うことにより、前記音線長を算出する音線解析部と、
前記音線解析部によって算出された前記音線長と、コンダクティビティ・テンプラチャ・デプス・プロファイラ(CTD)により計測される水温及び塩分から求まる基準時刻での前記断面平均音速とに基づいて、前記音波送受信器間における音波の平均伝播時間を算出し、算出された前記平均伝播時間と、前記下流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりと、前記上流の音波送受信器で受信される前記音波の相関波形のピークの立ち下がりとの関係に基づいて、前記閾値を決定する閾値算出部をさらに備え、
前記時間情報算出部は、
前記閾値算出部によって決定された前記閾値を用いて、前記第1、第2の時間を算出する、
ことを特徴とする請求項6に記載の音響トモグラフィー計測システム。
【請求項8】
前記物理量算出部によって算出された前記物理量の時系列データに含まれる異常値データをウェーブレット変換を行って検出し、
前記物理量の時系列データのうち、異常値データを除く他のデータに基づく統計モデルに対する前記異常値データの残差に基づいて前記異常値データを修正するスパイク修正部をさらに備える、
ことを特徴とする請求項1乃至7のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項9】
前記時間情報算出部は、GPS衛星からの時刻信号に基づいて、前記音波の伝播時間に関する情報を算出する、
ことを特徴とする請求項1乃至8のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項10】
前記所定のデータは、M系列法により生成された疑似乱数から構成されている、
ことを特徴とする請求項1乃至9のいずれか一項に記載の音響トモグラフィー計測システム。
【請求項11】
一方の河岸側に設置された上流の音波送受信器と他方の河岸側に設置された下流の音波送受信器との間で送受信される指向性のない音波によって河川の横断面に関する物理量を計測する音響トモグラフィー計測方法であって、
前記各音波送受信器から前記音波が送信された時刻を基準とする、前記各音波送受信器において受信される前記音波から復調されるデータと所定のデータとの相関波形に基づいて、前記音波送受信器間における前記音波の伝播時間に関する情報を算出する時間情報算出工程と、
前記上流の音波送受信器と前記下流の音波送受信器との間で送受信される前記音波の経路の長さである音線長と、前記時間情報算出工程において算出された前記伝播時間に関する情報とに基づいて、前記河川の横断面における断面平均音速及び断面平均流速との少なくとも一方の情報を算出する平均情報算出工程と、
前記平均情報算出工程において算出された情報に基づいて、前記河川の横断面に関する物理量を算出する物理量算出工程と、
を含む音響トモグラフィー計測方法。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate

【図11】
image rotate


【公開番号】特開2012−98093(P2012−98093A)
【公開日】平成24年5月24日(2012.5.24)
【国際特許分類】
【出願番号】特願2010−244673(P2010−244673)
【出願日】平成22年10月29日(2010.10.29)
【出願人】(504136568)国立大学法人広島大学 (924)
【Fターム(参考)】