(57)【要約】【課題】スモールループコイルで地下を電磁探査するとき、オフセット誤差を容易且つ高精度に処理できる電磁探査装置、方法、プログラムを提供することを目的とする。【解決手段】本発明は、地下の比抵抗分布モデルを仮定してこれを基に磁場強度を計算し、この計算された磁場強度と実測した磁場強度とを比較し、差が所定量以下に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推定する電磁探査装置であって、更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対してオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して更新された比抵抗分布モデルの磁場強度とオフセット量との差が所定量内に収まるまで比抵抗分布モデルの更新と計算とを行う逆解析手段を備えたことを主要な特徴とする。【選択図】図１
2009074953.tif 000002
【特許請求の範囲】
【請求項１】
１次磁場を形成するための送信コイルが設けられた磁場発生部と、前記１次磁場により地中に形成される２次磁場の磁場強度を検出する受信コイルが設けられた磁場検出部と、前記受信コイルが検出する磁場強度から前記１次磁場の寄与分をキャンセルするためのバッキングコイルと、を備え、
地下の比抵抗分布モデルを仮定してこれを基に磁場強度を計算し、この計算された磁場強度と前記受信コイルによって実測した磁場強度とを比較し、この２つの磁場強度の差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推定する電磁探査装置であって、
更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して前記寄与分を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して前記更新された比抵抗分布モデルの磁場強度と前記オフセット量との差が所定量内に収まるまで比抵抗分布モデルの更新と計算とを行う逆解析手段を備え、予めキャリブレーションを行わずに比抵抗分布を逆解析することを特徴とする電磁探査装置。
【請求項２】
前記受信コイルの高度を変えて二位置以上で磁場強度を検出したとき、前記逆解析手段が、この検出した二位置以上の磁場強度を磁場強度として統合すると共に、該磁場強度に対して独立したオフセット量として関係付けることを特徴とする請求項１記載の電磁探査装置。
【請求項３】
１次磁場を形成して、該１次磁場の寄与分をキャンセルしながら地中に形成される２次磁場の磁場強度を実測すると共に、地下の比抵抗分布モデルを仮定し、これに基づいて磁場強度を計算し、この計算した磁場強度と実測した磁場強度とを比較し、この２つの磁場強度の差が所定量以下に収斂するまで前記比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推測する電磁探査方法であって、
更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して前記１次磁場の影響を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して前記更新された磁場強度と前記オフセット量との差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とを行い、予めキャリブレーションを行わずに比抵抗分布を逆解析することを特徴とする電磁探査方法。
【請求項４】
高度を変えて２位置以上で磁場強度を検出し、この２位置以上の磁場強度を独立した対等な磁場強度として統合すると共に、該磁場強度に対して独立したオフセット量を関係付けることを特徴とする請求項３の電磁探査方法。
【請求項５】
コンピュータに、
実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算させる第１手順と、
更新された比抵抗分布のラアプラシアンを計算し前記１次磁場の影響を補正するためのオフセット量を関係付ける第２手順と、
更新された比抵抗分布モデルの磁場強度に対してヤコビアンを計算し、前記オフセット量の影響度を調整係数で調整して前記差との関係付けを行う行列を取得する第３手順と、
前記行列をＱＲ分解して前記比抵抗分布の変化分と前記オフセット量を得る第４手順と、
前記比抵抗分布の変化分と前記オフセット量がそれぞれ所定量内に収斂するまで比抵抗分布モデルを更新し、計算を繰り返す第５手順と、
を実行させるためのプログラム。
【請求項６】
前記第１手順が、高度を変えて２位置以上で磁場強度を検出し、この２位置以上の磁場強度を独立した対等な磁場強度として統合した実測の磁場強度とし、
前記第２手順が、更新された比抵抗分布モデルの磁場強度に対して前記１次磁場のほかに独立したオフセット量を関係付けることを特徴とする請求項５のプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、本発明は、電磁誘導作用を利用して地中探査を行い、２次磁場の測定値から地下の比抵抗分布を高精度に推定する電磁探査装置であって、とくに１次磁場の影響を除いたとき残留するオフセット誤差等を演算処理によって容易に処理できる電磁探査装置、及びその電磁探査方法、及びコンピュータに地下の比抵抗分布を逆解析させるためのプログラムに関する。
【背景技術】
【０００２】
従来、地中探査法にはＥＭ法（Electromagnetic Method）やＭＴ法（Magneto telluric Method）などが知られている。ＥＭ法とは、地盤に対向させた送信側のループコイルに交番電流を流して１次磁場を発生させ、これに伴って地中に発生する渦電流によって形成される２次磁場を受信側で観測し、これによって地盤の比抵抗等を測定する地中探査方法である。
【０００３】
これに対し、ＭＴ法は、地中で起こる電磁誘導現象により変動する電場と磁場を地面に接触させた電極とアンテナで観測して地盤の比抵抗等を測定する地中探査方法である。いずれも、電磁誘導を利用した地中探査法であるが、ＥＭ法は地面に電極を打ち込んだりしない。なお、ＭＴ法と同様に、電極を打ち込んで比抵抗分布を推定するものに電気探査法がある。従って、測定が容易になことから、ＥＭ法は地中の比抵抗の分布を推定する場合に多用されている。
【０００４】
このようにＥＭ法によれば、人間がループコイル持って移動すると低コストで手軽に空中で地中探査が行える。地形の影響をあまり受けず、比抵抗分布の平面的な広がりと深度方向の変化を容易に把握できるために、３次元な探査が可能な優れた方法である。
【０００５】
しかし、人間によって手軽に地中探査が行えると言うメリットがある反面、広範囲の地中探査を行うのにはあまり馴染まない。このため、ヘリコプターでループコイルを曳航する多周波数型の空中電磁法（helicopterboren electromagnetic method）が提案されている。これは空中でループコイルに交流電流を流し、交流磁場が地中を通過する際に生ずる電磁誘導現象を利用して地下の比抵抗分布を求めるものである。
【０００６】
但し、この空中電磁法によるデータ処理と解析は、多数の各種測定値のデータベースの作成、測定データと位置データとの統合、測定ＥＭ成分（各周波数における同相成分Ｉ、離相成分Ｑ）のデータの相互のレベルを調整するレベリング、各測定値の確定、比抵抗計算、最後に解析の順に行われる。この中で、レベリングは比抵抗計算解析の中で最も困難で時間の要する煩わしいものである。
【０００７】
すなわち、比抵抗を計算するには２次磁場の正しい同相成分Ｉ、及び離相成分Ｑが必要であるが、従来からこのデータの０レベルを調整するための決定的な方法はない。その理由を説明すると、ＥＭ法においては、１次磁場（ＨＰ）に対する２次磁場（ＨＳ）の割合ＨＳ／ＨＰをＩ、Ｑに分離してｐｐｍ単位で測定している。このため、図２に示すようにＨＰの大きさを既知として、受信側でＨＰをキャンセル（相殺）してＨＳのみを測定し、間接的にＨＳ／ＨＰを求める方法をとっている。
【０００８】
図２を参照して説明すると（詳細は後述する）、送信コイルの近くに小さなバッキングコイル（受信コイルと逆巻きに直列に接続）を設けて１次磁場の影響を相殺し、２次磁場を増幅して同相成分と離相成分とに分けて、１次磁場に対する同相成分比Ｉ（ｐｐｍ）と離相成分比Ｑ（ｐｐｍ）とを測定している。
【０００９】
従って、ＨＳ／ＨＰの値を直接測定しているわけではないので、残留ＨＰによるオフセット誤差等のΔＨＰがどうしても残り問題となる。そこで受信側の１次磁場ＨＰをバッキングコイルで相殺した後の電位を見掛け上の０レベル（相対０電位）としている。図３においてはΔＨＰ＝（ＨＰ−Ｈ’’Ｐ）が相当する。ここで−Ｈ’’は相殺のための磁場である。
【００１０】
このオフセット誤差等は工場出荷時に初期設定されているが、例えば温度や地形など、場所が変われば多くの環境要因が影響し、一定値にはならない。現場（温度、地表からの距離などのパラメータ）で補正される必要がある。送信コイル、バッキングコイル、受信コイルの関係位置は決まっているため、地表面からの２次磁場の影響がなくなる高い高度、通常例えば１０ｍ以上の樹木などにループコイルを設置して補正するのが最良の方法とされている。ただし、近くにこうした樹木等が存在するときのみ可能な方法である。しかし、この方法でも完全な誤差は残る。要するに、結果として、現場で補正可能な誤差と、避けることのできない誤差があり、ＥＭ方法の問題となっている。
【００１１】
こうした問題に対処するため、次の空中電磁法のレベリングにおけるドリフト補正方法が提案された（特許文献１参照）。特許部文献１のドリフト補正方法は、測定飛行前に通常のキャリブレーション（オフセット誤差やその他の誤差）を行い、高高度で０レベルの設定を行ってから、調査対象地域等の平坦地を選んで垂直多点測定によるドリフト補正のためのデータを取得する。
【００１２】
次に通常の測定飛行を行って調査終了後に着陸する直前に、垂直多点測定を行った同一箇所で、同一内容の垂直多点測定を行い、データ処理時に通常のレベリング手法を用いてそれぞれ０レベルを設定し、実測値の大きさを求め、測定飛行前後の垂直多点測定値それぞれについて、インバージョン（逆解析）によりキャリブレーション時の利得及び位相誤差によるドリフトを求め、ドリフト補正を施すことにより正確な応答値を得る方法である。ドリフトには測定システムの誤差や残留磁場の影響、ノイズ部分がある。
【００１３】
インバージョンは、垂直多点測定を行った地点の地下が水平多層構造からなると仮定し、各周波数毎のドリフトの数と水平多層構造パラメータの数の総計を未知数とし、垂直多点測定によるデータ数が前記未知数の総数を上回るようにして、ドリフトを含む全未知数を最小２乗法で求めている。
【００１４】
【特許文献１】特開２００３−４３１５７号公報
【発明の開示】
【発明が解決しようとする課題】
【００１５】
さて、以上説明した特許部文献１のドリフト補正方法は、ヘリコプターによる対地高度２０ｍ〜１５０ｍの連続的測定で、データ数が未知数の総数を上回ったとき全未知数を最小二乗法で解く多周波数型の空中電磁法（ＨＥＭ）の補正方法であり、未知数を線形の最適化問題として解いている。具体的には多数の未知変量に初期値を与え、演算を繰り返し、最小解を求めている。
【００１６】
しかし、対地高度誤差を除いたドリフトだけですらα、β、γ、δなど４個の未知数を有し、情報不足であるため、ヘリコプターによる対地高度２０ｍ〜１５０ｍの連続的測定を行い、データ数が未知数の総数を上回ったとき全未知数を最小二乗法で解くため、演算量が膨大で、高能力の計算機と、またコストが度外視されてもよいとき、のような探査方法である。なお、補正時に高度が上がると応答値は３乗に逆比例して悪化し、データ量が増加したからと言って、データの誤差は相対的に大きくなる。
【００１７】
これに対し、人間がループコイルを移動させるスモールループコイル式の電磁探査方法においては、簡易に地中探査が可能になる。もし、このスモールループコイル式の電磁探査方法でオフセット誤差等が問題なく処理され、分析結果が高精度であれば、今後新たな電磁探査方法の途が開ける。
【００１８】
そこで本発明は、スモールループコイルで地下を電磁探査するとき、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる電磁探査装置、また、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる電磁探査方法、及び、比抵抗分布を高精度に推定でき、オフセット量の処理を容易に実行できるプログラムを提供することを目的とする。
【課題を解決するための手段】
【００１９】
本発明は、上記課題を解決するために、１次磁場を形成するための送信コイルが設けられた磁場発生部と、１次磁場により地中に形成される２次磁場の磁場強度を検出する受信コイルが設けられた磁場検出部と、受信コイルが検出する磁場強度から１次磁場の寄与分をキャンセルするためのバッキングコイルと、を備え、地下の比抵抗分布モデルを仮定してこれを基に磁場強度を計算し、この計算された磁場強度と受信コイルによって実測した磁場強度とを比較し、この２つの磁場強度の差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推定する電磁探査装置であって、更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して寄与分を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して更新された比抵抗分布モデルの磁場強度とオフセット量との差が所定量内に収まるまで計算とモデル更新を行う逆解析手段を備え、予めキャリブレーションを行わずに比抵抗分布を逆解析することを主要な特徴とする。
【００２０】
また同様に、本発明は、１次磁場を形成して、該１次磁場の寄与分をキャンセルしながら地中に形成される２次磁場の磁場強度を実測すると共に、地下の比抵抗分布モデルを仮定し、これに基づいて磁場強度を計算し、この計算した磁場強度と実測した磁場強度とを比較し、この２つの磁場強度の差が所定量以下に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推測する電磁探査方法であって、更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して１次磁場の影響を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して更新された磁場強度とオフセット量との差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とをおこない、予めキャリブレーションを行わずに比抵抗分布を逆解析することを主要な特徴とする。
【００２１】
さらに、本発明のプログラムは、コンピュータに、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算させる第１手順と、更新された比抵抗分布モデルのラアプラシアンを計算し１次磁場の影響を補正するためのオフセット量を関係付ける第２手順と、更新された比抵抗分布モデルの磁場強度に対してヤコビアンを計算し、オフセット量の影響度を調整係数で調整しながら関係付けを行う行列を取得する第３手順と、行列をＱＲ分解して比抵抗分布の変化分とオフセット量、係数を得る第４手順と、比抵抗分布の変化分とオフセット量、係数がそれぞれ所定量内に収斂するまで比抵抗分布モデルを更新し、計算を繰り返す第５手順と、を実行させることを主要な特徴とする。
【発明の効果】
【００２２】
本発明の電磁探査装置によれば、スモールループコイルで地下を電磁探査するとき、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる。
【００２３】
また、本発明の電磁探査方法によれば、スモールループコイルで地下を電磁探査するとき、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる。
【００２４】
さらに、本発明のプログラムによれば、地下の比抵抗分布を高精度に推定でき、オフセット量の処理を容易に実行できる。
【発明を実施するための最良の形態】
【００２５】
本発明の第１の形態は、１次磁場を形成するための送信コイルが設けられた磁場発生部と、１次磁場により地中に形成される２次磁場の磁場強度を検出する受信コイルが設けられた磁場検出部と、受信コイルが検出する磁場強度から１次磁場の寄与分をキャンセルするためのバッキングコイルと、を備え、地下の比抵抗分布モデルを仮定してこれを基に磁場強度を計算し、この計算された磁場強度と受信コイルによって実測した磁場強度とを比較し、この２つの磁場強度の差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推定する電磁探査装置であって、更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して寄与分を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して更新された比抵抗分布モデルの磁場強度とオフセット量との差が所定量内に収まるまで比抵抗分布モデルの更新と計算とを行う逆解析手段を備え、予めキャリブレーションを行わずに比抵抗分布を逆解析することを特徴とする電磁探査装置である。この構成によって、スモールループコイルで地下を電磁探査するとき、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる。
【００２６】
本発明の第２の形態は、第１の形態に従属する形態であって、受信コイルの高度を変えて二位置以上で磁場強度を検出したとき、逆解析手段が、この検出した二位置以上の磁場強度を磁場強度として統合すると共に、該磁場強度に対して独立したオフセット量として関係付けることを特徴とする電磁探査装置である。この構成によって、オフセット誤差の調整や測定高さを気にすることなく、二位置以上の磁場強度の情報を統合できる。
【００２７】
本発明の第３の形態は、１次磁場を形成して、該１次磁場の寄与分をキャンセルしながら地中に形成される２次磁場の磁場強度を実測すると共に、地下の比抵抗分布モデルを仮定し、これに基づいて磁場強度を計算し、この計算した磁場強度と実測した磁場強度とを比較し、この２つの磁場強度の差が所定量以下に収斂するまで比抵抗分布モデルの更新と計算とを続けて地下の比抵抗分布を推測する電磁探査方法であって、更新された比抵抗分布モデルに基づいて初期値付近で線形化された磁場強度を計算し、この計算された磁場強度に対して１次磁場の影響を補正するためのオフセット量を関係付け、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算し、この差に対して更新された磁場強度とオフセット量との差が所定量内に収斂するまで比抵抗分布モデルの更新と計算とを行い、予めキャリブレーションを行わずに比抵抗分布を逆解析することを特徴とする電磁探査方法である。この構成によって、スモールループコイルで地下を電磁探査するとき、２次磁場に基づいて地下の比抵抗分布を高精度に推定でき、１次磁場の影響を除くとき残留するオフセット誤差等を演算処理によって容易に処理できる。
【００２８】
本発明の第４の形態は、第３の形態に従属する形態であって、高度を変えて２位置以上で磁場強度を検出し、この２位置以上の磁場強度を独立した対等な磁場強度として統合すると共に、該磁場強度に対して独立したオフセット量を関係付けることを特徴とする電磁探査方法である。この構成によって、オフセット誤差の調整や測定高さを気にすることなく、二位置以上の磁場強度の情報を統合できる。
【００２９】
本発明の第５の形態は、コンピュータに、実測した磁場強度と更新前の比抵抗分布モデルの磁場強度との差を計算させる第１手順と、更新された比抵抗分布のラアプラシアンを計算し１次磁場の影響を補正するためのオフセット量を関係付ける第２手順と、更新された比抵抗分布モデルの磁場強度に対してヤコビアンを計算し、オフセット量の影響度を調整係数で調整して上記差との関係付けを行う行列を取得する第３手順と、行列をＱＲ分解して比抵抗分布の変化分とオフセット量を得る第４手順と、比抵抗分布の変化分とオフセット量がそれぞれ所定量内に収斂するまで比抵抗分布モデルを更新し、計算を繰り返す第５手順と、を実行させるためのプログラムである。この構成によって、地下の比抵抗分布を高精度に推定でき、オフセット量の処理を容易に実行できる。この第４手順において、比抵抗分布の変化分とオフセット量を得るだけでなく、調整係数の値も取得し、また、第５手順で比抵抗分布の変化分とオフセット量を得るだけでなく、調整係数の値が所定量内に収斂するようにするのも好適である。
【００３０】
本発明の第６の形態は、第５の形態に従属する形態であって、第１手順が、高度を変えて２位置以上で磁場強度を検出し、この２位置以上の磁場強度を独立した対等な磁場強度として統合した実測の磁場強度とし、第２手順が、更新された比抵抗分布モデルの磁場強度に対して１次磁場のほかに独立したオフセット量を関係付けることを特徴とするプログラムである。この構成によって、オフセット誤差の調整や測定高さを気にすることなく、二位置以上の磁場強度の情報を統合できる。
【００３１】
（実施の形態１）
以下、本発明の実施の形態１における電磁探査装置、電磁探査方法及びプログラムについて説明をする。図１は本発明の実施の形態１における電磁探査装置のブロック構造図、図２は測定の原理の説明図、図３は１次磁場と２次磁場、及びキャンセルのための磁場とオフセット量との関係を示す説明図である。
【００３２】
図１において、１は実施の形態１における電磁探査装置の１次磁場を発生するための磁場発生部であり、２は所定周波数の電流を流すスモールループの送信コイルである。本発明の実施の形態１では、１ｋＨｚ〜１０ｋＨｚ程度から選択された周波数の交番電流（交互に極性が変わる所定の電流であればよい）で１次磁場を形成し、この１次磁場が形成する渦電流により２次磁場を生成している。
【００３３】
次に、図１に示す３は、この電磁探査装置によって生成された磁場を検出するための磁場検出部であり、４は磁場強度を測定するためのスモールループの受信コイルである。そして、５は受信コイルと巻回方向は逆に巻かれているが、受信コイルと直列に接続されたバッキングコイル、６は制御部、７は得られた情報を分析する分析装置である。
【００３４】
そして、図２は送信コイル２、受信コイル４、バッキングコイル５と１次磁場ＨＰ（その磁場強度をＨＰとする）、２次磁場Ｈｓ（その磁場強度をＨｓとする）の関係を示している。送信コイル２に磁界発生部１から交番電流、例えば複素数表現でｉ＝ｉ０ｅｊωｔが流されると、１次磁場強度ＨＰ＝ＨＰ０ｅｊωｔが発生する。このように磁束が変化することで、これを打ち消すように渦電流ｉｗが発生する。さらにこの渦電流Ｉｗが２次磁場Ｈｓ＝Ｈｓ０ｅｊωｔを発生する。しかし、受信コイル４には、上述した１次磁場とこの２次磁場が加わったＨＰ＋Ｈｓ＝ＨＰ０ｅｊωｔ＋Ｈｓ０ｅｊωｔの磁束が測定される。１次磁場ＨＰに対する２次磁場Ｈｓのオーダはｐｐｍオーダである。
【００３５】
そこで、受信コイル４と直列に接続されたバッキングコイル５を利用して、受信コイル４で受信する磁束から１次磁場ＨＰの影響をキャンセル（相殺）する。すなわち、電磁誘導でＨＰ’＝ＨＰ０’ｅｊωｔ≒ＨＰにより、相殺のための電流を発生させ、受信コイル４で（−ＨＰ’’）＝（−ＨＰ０’’ｅｊωｔ）≒（−ＨＰ）の相殺磁界を生成している。これにより、受信コイル４では（ＨＰ＋ＨＳ−ＨＰ０’’）＝（ＨＰ０ｅｊωｔ＋Ｈｓ０ｅｊωｔ−ＨＰ０’’ｅｊωｔ）の磁場となり、誤差（オフセット誤差）をΔＨＰ＝（ＨＰ−Ｈ’’Ｐ）とすると、（Ｈｓ＋ΔＨＰ）＝（Ｈｓ０ｅｊωｔ＋ΔＨＰ）となり、キャリブレーションによってΔＨＰ≒０にすることが可能であるなら、受信コイルでＨｓ＝Ｈｓ０ｅｊωｔを検出できることが可能になる。
【００３６】
以上の構成を基にして、以下、実施の形態１における電磁探査装置、電磁探査方法がどのようにして比抵抗分布を求めるのか、について以下簡単に説明する。地盤中に電気伝導率σ（Ｓ／ｍ）の物質（ここでσ＝１／ρ（Ωｍ）、ρは比抵抗）が存在するとし、空気と地盤という２つの媒体が地表を境に接しており、空気中で地盤に向けて電磁波を送信するため、境界条件を満たしつつ、マクスウェルの方程式を解くことにより、地盤からの応答（磁場強度）の応答値（理論値）を計算できる。
【００３７】
しかし、この方程式は数値計算により差分や有限要素法を用いて解くしかない。そこで、地下の比抵抗分布（地下モデル）を仮定し（フォワードモデリング）、電流ｉを送信コイル４に流して１次磁界ＨＰを形成し、この１次磁場ＨＰに対する地下モデルでの渦電流ｉｗの発生、さらにはこれに由来して２次的に派生する２次磁場Ｈｓが計算できる。これに対して、この２次磁界のＨｓの実測値が得られているので、この実測値と応答値（理論値）との差がなくなるように比抵抗分布モデルを修正し、これを収斂させれば地中の比抵抗分布の可能性が大である。言い換えれば地下モデルを仮定して逆解析することが実施の形態１における電磁探査方法と言うことになる。
【００３８】
さて、２次磁界ＨＳを測定するとき、１次磁場ＨＰの影響が十分キャンセルできているか、またその他の誤差はどれほどかが問題になる。すなわち、２次磁場ＨＳの上記誤差ΔＨＰが問題となる。バッキングコイル５を利用して１次磁場ＨＰの影響はほぼ解消できるが、これは完全でなく、誤差が残る。システム誤差等のようにキャリブレーションで解決できる誤差もあるが、現地環境の影響（温度、地磁気、その他の自然環境など）による誤差等を完全に取り除くことはできない。
【００３９】
そこで、本発明においては、こうした誤差を分析的に追求するのではなく、このような誤差をすべて含めて、地下モデルの比抵抗分布のデータと同様のデータと位置付け、地下モデルの要素として処理を行い、すべてに矛盾がないように逆解析するものである。
【００４０】
以下この説明を行う。まず、マクスウェルの方程式を解くことに関して、であるが、地盤中での磁場と比抵抗が関係していることを具体的に示すと、均質な大地（半無限媒質）で周波数ｆが低いと仮定したとき、磁場強度Ｈｓ、ＨＰ、比抵抗ρ（Ωｍ）との間に次の（数１）の関係がある。ここでμ０は真空の透磁率（４π×１０−７Ｈ／ｍ）であり、ｓは送信コイル２と受信コイル４のコイル間隔（ｍ）である。
【００４１】
【数１】
2009074953.tif 000003
この（数１）が示唆するように、それぞれが所定の比抵抗を有する複数のブロックから構成される地下モデルを仮定し、これに初期値を与え、周波数ｆ、磁場強度ＨＰを基に地表での２次磁場Ｈｓの理論値を計算すれば、実際に２次磁場Ｈｓを測定することで、その差が判明し、この差がなくなるように地下モデルの修正を行い、これを繰り返せば目的とする比抵抗分布を逆解析することができる。
【００４２】
さて、再び図１に戻って、実施の形態１における電磁探査装置の構成について説明を続ける。６は電磁探査装置の制御部である。制御部６は送信コイル２に供給する電流ｉ＝ｉ０ｅｊωｔを制御し、受信コイル４で受信した磁場情報（ＨＳ＋ΔＨＰ）≒Ｈｓから２次磁場Ｈｓの同相成分Ｉ（実数部）及び離相成分Ｑ（虚数部）の測定値を求めることができる。周波数を様々変更して測定することもできる。
【００４３】
なお、地盤の比抵抗ρは（数２）によって演算する。（数１）を変形したものである。
【００４４】
【数２】
2009074953.tif 000004
１次磁場ＨＰと２次磁場Ｈｓとの関係、ΔＨＰ＝（ＨＰ−Ｈ’’Ｐ）の関係を示したのが図３である。合成磁場Ｈは１次磁場ＨＰと交角θであり、１次磁場ＨＰと２次磁場Ｈｓの位相差（π／２＋φ）のφは（数３）となる。これからも分かるように−ΔＨＰが地下モデルの比抵抗分布の正確さに大きく影響することが分かる。
【００４５】
【数３】
2009074953.tif 000005
さて、実施の形態１の制御部６は、１次磁界ＨＰを形成するために送信コイル２の電流Ｉを制御し、磁場検出部３で検出されたＨＳの同相成分Ｉ及び離相成分Ｑを取り出し、処理結果を分析装置７に渡す。
【００４６】
この分析装置７は、実施の形態１の場合はパーソナルコンピュータであって、ＣＰＵ（中央処理演算装置）７ａと記憶装置７ｂを備えている。記憶装置７ｂには、ＲＯＭ（Read only Memory）と一時記憶のためのＲＡＭ（Random
access Memory）、さらに実施の形態１のプログラムを格納するための不揮発性メモリが設けられている。ＣＰＵ７ａはＲＯＭや不揮発性メモリに格納された各機能を実行するプログラムを読み出して、それぞれの機能を実行する機能実現手段として機能する。
【００４７】
従って、実施の形態１の分析装置７のＣＰＵ７ａには、後で詳述する手順で磁場の計算を行うと共に、所定の比抵抗分布を有する地下モデルを逐次更新しながら目的の比抵抗分布を有する地下モデルにまで到達させるプログラムがデータと共に読み込まれ、機能実現手段（本発明の逆解析手段）として地盤を逆解析する機能を実行する。なお、このプログラムを外付けの記憶媒体に格納しておくのもよい。
【００４８】
以下、実施の形態１の電磁探査装置、電磁探査方法で行われる処理の詳細について説明する。地下モデルは地盤を複数のブロックで分割して構成する。各ブロックの比抵抗ρの組をベクトル
2009074953.tif 000006
（外字１、以下ｍ）とし、電磁探査データ（実測値）の組をベクトル
2009074953.tif 000007
（外字２、以下ｄ）と表す。地下モデルはマクスウェルの方程式を差分または有限要素法で解くものとする。
【００４９】
さらに、各ブロックにおける比抵抗ρの組ｍを入力すると、地下モデルに対応した電磁探査における応答値（各ブロックが形成する２次磁場の総和としてのモデル計算値）を出力する関数を
2009074953.tif 000008
（外字３、以下ｆ（ｍ））とする。
【００５０】
すなわち、ブロック単位でそれぞれの比抵抗ρをもつ地下モデルに対して、地表において、この比抵抗分布のｍに対応した２次磁場（モデル計算値）を生成する関数である。電磁探査データは、地表上の二位置の高さあるいは三位置以上の高さで実測した電磁探査データｄを統合して構成される。このとき、電磁探査データとの差を示すオフセット量（システムのオフセット誤差、その他の誤差）をベクトル
2009074953.tif 000009
（外字４、以下ｓ）で表されるとすると、（数４）の関係になる。
【００５１】
【数４】
2009074953.tif 000010
ここで、
2009074953.tif 000011
（外字５、以下Ｇ）はオフセット量ｓを探査データｄの各要素に関係付けるマトリックスであり、各列の適当な位置だけが１であり、他の要素はすべて０の行列である。ｓは比抵抗ρの組ｍと同様に求めるべき未知数のベクトルである。
【００５２】
さて、（数４）を直接解くことはできないので、最小二乗法でこれを解くことになる。すなわち、本発明の実施の形態１における電磁探査方法は、地下の比抵抗分布のモデルを（数５）に示す目的関数φを最小化する比抵抗ρの組ｍとｓを決定する問題に帰着する。ここで
2009074953.tif 000012
（外字６、以下ノルム）はベクトルのノルム（要素の二乗和）であり、
2009074953.tif 000013
（外字７、以下Ｃ）はラプラス（二次微分）線形操作を表している。
【００５３】
【数５】
2009074953.tif 000014
（数５）の第１項は探査データとモデル計算値及び誤差との差の二乗和を表し、第２項は比抵抗分布の粗さを表す。このλはラグランジェの未定係数（本発明の調整係数）であって、本来比抵抗分布の粗さは最小の０でなければならないという条件付きにしたものである。この第２項を加えたことによって第１項の残差と第２項のバランスが調整され、探査データの地下モデルへの当てはまりが最小化されるだけでなく、オフセット量が明らかになり、比抵抗分布が滑らかに変化する地下モデルを求めることができる。ただ、λについても最小値（データ残差）を与えるλを求めるように同時計算するのがよいが、計算時間がかかる上に、経験上λを固定した場合と大差ないので、複数の、例えば３値の中で最小値を与えるλを選択する方法を選んだり、経験上良好な結果を与える所定値λ０に固定した値を使用したりするのでもよい。
【００５４】
なお、実施の形態１における電磁探査装置、電磁探査方法で実測した電磁探査データは、地表で測定した１回目の電磁探査データと、地表面から１ｍ〜２ｍの高さで測定した２回目の電磁探査データとを統合したデータである。上述したように、更に３回以上の探査データまで統合した電磁探査データであってもよい。オフセット量ｓは１回目，２回目・・の電磁探査データ間で独立した値である。
【００５５】
さて、ｆ（ｍ）はｍに関して非線形であるため、初期値
2009074953.tif 000015
（外字８、以下ｍｋ−１）を与えて、初期値周りでｆ（ｍ）のテイラー展開を行って線形化する。すると、（数５）は、初期値ｍｋ−１からの変化分
2009074953.tif 000016
（外字９、以下Δｍｋ）、すなわち、（数６）を最小にするΔｍｋとオフセット量ｓを求めることに帰着する。オフセット量ｓはランダム誤差を除いた誤差を吸収する。
【００５６】
【数６】
2009074953.tif 000017
ここで、
2009074953.tif 000018
（外字１０、以下Ｊ）はヤコビアン（微係数行列）である。ここで、
2009074953.tif 000019
（外字１１、以下Δｄ）は（数７）の関係を有する。
【００５７】
【数７】
2009074953.tif 000020
以上を言い換えると、（数６）を最小化することは次の観測方程式（数８）の最小二乗解を求めることと等価である。
【００５８】
【数８】
2009074953.tif 000021
（数８）を簡潔な形に形式的に書き換えると（数９）となる。ｎ×ｍ行列
2009074953.tif 000022
（外字１２、以下Ａ）は修正グラム・シュミット法やハウスホルダー法により（数１０）のように直交分解、ここではＱＲ分解が行える。
【００５９】
【数９】
2009074953.tif 000023
【００６０】
【数１０】
2009074953.tif 000024
ここで、
2009074953.tif 000025
（外字１３、以下Ｑ）は列ベクトルが正規直交するｎ×ｍ行列（すなわち、ＱＴＱ＝Ｉ、ここでＱＴはＱの転置行列）であり、
2009074953.tif 000026
（外字１４、以下Ｒ）は列ベクトルが正規直交するｍ×ｍ行列の右上三角行列である。言い換えれば、ＱＲ分解できると言うことは、このようなＱとＲが存在することを意味する。このときの最小二乗解は（数１１）となる。
【００６１】
【数１１】
2009074953.tif 000027
この（数１１）は後退代入していくだけで解くことができる。これによって、初期値ｍｋ−１からの変化分Δｍｋが求められ、（数１２）によって地下モデルを逐次修正して、地下モデルの更新を行いながら、ｍとオフセット量ｓが所定の値の範囲内にそれぞれ収斂するようにする。これにより変化分Δｍｋが０に近づく。あるいは所定回数修正を繰り返してもよい。このような処理によって、比抵抗ρの組ｍから構成された目的とする地下モデルと、オフセット量ｓが求まる。
【００６２】
【数１２】
2009074953.tif 000028
このように実施の形態１の電磁探査装置、電磁探査方法によれば、従来のように高高度でのレベリングを行ってオフセット誤差（０レベルの設定）を調整することなく、１次磁場の影響として残留する誤差や測定位置、環境に由来するオフセット量を一括して処理でき、良好な比抵抗分布をもつ地下モデルを得ることができる。比抵抗分布と共にオフセット量を計算するので、高い場所でのレベリング作業が不要になり、測定段階でオフセット誤差を考慮する必要がなくなる。
【００６３】
次に、図４に基づいて実施の形態１のプログラムの手順について説明する。図４は本発明の実施の形態１におけるプログラムのフローチャートである。まず、実測の探査データｄとモデル計算値を出力するｆ（ｍ）を読み出す（ｓｔｅｐ１）。次いで、比抵抗の組ｍｋの初期値ｍ０とラグランジェ係数λの初期値λ０（固定したλの場合は固定値）を設定する（ｓｔｅｐ２）。
【００６４】
続いて、ｍｋ−１＝ｍ０として、Δｄ＝ｄ−ｆ（ｍｋ−１）を計算し（ｓｔｅｐ３）、比抵抗分布のラプラシアンＣｍｋ−１を計算し、これを基に−λＣｍｋ−１を計算する（ｓｔｅｐ４）。さらに、ｆ（ｍｋ）のヤコビアンＪを計算してｎ×ｍ行列Ａを計算する（ｓｔｅｐ５）。
【００６５】
次に、得られた行列ＡをＱＲ分解して（ｓｔｅｐ６）、Δｍｋ，ｓｋ，λｋ（固定したλの場合は除く）を求める（ｓｔｅｐ７）。最小二乗解を与えるΔｍｋ，ｓｋ，λｋ（固定したλの場合は除く）を用いて、Δｄ−ＪΔｍｋ−１−Ｇｓのノルムが所定量εより小さくなって収斂したか否かを判定する（ｓｔｅｐ８）。
【００６６】
ｓｔｅｐ８において収斂していない場合は、ｍｋ−１＋Δｍｋを計算して（ｓｔｅｐ９）、ｓｔｅｐ３のｍｋ−１に代えてｍｋを代入する（ｓｔｅｐ１０）。ｓｔｅｐ８において、Δｄ−ＪΔｍｋ−１−Ｇｓのノルムが所定量εより小さくなって収斂した場合、計算を終了する。
【００６７】
このように実施の形態１の電磁探査方法を実行するプログラムによれば、高高度でのレベリングを行ってオフセット誤差（０レベルの設定）を調整することなく、１次磁場の影響として残留する誤差や測定位置、環境に由来するオフセット量を容易且つ高精度に処理でき、良好な比抵抗分布をもつ地下モデルを得ることができる。パーソナルコンピュータ等のコンピュータの記憶装置に実施の形態１のプログラムを格納しておき、上述した手順で、比抵抗分布と共にオフセット量を計算すすればよい。従って、高い場所でのレベリング作業が不要になり、予め（測定段階で）オフセット誤差を補正する必要がなくなる。
【実施例】
【００６８】
（実施例１）
本発明の実施例の作用効果をみるために、図５のような２次元の地下モデルで計算によるモデル実験を行った。図５は本発明の実施例１におけるモデル実験の構造の説明図、図６（ａ）は比較例としてのモデル実験のオフセット誤差がない場合の説明図、図６（ｂ）は比較例としてのモデル実験のオフセット誤差がある場合の説明図、７（ａ）は本発明の実施例１におけるモデル実験の初期値１００Ωｍの場合の説明図、図７（ｂ）は本発明の実施例１におけるモデル実験の初期値１０Ωｍの場合の説明図である。
【００６９】
図５に示すように、地下モデルには比抵抗ρが２Ωｍのブロックで形成された階段状になった山の部分が３つある。深さ４ｍから２ｍまで膨出し、これが再び６ｍの深さにまで凹み、再度深さ２ｍにまで膨出し、これをもう一度繰り返した連山状の分布である。地表側はこれを裏返した形状で比抵抗ρが２０Ωｍとした。
【００７０】
図６（ａ）はシステム上のオフセット誤差がない（補正済み）とした場合の理想的な電磁探査装置、電磁探査方法による解析結果である。等高線は比抵抗ρが一定値を示す。比較的正しく比抵抗ρが２Ωｍの山を検出している。しかし、２つの山（膨出部分）しか検出できていない。そして、図６（ｂ）はオフセット誤差を残し、これを考慮していない従来の電磁探査装置による解析結果である。実験モデルとはまったく別の地下モデルになっている。逆に言えば従来の電磁探査装置でオフセット誤差等を補正しなければ、この程度の結果しか得られない、と言うことである。
【００７１】
これに対して、図７（ａ）（ｂ）は本発明のオフセット量を考慮して処理する場合の解析結果である。図７（ａ）の初期値は１００Ωｍである。この結果をみれば、比抵抗ρが２Ωｍの山が３つほぼ正確な位置で捉えられている。また、図７（ｂ）の初期値は１０Ωｍである。この結果も、図７（ａ）と同様に、比抵抗ρが２Ωｍの山が３つほぼ正確な位置で捉えられている。これらの解析結果によれば、本発明のオフセット量を考慮した処理が、オフセット誤差がない（補正済みの）場合よりも、はるかに正確に地盤構造を捉えていることが分かる。なお、データには標準偏差が１００ｐｐｍの正規分布誤差（ノイズ）を加えて計算を行った。
【００７２】
本発明の電磁探査方法で、実際に現地に赴いて電磁探査を行った解析結果を示す。探査は韓国Ｙｏｅｎｋｗａｎｇで行われた。図８（ａ）は本発明の実施例１における現地での電磁探査を初期値１００Ωｍで計算した場合の解析結果の説明図である。これによればなだらかな地層が形成されているのが分かる。また、図８（ｂ）は本発明の実施例１における現地での電磁探査を初期値２０Ωｍで計算した場合の解析結果の説明図である。図８（ａ）とほとんど同様の解析結果を示している。
【００７３】
図９は比較例としての電気探査方法で本発明の実施例１の電磁探査と同時に行った解析結果の説明図である。電気探査法は浅い地盤で正確だと言われているが、本実験では、実施例１の電磁探査と同時に、この電気探査による解析を行った。その結果を示している。電気探査法によって解析された地盤は、本発明の電磁探査方法で解析した地盤に似ているが、実施例１と比べて地層の形態がそれほどよくない。
【００７４】
図１０（ａ）はモデル実験で計算を繰り返した回数とデータ残差の関係を示した説明図である。データ残差とは実測値と計算値の差のことを意味する。すなわち、この１回の計算でそれぞれ比抵抗分布とオフセット量が同時に得られたことを意味している。横軸は反復回数であり、縦軸がデータ残差である。これによれば、ケース１は初期値を１００Ωｍにした場合、４回の反復計算で標準偏差誤差のほぼ１００ｐｐｍに収斂しており、ケース２は初期値を１０Ωｍにした場合は３回の計算だけで標準偏差誤差のほぼ１００ｐｐｍに収斂している。
【００７５】
同様に、図１０（ｂ）は韓国Ｙｏｅｎｋｗａｎｇで測定したデータを使って繰り返して計算した回数とデータ残差の関係を示した説明図である。同じく、横軸は反復回数であり、縦軸がデータ残差である。ケース１は初期値を１００Ωｍとしたものであるが、４回ほどの繰り返し計算によってほぼ１００ｐｐｍに収斂している。また、初期値を１０Ωｍとした場合は２回の繰り返し計算によってほぼ１００ｐｐｍに収斂した。なお、この１００ｐｐｍは上記標準偏差誤差とは異なるが、この実験の場合一致した値に収斂した。
【００７６】
さらに、本発明の電磁探査方法による実験結果と電気探査法による実験結果の比較を行ったものを図１１として示す。図１１は本発明の実施例１における電磁探査方法の実験結果と電気探査法による実験結果の比較を行った比較図である。横軸は本発明の電磁探査から求められた比抵抗であり、縦軸は電気探査から求められた比抵抗である。対数表示している。比較的良い相関をしている。
【００７７】
従って、本発明の電磁探査方法によって地中探査を行えば、１０ｍ以上の高度でレベリングを行ってオフセット誤差を補正することなく、たった２回〜４回の繰り返し計算により簡単にオフセット量を高精度に処理できる。このため短時間での地下探査、低コストの地下探査が可能となる。演算時間が少ない上に、きわめて良好な比抵抗分布を推定できる。比抵抗分布と共にオフセット量を計算するので、測定段階では（予め）オフセット誤差を考慮する必要がなくなる。
【産業上の利用可能性】
【００７８】
本発明は、電磁誘導作用を利用して地中探査を行う電磁探査装置、コンピュータに電磁探査の比抵抗分布を推定させるためのプログラムに適用できる。
【図面の簡単な説明】
【００７９】
【図１】本発明の実施の形態１における電磁探査装置のブロック構造図
【図２】測定の原理の説明図
【図３】１次磁場と２次磁場、及びキャンセルのための磁場とオフセット量との関係を示す説明図
【図４】本発明の実施の形態１におけるプログラムのフローチャート
【図５】本発明の実施例１におけるモデル実験の構造の説明図
【図６】（ａ）比較例としてのモデル実験のオフセット誤差がない場合の説明図、（ｂ）比較例としてのモデル実験のオフセット誤差がある場合の説明図
【図７】（ａ）本発明の実施例１におけるモデル実験の初期値１００Ωｍの場合の説明図、（ｂ）本発明の実施例１におけるモデル実験の初期値１０Ωｍの場合の説明図
【図８】（ａ）本発明の実施例１における現地での電磁探査を初期値１００Ωｍで計算した場合の解析結果の説明図、（ｂ）本発明の実施例１における現地での電磁探査を初期値２０Ωｍで計算した場合の解析結果の説明図
【図９】比較例としての電気探査方法で本発明の実施例１の電磁探査と同時に行った解析結果の説明図
【図１０】（ａ）モデル実験で計算を繰り返した回数とデータ残差の関係を示した説明図、（ｂ）韓国Ｙｏｅｎｋｗａｎｇで測定したデータを使って繰り返して計算した回数とデータ残差の関係を示した説明図
【図１１】本発明の実施例１における電磁探査方法の実験結果と電気探査法による実験結果の比較を行った比較図
【符号の説明】
【００８０】
１    磁場発生部
２    送信コイル
３    磁場検出部
４    受信コイル
５    バッキングコイル
６    制御部
７    分析装置
７ａ    ＣＰＵ
７ｂ    記憶装置
【図１】
2009074953.tif 000029
【図２】
2009074953.tif 000030
【図３】
2009074953.tif 000031
【図４】
2009074953.tif 000032
【図５】
2009074953.tif 000033
【図６】
2009074953.tif 000034
【図７】
2009074953.tif 000035
【図８】
2009074953.tif 000036
【図９】
2009074953.tif 000037
【図１０】
2009074953.tif 000038
【図１１】
2009074953.tif 000039