第 17 回 MAGDA コンファレンス in 日立、2008.11.20~21P47
電流分布の赤外線による可視化とその応用
Visualization of Current Distributions by Infrared Imaging and Its Application-赤外線による過渡熱分布測定-Transient Temperature Distribution Measurement by Infrared Imaging鈴木剛*1(学生員),齋藤

兆古*2（正員）

Goh SUZUKI (Stu.), Yoshifuru SAITO (Mem.)
In order to realize the reliable nondestructive inspection systems, many infrared image sensor systems aredeveloped and utilized because of their high visual capability. Recently, nuclear electric power plant is increasing its
usefulness not warm up earth atmosphere accompanying with human life activities. We are now developing toenhance an infrared inspecting system to carry out the condition based maintenance methodology mainly for nuclear
electric power plant use. At first, this paper clarifies a relationship between the pixel value constructing infrared image
and absolute temperature even though automatic gain control function equipped in infrared CCD camera is activating
to get the highest contrast image. Second we demonstrate the absolute temperature visualization when heating thetargets by magnetic induction. This makes it possible to clarify the versatile capability of our scheme for various
condition based maintenance applications.
Keywords: Infrared image, Nondestructive testing, induction heating

１ 緒論

２ 基準温度同時撮影法

映像情報による診断は,医学の世界などでも用いら本実験において使用した赤外線カメラは「三菱サれているように極めて具体的な状態把握を可能にする.

ーマルイメージャ（形式 IR-SC１

たとえ言語の異なる民族間でも,映像情報を用いるこ

社）
」である.

三菱電機株式会とで意思疎通が可能となる.このように人間の視覚情

赤外線画像は対象物の温度分布をモノクロ濃淡情

報処理は全地球人類共通の強力な知的機能である.人報として表示するものであり,温度の高い部分を白色,

間は外部から得る情報のなかで 80 パーセント以上を温度が低い部分を黒色で可視化する．温度の高低に伴視覚から取り入れている[1].

う白黒濃淡情報の範囲（ダイナミックレンジ）は有限

近年,商品の品質検査・管理のために,赤外線カメラをであるため,対象の温度情報を高コントラストで取得用いた非破壊検査技術が広範に用いられている[2-4]．

する方途として,AGC(Automatic Gain Control, 対象のこの背景として,赤外線カメラの広汎な普及と監視技

温度分布で最も支配的な温度を中心としたダイナミッ術の高度化が挙げられる.

クレンジへ自動制御する機能)が組み込まれている．
こ

本研究は,第 2 世代赤外線非破壊検査として,コンピれは相対的な温度分布可視化に有効であるが,絶対的ュータを前提とする画像処理技術,すなわち,画像認

な温度分布可視化を不可能にする．可動部分や動力源

識・識別・監視技術を含めた,赤外線画像情報による知を含むシステムにおいては,相対的な温度分布よりも的非破壊検査・監視技術の開発を目指すものである.

遙かに絶対的温度分布が有意義である場合もある.こ

本論文では,従来から我々が提唱する“基準温度同時撮影法”を用いて検査対象である金属の過渡温度上昇

れは,多くの機器が特定の許容絶対温度内で正常に機能すべく設計されていることから自明である．

分布動画像を赤外線ビデオカメラの熱画像から抽出し,以上の事から,単純な赤外線可視化画像から各種機

金属中の欠損やシステム欠陥箇所,また非金属体へ混器の適正動作温度を前提とする高精度連続監視シス

入した金属片可視化法を提案する.

テム構築は困難である. しかし,AGC は広範囲な温度幅

連絡先： 鈴木 剛，〒184-8584 東京都小金井市梶野町
3-7-2 ，法政大学 大学院工学研究科 電気工学専攻 齋
藤兆古研究室e-mail: suzuki@ysaitoh.k.hosei.ac.jpの撮影を可能とし,赤外線カメラの適用範囲をより広*1

法政大学大学院*2

法政大学大学院

汎なものとする.換言すれば AGC が備わっていない場
合,温度感知範囲が限定されるため,不測の温度に対する熱画像が得られない．
-479-

第 17 回 MAGDA コンファレンス in 日立、2008.11.20~21

Fig. 2 Temperature versus pixel value
Fig.1 Infrared cup images containing different
temperature water.

的情報を求める航空写真に拠る地図作成法と同様であ赤外線 CCD に拠る熱画像は対象の赤外線反射率に依既知である対象の赤外線画像を構成する画素値と温度る．このため,“基準温度同時撮影法”では,温度が存する．
すなわち,対象の赤外線反射率が厳密に既知での関係を把握しなければならない．
多くの自然科学の問題では,温度上昇などの拡散現無い限り厳密な熱画像・温度画像は得られない．
AGC と赤外線反射率問題を解決する一方法として,

象は時間・空間に対しては指数関数的に変化する.この

本論文では,赤外線画像中に既知の異なる温度を持つ

ため,ここでは赤外線画像を構成する画素値と絶対温

複数の温度基準画像を被温度測定対象と同時に撮影し,度の関係も指数関数的に関係すると仮定し, Fig.2 に示

既知温度を有する複数個の画像を構成するそれぞれのすように両者を両対数グラフへプロットする．Fig.2 か

画素値とそれらの温度間の関係を用いて被温度測定対

ら画素値と温度が比例するため,明らかに両者の関係

象の絶対温度を測定する“基準温度同時撮影法”を提は指数関数的であり,Fig.2 を用いて画素値を温度へ換唱する[2,3,4].さらに,この基準温度同時撮影法と誘算すれば絶対温度分布の可視化画像を得ることが可能

導加熱を併用した赤外線非破壊検査の幾つかの例を示

となる.

す．

３ 誘導加熱に拠る非破壊検査

Fig. 1 は異なる温度の湯を満たしたコップの赤外線画像である.それぞれのお湯の温度は,上段の左から46．4℃,29．0℃,25．2℃,下段の左から 15．7℃,
12.8℃である．赤外線カメラの ACG 機能が有効に機能し,高温から低温へ至るカップのお湯を最大コントラ
ストで可視化している．お湯の表面画像を構成する画
素値は,高温から順に,180,101,87,52,43 となり,絶対Fig. 3 Two Iron plates to simulate a straight
crack.

温度と画素値は比例しないことが判る．
従って,赤外線
画像から対象物の相対的な温度分布は可視化可能であるが,高精度な絶対温度分布の可視化は困難である．
通常のデジタルカメラ用 CCD であっても,赤外線波長領域に於ける光エネルギーの入射によって電荷移動が励起されることから,赤外線 CCD と同じ機能を有する．
このため,仮に通常の CCD を用いて正確な絶対温度画像

Fig. 4 Thermal Image of two iron plates when
heating by alternating magnetic fields.

が取得可能であれば,遙かに低コストで温度センシン
グシステムの構築が可能となる．
“基準温度同時撮影法”の基本的な着想は,長さな3.1 鉄板中の亀裂探査

ど幾何学的情報が既知の対象を画像中に同時に写し込

鉄板中のクラックを模擬するため 2 枚の鉄板を接合

み,全体の画像中で,幾何学的情報が既知の対象画像をした実験モデルを Fig.3 に示す．平面状スパイラル励

基準として幾何学的情報が未知である対象物の幾何学

磁コイルへ通電し Fig.3 のモデルを誘導加熱した．定常的な加熱状態を表す赤外線画像を Fig.4 に示す．
Fig.5

-480-

第 17 回 MAGDA コンファレンス in 日立、2008.11.20~21

Fig. 7 Exciting coil used for induction heating

Fig.5 Absolute temperature distribution of two ironplates when heating by alternating magnetic fields.
は基準温度同時撮影法によって得られた絶対温度分布である．赤部分が最も高温部を示し,白色に近づくにつれ低温部を示す．この結果から鉄板接合部分は周囲に比べて摂氏 30～40 度高温であることがわかる.すなわ
ち，加熱電流である渦電流密度は鉄板の接合部で最大となる．
3.2 非金属材料中の金属探査

Fig. 8 Absolute temperature rise distribution of thecloth folding in an iron plate
Fig.9 はノイズの状況を把握するために Fig.8 の温度上昇を縦軸に取った Fig.8 の 3 次元表示である．
Fig.8, 9 より金属片が位置する部分の温度は,誘導加熱前と比較して約 3℃温度上昇することが確認できる.
すなわち,金属片は布の間に挟まれているため目視が不可能であるが,誘導加熱によって金属片が加熱され，結果として,非金属中に混入した金属が可視化される.

Fig. 6 A cloth folding in a 1.0×1.0cm2 iron plate
衣料品中の縫製針や魚介類中の金属網破片探査は
直接人体の内外へ障害を与えるため,極めて精密で正確な検査を必要とする．
Fig.6 に示すように非金属中に混入した金属片を模擬し,二つ折りにした布間に 1.0×1.0cm2 の鋼板を挟んだ．Fig.7 に示す励磁コイルへ通電して生成される交流磁界中へ Fig.6 の金属片を挟んだ布を配置する．誘導
加熱開始２分後に赤外線画像を撮影した．得られた赤Fig. 9外線画像へ基準温度同時撮影法を適用して絶対温度上

3D representation of the absolute temperature rise
distribution

昇分布を可視化した．Fig.8 は温度上昇画像を示す．
次に,加熱対象となる金属と,温度上昇の関係について吟味する．Fig. 10 はそれぞれ 0.5×0.5, 1.0×1.0, 1.5-481-

第 17 回 MAGDA コンファレンス in 日立、2008.11.20~21

1) 同じ条件下で誘導加熱を行った場合,励磁コイル
面と平行な加熱対象金属の表面積が大きいほど発熱量が大きい,2) 最大温度上昇値から金属片の励磁コイルに平行な面の相対的な寸法の概略が把握できる.
４ 結言
Fig. 10Six tested metallic plates having different
surface area.

本論文では,赤外線 CCD カメラを用いて絶対温度測定を可能とする一方法として“基準温度同時撮影法”を提唱し,その有効性を検証した.その結果,赤外線カメラを用いた絶対温度状態監視システム構築の一端が
開けた．

参考文献
[1] 中島 龍興,近田 玲子,面出 薫:,”照明デザイン入門”, 彰国社(1955)[2] 鈴木剛, 齋藤兆古，堀井清之, “過渡熱伝達分布測定によFig. 11 Temperature rise versus surface area of thetested metallic plates
×1.5,2.0×2.0, 2.5×2.5, 3.0×3.0cm2 の 6 種類の金属片を示す．これらの金属片を 3.1 や 3.2 節と同様にそれぞる金属欠損の可視化”, 第35回可視化情報シンポジウム,
可視化情報 C211, Vol.27, Suppl. No.1 (2007年7月)[3] 鈴木 剛，齋藤 兆古, 堀井 清之,”:電気抵抗の赤外線画像
による可視化”, 第36回可視化情報シンポジウム, 可視化

れを布に挟んだ状態で誘導加熱を行い,絶対温度上昇
分布画像を基準温度同時撮影法によって得た．
6 種類の金属片で,それぞれの最高温度上昇は 1.76℃,

情 P01-001, Vol. 28, Suppl. No.1 (2008年7月）[4] 鈴木 剛, 齋藤 兆古, “過渡熱伝達分布測定による金属欠2.95℃, 5.66℃, 6.39℃, 7.59℃, 10.63℃であった．金属片の表面積と,温度上昇の関係を Fig. 11 に示す．Fig.11 の結果から以下のことがわかる．

-482-

損の可視化”, 第16 回MAGDA コンフェレンス in 京都,(2007.11.22-23)

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-004

セルオートマトンによる磁化特性の可視化‐複雑系の物理とプライザッハモデル‐
宮坂総 ○，齋藤兆古(法政大学大学院)，堀井清之(白百合女子大学)

Visualization of Magnetization Characteristics by Cellular Automaton- Physics of complex systems and Preisach model Sou MIYASAKA, Yoshifuru SAITO, Kiyoshi HORII
ABSTRACTThe cellular automaton is being widely used for analyzing complex physical systems
as a deterministic methodology. Although theory of cellar automaton innovated in 1940
by Neumann, German scientist Prisach had worked out the same to the cellar
automaton model to represent the nonlinear magnetization characteristics offerromagnetic materials in 1935. In 1984, S. Wolfram pointed out that method of cellar
automaton was a reasonable methodology to represent the complex physical systems but
Preisach had derived as is called Preisach model by cellar automaton for typicalnonlinear complex ferromagnetic physical system.
In this paper, we try to derive the Preisach type, i.e. cellar automaton, model byvisualizing magnetic domain movements by Bitter method.
Keywords: Cellular Automaton, Preisach model, Hysteresis
１．緒

論

値を持つ不連続な系からなり、これらのサイトはそれぞれのサイトがいくつかの局所的な隣接サイトの値に基づ

セルオートマトンは、1940 年代にノイマンによって提いた新しい値と有限な数の過去の時間ステップを仮定す唱された．その後、1984 年にウルフラムによって複雑系

ると、離散的な時間ステップで状態が変化するものであのシミュレーション行う決定的な方法として紹介され、る。2)

近年、多くの分野で新世代シミュレーション技術として

また、セルオートマトンには、空間の広がりに応じて研究・開発が行われている．特に、従来の方法ではモデ1 次元、2 次元、さらに 3 次元のセルオートマトンがある．

ル化が困難であった交通渋滞、雪崩、さらに火事などの
自然災害のシミュレーションにも用いられており、今後さらなる研究開発と実用化が期待されている．

2.2

セルオートマトンの例

ある一定のルールを決め、実行した 1 次元のセルオー

一方、磁性体は極めて複雑な磁化特性を呈するが、1930

トマトン例を示す。Fig.1 は初期値(0 1 0 0 1 1 0 1 0 1

年代にこの磁化特性を表現するモデルがプライザッハに

1 0 0 0 1 0)でセルオートマトンのルール 30((左近傍 tよって提唱された．このプライザッハモデルは極めて具の値 右近傍)→t+1 の値：(1 1 1)→0,(1 1 0) →0,(1 0体的な磁区理論に基づくため、他に代替えとなるモデル

1) →0,(1 0 0) →1,(0 1 1) →1,(0 1 0) →1,(0 0 1) →が無く、現代でも磁気記録理論や磁性体を含む磁界計算

1,(0 0 0) →0)を適用した場合を Fig.1 に示す．

に広範に使われている．

体的には、ビッター法を用いた可視化磁区動画像からセ

t=0
t=1
t=2
t=3
t=4
t=5

ルオートマトンモデル導き、このセルオートマトンモデ

Fig.1 Example of Cellar Automaton based on Rule 30.

本稿は、セルオートマトンを用いてプライザッハモデ
ルを導き、プライザッハモデルはセルオートマトン理論が提唱される以前に既に提案され実用化されたセルオー
トマトンによる磁化特性モデルであることを述べる。具ルから強磁性体の磁化特性を導く．1)

３．磁化特性のセルオートマトンモデル
２．セルオートマトン
3.1 プライザッハモデル
2.1

セルオートマトンとは

セルオートマトンとは、格子サイトがいろいろな初期

セルオートマトンを用いたプライザッハモデルのシミュレーションを行った結果を Fig.2 に示す．

1

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-0043)プライザッハ線図からレイリーループを描き、4)ビッmagnetic flux density

ター法によって可視化された磁区画像を、磁化過程の相違で 3 領域へ分類し、5)分類された領域の中で可逆的磁壁移動範囲に対する磁区画像を画素値に応じて１かゼロ
へ 2 値化してセルオートマトンモデルを作成した．6)磁区画像から得られたセルオートマトンモデルから磁化特性へ対応する線図を得ることに成功した．
100

magnetic

Fig.2

field

Magnetization Curve computed by Cellar
Automaton, i.e. Preisach Model.

200

300

100

200

300

50

50

50

50

100

100

100

100

150

150

150

150

200

200

200

Fig.2 のシミュレーションに使ったセルオートマトン

100

200

200

100

300

0.00[s]200300

0.167[s]を Fig.3 に示す． Fig.3 の図の黒色部分は１の値を、灰
100

200

100

300

200

300

色部分は−１の値を持つ。Fig.2 はそれぞれの画像の値の総和を各フレームの値とし黒に転移する場合はプラスに磁化され灰色に転移する場合はマイナスに磁化されると考える．Fig.3 はプライザッハが提案したプライザッ

50

50

50

50

100

100

100

100

150

150

150

150

200

200

200

ハ線図そのものに他ならない．

100

Frame21

Frame30
Frame41
Frame50
Fig.3 Cellar Automaton Images representing theMagnetization Processes.

得られた磁区挙動画像を、磁化過程によって 3 過程に分類した．磁性体の磁化過程は、可逆的磁壁移動範囲、非可逆的磁壁移動範囲さらに可逆的回転磁化領域分けられる．簡単のため、ここでは可逆的磁壁移動範囲について

300

0.50[s]Binary Images of Magnetic Domains.

field magnetic

Fig.5 Diagram corresponding to the Magnetization
Characteristics.

3.2 磁区挙動からの解析
磁性体表面磁区挙動をビッター法によって可視化した．

200

density flux magnetic

Frame10

200

100

300

0.333[s]Fig.4

Frame1200

参 考 文 献1) 須永高志、齋藤兆古、堀井清之：ビッター法による磁性材
料の鉄損分布の可視化、可視化情報学会シンポジウム(2007)

C207

2) R.J.ゲイロード/P.R.ウェリン、荒井隆＝訳：Mathematica
複雑系のシミュレーション物理学と生物学の探求

考える．可逆的磁壁移動範囲で、正に単位磁化された部

3) R.J.ゲイロード/西舘数芽＝著、西舘数芽/西川清＝訳：

分を白色（数値１）、磁化されてない部分を黒色（数値ゼMathematica 自然現象の計算モデル化セルラーオートマロ）として 2 値化した磁区画像 Fig.4 は、Fig.3 の白色部がランダムに存在する、セルオートマトン画像に対応

タ・シミュレーション
4) 西成活裕：セルラー・オートマトンによる複雑現象のモデ
ル化、東京大学ホームページ、する．Fig.4 から Fig.3 と同様にして Fig.5 の磁化特性

http://soliton.t.u-tokyo.ac.jp/nishilab/mypapers/rikに対応する線図が得られる．

ouJ.pdf

4. まとめ
本稿では、１）簡単なセルオートマトンの例を述べ、2）2 次元セルオートマトンとしてプライザッハ線図を示し、

2

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006

1/f
齊藤

ゆらぎの可視化

兆古 ○ (法政大学），宮坂総(法政大学）、加藤千恵子(東洋大学）

Visualization of 1/f Fluctuation

Yoshifuru SAITO, Soh MIYASAKA and Chieko KATO
ABSTRACT
1/f fluctuation analysis of the non-life-bearing, life-bearing and self driven particle
targets is carried out. As a result, it is clarified the followings. State transition, e.g., melting an
ice to water, of all the non-life-bearing targets exhibits 1/f fluctuation frequency characteristics.
Namely, all the non-life-bearing materials exhibit 1/f fluctuations when changing their state
aspect such as solid, liquid and gas. On the other side, 1/f fluctuations can be observed thecreative works, e.g., artistic paints when transitioning mental situation from depression tomania of life-bearing target, i.e., his or her. In the other words, a life-bearing target typicallyhuman exhibits 1/f fluctuation characteristics via his own creative work when changing his
mental situation.

Keywords: Visualization, 1/f Fluctuation , Healing, Artistic Paint Works
１．緒論

約 30 年近く前になるであろうか、毎日々計算機で非線形磁化特性問題を解析していた．この問題は依然として想される．これは「渋滞学」で使われる「自己駆動粒子（Self-driven particle）」モデルが磁性体を構成する磁区モ
デルに対応し、Cellar Automaton モデルそのものが異方性を前提としているためである 3)．

計算電磁気学で解決できていない問題であり、汎用電磁非線形系解析の理論的背景は 21 世紀に入ってかなり

界解析パッケージが持つ最大のウイークポイントである．

体系的に系統化され、その結果、非線形系に特有な「1/f1981 年 9 月に米国シカゴで開催された最初の IEEE 主催

ゆらぎ」は化学反応系の相転移（Phase Transition）時に

COMPUMAG（計算電磁気学に関する国際会議）で自分生ずることが解明された 4)．本稿では、非生物（物質）

なりの非線形磁化特性モデルの提唱と準 3 次元非線形磁に於ける相転移（Phase transition、例えば氷が溶けて水へ界解析結果を発表して以来、非線形問題は「力仕事」と変化やガスが燃焼によって酸化化合物へ変化）時に「1/f言う先入観が筆者(齋藤)のトラウマであった．

ゆらぎ」現象が伴う実験的検証例を述べる.

非線形問題のイメージが変わったのは 2007 年の夏休

次に、生物系に於ける「1/f ゆらぎ」現象は精神・心理みであった．某工務店からの依頼で居住空間に於ける直的な側面が投影された行動や行為の結果に表れることを火使用に拠る環境効果の調査をきっかけとして始めた

示す．具体的な例として、精神・心理的変化（Mental「1/f ゆらぎ」抽出に関する研究から「1/f ゆらぎ」生成/Psychological transition：例えば鬱状態から回復へ精神・の研究へ移行するために、チューリング・パターンに関

心理的状態が変化）時の人間が描いた絵画などに「1/fする状態反応拡散系、「箱玉理論」によるソリトン波、ゆらぎ」現象が生ずることを述べる.

「 超 離 散 化 」、 セ ル ラ ー ・ オ ー ト マ ト ン （ Cellular

非生物と生物の中間的な物として自己駆動粒子の「1/fAutomaton）、「渋滞学」を少しかじり始めてから非線形系

ゆらぎ」について述べる．自己駆動粒子とは、例えば高解析のイメージが変わった 1,2)．

速道路を走行中の自動車などが具体的な例であり、自己非線形系の解析は計算機支援の力仕事または超幾何級駆動機能を持つがその行動・動作が制限されるものを言数で代表される複雑な関数を用いた解析解から有限差分う．自己駆動粒子モデルを磁性体の磁区挙動に適用し、を基礎とする「超離散化」などの離散数学と 1984 年に磁性体が磁気飽和へ至る状態遷移（State transition）時に

S.Wolfram によって提唱されたに「セルラー・オートマ「1/f ゆらぎ」現象が生ずることを示す．同様に自己駆動トン」よる複雑系の物理へ繋がったと考える．異方性を粒子モデルの例として、一様な層流からカルマン渦へ変含めた非線形磁化特性問題などは案外「渋滞学」的なセ化する状態遷移時に「1/f ゆらぎ」現象が生ずることを示

ルラー・オートマトンモデルが解決の決め手となると予す．

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006
纏めれば、非生物では相転移時、生物では精神・心理2.2 相転移時における 1/f ゆらぎ
Fig.2 に示すように赤外線カメラを用いて，氷的変化時の行為や行動の結果、そして、非生物と生物の(a)融解中間である自己駆動粒子モデルでは状態遷移時にそれぞの融解状態を撮影した．フレームレートは 1 秒間に 30

れ「1/f ゆらぎ」現象を生ずることを述べる．

フレームである.

最後に、「1/f ゆらぎ」は何故に「癒し（Healing）」効果を与えるのかを考える。癒し効果の実験的検証は脳波のアルファ波が増加することで検証されるが、アルファ波が増加する理由は、「1/f ゆらぎ」を含有する自然現象、
例えば灯された蝋燭などを観察することで人間の「1/fゆらぎ」状態が促され（共振し）、結果として「癒し」へ(a) Experimental equipment

繋がる精神・心理的変化を喚起する点にあると考えられ

Fig.2(b) Infrared Images

Ice Melting by Boiled Water

る．
Fig.2(a)は沸騰した熱湯をビーカーに 150cc 入れ，その２．1/f ゆらぎ解析

熱湯中に氷を入れて，氷の融解状況を 30 秒間撮影した実

2.1

1/f ゆらぎ

検の模式図である．また，赤外線カメラのダイナミック

任意の周期関数 f(t)をフーリエ級数で表現すると、レンジを固定するため，基準温度としてコップに入れた氷をビーカーの横に置いて撮影した.

∞

f (t ) = a0 + ∑ [ai cos(i 2πft ) + bi sin (i 2πft )]i =1⎡⎛ b ⎞⎤= a0 + ∑ ai2 + bi2 cos ⎢2πift − tan −1 ⎜⎜ i ⎟⎟⎥i =1⎝ ai ⎠ ⎦⎣∞∞⎡⎛ b ⎞⎤= a0 + ∑ ci cos ⎢2πift − tan −1 ⎜⎜ i ⎟⎟⎥i =1⎝ ai ⎠ ⎦⎣(1)
White : 1/f Frequency(a)となる．離散値系フーリエ変換では、定数項と余弦波の項が実部、正弦波の項が虚部へそれぞれ対応する。スペクトラムの振幅は実部と虚部のノルムとなる．
（1）式のパワースペクトラムの振幅 ci が周波数 if に(b)
Ice Melting by Boiled Water

Left: 1/f-fluctuation Frequency,
Right: Frequency Fluctuation Distribution

Fig.3

1/f-fluctuation accompanying ice melting

対して反比例する信号を「1/f ゆらぎ」と呼ぶ．フーリ

Fig.2(b)は撮影したフレーム画像例を示す．黒い部分がエ・パワースペクトラム対周波数の関係を両対数グラフ

氷である．Fig.3(a) は Fig. 2(b)に示した氷の融解時の動画に描き，描かれた線の傾きによってゆらぎの種類を大別

像を構成する各フレームの画素値を時系列（フレーム）する．Fig.1 にフーリエ・パワースペクトル対周波数のグ方向に離散値系フーリエ変換し，各周波数に対するフーラフの一例を示す．

リエ・パワースペクトラムから求めた周波数対振幅の傾10

き を 濃 淡 画 素 値 と す る 可 視 化 画 像 で あ る ．Fig.3(b)は0

パワースペクトル

ホワイトノイズ
10

-1

10

-2

Fig.3(a)で求めたフーリエ・パワースペクトラムの傾きが-1.05 から-0.95 の画素値を 1 とし，その他の画素値をゼロとして得られた「１/ｆゆらぎ」周波数の 2 値化画像で1/f ゆらぎある．すなわち，Fig.3(b)は相転移状態における「１/ｆ

単調ゆらぎ」可視化画像であり、融解時の相転移時において、
10

-3

10

-4

10

「１/ｆゆらぎ」成分が存在することを示している 4)．
0

10

1

10

2

10

3

10

4

周波数

Fig.1

Definition of 1/f Fluctuation

Fig.1 で，直線の傾きが 0 の場合は主にホワイトノイズ
である．また、直線の傾きが急なほど単調な信号となる.
そして，ホワイトノイズと単調な信号の中間で直線の傾
きが-1 となる場合を「1/f ゆらぎ」と呼ぶ 5)．

Yellow shows 1/f fluctuation points
Fig.4 1/f fluctuation accompanying gas burning

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006(b)燃焼Fig.4 はガスバーナの燃焼時において空気とガスが化合し酸化化合物へ相転移する黄色部分に「1/f ゆらぎ」が生じていることを示している 4)．
2.3

精神・心理的転移における 1/f ゆらぎ(a)失恋

Fig.5 はゴッホが失恋した時期に描いた絵画の1 例を示す。デジタル画像は光の三原色(赤 R，緑 G，青B)からなっている．1 枚のカラー画像は R,G,B 成分が構成する 3 枚のモノクロ画像からなっている．RGB 各成分
画像はそれぞれ各成分の強度分布が 0 から 255 の 8 ビットのレンジを持つ画素値で描かれている．それぞれのRGB 成分画像の RGB 強度（画素値）のヒストグラムを作成する．RGB 成分画像それぞれの強度分布ヒストグラ
ムを RGB 固有パターンと定義する．この RGB 固有パターンはそれぞれの絵画によって幾何学的形状が変化するにも拘わらず、色彩情報で絵画の特徴変化を吟味可能とFig.6 Blue color component 1/f fluctuation derived fromtime sequentially arranged Blue color histograms of the
Paint Works by Vincent van Gogh
真実の愛」である．このように青には、マイナスのイ
メージや失恋のイメージがある．
ゴッホはこの作品を描いた時期に失恋している．好きな人の両親に自分の気持ちが真剣であることを伝えるた

する．

め，ゴッホは左手を蝋燭の炎にかざすという自傷行為を行っている[5]．このため、この時期に描いたゴッホの作
品に失恋時の心理状態が反映されていると考えられる．
また、弟と画家としての将来について話し合い、自分が画家として生計を立てられるか否かを悩んでいる様子が伺える時期であり，さらに，水彩画の道具一式を尊敬する画家から与えられ，絵画の方向性を変えようとしている時期である．すなわち、画家としての原点に帰り，創造するクリエイティブなエネルギーを生み出そうとしている時期であり、当に精神・心理的な側面が変化状態に
Fig.5 Example of Vincent van Goghʼs paint workwhen loosing sweetheart.
2.3

精神・心理的転移における 1/f ゆらぎ(a)失恋

Fig.5 はゴッホが失恋した時期に描いた絵画の1 例を示す．
ゴッホの個々の作品それぞれに対する RGB 固有パタあると言えよう 5)．
Fig.７はゴッホが精神的に深く病んでいる時期に描いた絵画の 1 例である．この時期、彼は友人のゴーギャン
へコップを投げつけたり、耳を切って売春宿へ送りつけたりなど奇矯な振る舞いをした．Fig.8 は Fig.6 と同様に
RGB 色情報の中で 1/f ゆらぎが観測された B(青)成分を示す．

ーンを時系列に並べ，同一強度（画素値）に対して時系
列方向に離散値系フーリエ変換を適用し，フーリエ・パ
ワースペクトラムを得る．フーリエ・パワースペクトラ
ムの周波数に対する変化からゆらぎ成分を抽出する．その結果 RGB 成分の 0 から 255 までの画素値の時系列の変化率に「1/f ゆらぎ」が存在することを示す．横軸にR,G,B 成分の画素値，縦軸にフーリエ・パワースペクト
ラムの傾きをとり、赤、緑、青成分画に対するゆらぎを描いた．Fig.6 は縦軸の値が最も-1 に近い値を取る B（青）
成分のゆらぎを示す．
青は「寒冷」「冷淡」「陰気」の象徴としてよく使われている。さらに，原点・原初の色であり，創造するクリ
エイティブなエネルギーであり，創造して物を互いにコ
ミュニケートするエネルギーである．また，青の中でも画素値が 50 程度の（淡い）部分に最も「1/f ゆらぎ」に近い揺らぎが見られる。淡い青の花として勿忘草が上げられる．この花の花言葉は「私を忘れないでください・

Fig.7 Example of Vincent van Goghʼs paint workwhen diseasing mental situation.

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006
0

傾き

-0.2
-0.4
-0.6
-0.8

0

50

100

150

200

250

Blue

Fig.8 Blue color component 1/f fluctuation derived fromtime sequentially arranged Blue color histograms of the
Paint Works by Vincent van Gogh
なお、ゴッホの絵画の RGB 色情報を作品が描かれた
時代順に単純に時系列で並べた可視化情報から周波数特性を計算しても 1/f ゆらぎ特性は可視化できなかった．
ゴッホの精神的葛藤が奇矯な振る舞いとして判明している時期の絵画を時系列データとして解析してはじめて

Arrows show the 1/f fluctuation points
Fig.10 RGB 1/f fluctuation parts derived from paintworks of a patient when recovering psychotic depression
disease.
ぎ」が存在し，沈静を求めている様子が判る 7)．
鬱病患者の描いた絵画に於いても単純に描かれた時期順に並べて時系列データとして周波数解析を行っても1/f ゆらぎは可視化できなかった．患者の症状の変化時に描かれた一連の絵画から 1/f ゆらぎが可視化できた．
2.4

状態転移における 1/f ゆらぎ(a)磁性体の磁気飽和ビッター法による磁壁可視化の1/f ゆらぎが可視化できた．

原理について述べる．試料表面の磁壁を観察するために，(b)鬱病

反射型の金属光学顕微鏡を用いる．試料台に薄板状強磁

Fig.9 は鬱病患者が回復期に描いた絵である。

この患者は極めて真面目で仕事熱心な地方公務員から外性体を載せ，強磁性体の表面歪みを電界研磨によって取資系の企業へ転職した方である．しかし、転職後、彼はり除いた研磨面にスポイトで磁性コロイド溶液をたらし，仕事に真面目すぎるために鬱病に侵された．Fig.9 は彼カバーガラスを載せて観察する．強磁性微粒子（γ 酸化が鬱病から回復期において描いた一連の絵を纏めて示し鉄）のコロイド液を強磁性体表面につけるとコロイド液ている．

内に分散していた酸化鉄の微粒子が，磁壁付近の急峻な磁化変化に起因する漏れ磁束の急激な勾配に反応し，表面近くの磁壁の観測ができる．この引き付けられた微粒子を光学顕微鏡で観察する．

Fig.9 Paint works by a psychotic depression patient.
Fig.9 の絵が描かれた順番（時系列）に並べて Figs.6,8と 同 様 に RGB 固 有 パ タ ー ン を 抽 出 し てR(赤),G(緑),B(青)色成分強度（画素値）ヒストグラムから 1/f ゆらぎ求めた結果を Fig.10 に示す．
明らかに精神・心理的変化時（回復初期）に描いた絵
は 1/f ゆらぎを赤、緑、青成分全てに呈していることがFig.11 Experimental device for magnetic wall observation.
Fig.11 は磁壁移動の可視化に使用した実験装置を示す．
励磁コイルを使って周波数 1[Hz]の交流磁界を観測対象に印加する．顕微鏡は KEYENCE 社の HD 中倍率ズーム
レ ン ズ VH-Z75 に 高 精 細 ク イ ッ ク マ イ ク ロ ス コ ー プ

判る．すなわち、回復期に入ると何れの色成分でも「1/f」

VH-5000 を接続したものである．軟鉄の動的磁壁移動をゆらぎが観察される．特に安らぎや癒しを示す緑におい動画像として記録した．Fig.12 は観測したフレーム画像て「1/f ゆらぎ」が見られ、青成分においても「1/f ゆらの例を示す．

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006(a) 0[s](b) 0.14[s](a) Magnetic Field vs. number of 1/f fluctuation parts(e) 0.56[s](f) 0.70[s]Fig.12 Frame images of magnetic wall. Sample: soft iron.
飽和磁束密度は最も高いが鉄損も大きい軟鉄における
1/f ゆらぎ周波数成分（厳密には周波数 f に対して-0.95
乗から-1.05 乗）分布を調べた結果を Fig.13 に示す．
Fig.13 は周波数 f に対して-0.95 乗から-1.05 乗を 1，その他を 0 と二値化して示している．白色の点が 1/f ゆら
ぎ周波数で運動する磁壁部分を示す。尚、1/f ゆらぎ部分の抽出は Figs.3,4 と同様な方法で行った．

White : 1/f frequency
Fig.13 Extracted 1/f frequency fluctuation parts in soft iron.
Fig.13 で、励磁電流の振幅を変更して白色部で示されている 1/f ゆらぎ部分の画素数をカウントした．この結(b) Typical magnetization curve of ferromagnetic
materials.
Fig. 14 Magnetic Field vs. number of 1/f fluctuation
parts and typical magnetization curve of
ferromagnetic materials.

Yellow shows 1/f fluctuation points
Fig.15 1/f fluctuations occur at the boundaries
between layer and eddy flows.
れる 1/f ゆらぎを示す 9)．
レイノルズ数 Re が 0〜約 3000 の状態では、流れは層果を励磁磁界 H 対個数として図 15(a)に示す．図 15(b)流となるが、レイノルズ数が約 3000 を越えると、流れは典型的な磁気飽和曲線である．図 15(a)と(b)を比較すの状態は急激に変化し、大小の渦が入り乱れながら流れることで以下のことが判明した．可逆的磁壁移動範囲でる乱流となる。この層流と乱流の境は明らかに状態遷移

は、1/f ゆらぎ部分は少ないが非可逆的磁壁移動範囲でであり、その結果、層流と乱流の境界で 1/f ゆらぎが観

1/f ゆらぎ部分が劇的に増加する．しかし、可逆的回転磁察される．

化範囲へ入ると 1/f ゆらぎの個数は飽和し、逆に磁界 Hの増加にも拘わらず減り始めることが判る 8)．これは磁３．まとめ性体の磁化過程が磁壁移動から各磁区内の磁化ベクトル

非線形系解析の理論的背景は 21 世紀に入ってかなりの回転へ状態遷移時に「1/f ゆらぎ」が最大となることを体系的に系統化され、その結果、多くの非線形特有な事意味している．この関係は「渋滞学」に於ける自由走行

象が明らかにされた．

対密度の関係に類似している 3)．
(b)流体Fig.15 は層流とカルマン渦間の境界に生成さ本稿では、物質の相転移（Phase transition）時に「1/fゆらぎ」現象が伴う実験例としてガスが酸素と化合して

可視化情報 Vol. 28 Suppl. No.1(2008年7月）
P01-006
酸化物へ変わる相転移時に生成される「1/f ゆらぎ」と氷

3) 西成活裕、渋滞学の新しい世界、東京大学ホームペが融解する相転移時に生ずる「1/f ゆらぎ」を示した．

−ジ、http://soliton.t.u-tokyo.ac.jp/nishilab/http://soli
ton.t.u-tokyo.ac.jp/.

次に、失恋時のゴッホの作品や鬱病患者が回復時に描いた絵の色情報から「1/f ゆらぎ」を抽出し、生物（人間）4) 寺西正晃，丸山和夫，早野誠治，齋藤兆古，堀井清之，自然界の画像が持つ 1/f 周波数成分の可視化，可視化

における「1/f ゆらぎ」現象は精神・心理的変化（Mental

情報シンポジウム,B108,2005.

/Psychological transition）時に人間の行動や行為の結果に表れることを示した．

5) 丸山和夫，早野誠治，斎藤兆古，堀井清之：色情報を利用した知的動画像認識，可視化情報学会誌，vol.23，磁性体の磁気飽和特性へ至る状態遷移時と流れが一様

No.1，pp.95-98，2003

な層流と渦が生ずる状態遷移境界に「1/f ゆらぎ」が生成された具体例を用いて、非生物と生物の中間的な物を意

6)

宮坂総，齋藤兆古，加藤千恵子,” 動画像の色彩情報

味する自己駆動粒子モデルにおいては状態遷移（State

可視化とその応用,” 第35回可視化情報シンポジウム 工学

transition）時に「1/f ゆらぎ」が生ずることを示した．

院大学2007年７月25日, C213, Vol.27, Suppl. No.1 (2007年7月)pp. 227-228

参 考 文 献7)加藤千恵子、法政大学学位論文、2007年3月

1) 東 京 大 学 ホ ー ム ペ ー ジ 、 http://soliton.t.u-tokyo.ac.jp

8)

須永高志，寺西正晃，齋藤兆古, “ビッター法による可視化
磁区画像から周波数特性の抽出,” 日本AEM学会誌/nishilab/mypapers/rikouJ.pdf.
2) 高橋大輔、超離散化された波、特集/<波の世界>の
魅 力 と 数 理 , 早 稲 田 大 学 ホ ー ム ペ ー ジ 、 http://takahashi.
math.sci.waseda.ac.jp/works/public/03-Surikagaku-478-35-prepri
nt.pdf.

Vol.15,No.2(2007) pp.195-200
9)古川裕之、カルマン渦の実験画像、名城大学ホーム
ページ、http://mech.meijo-u.ac.jp/prof/furukawa//mj/gazoupage1.htm

zemi

213

第 51 回自動制御連合講演会
2008 年 11 月 22 日，23 日 山形大学工学部

1/f ゆらぎと複雑系に関する考察
○齊籐兆古

宮坂総

菅井桂子（法政大学）

Study of 1/f Fluctuation Frequency and Complex System* Y.Saito, S.Miyasaka and K.Sugai (Hosei University)
Abstract－－ Based the Newton mechanics and continuous mathematics, recent engineering sciences havedeveloped various tools for modern human society such as high speed traffic/transport, high quality consumer
electronics, intelligent housing and enormous information network systems. One of the most effective methodologies for next generation of the engineering sciences may be a nonlinear complex system methodology, because the human is one of the complex systems. This paper describes an extracting method of the 1/f fluctuation
frequencies intrinsically accompanying non-linear physical system operation such as non-life being, life-being
and self-driven particle systems.
Key Words: 1/f fluctuation frequency, Complex system

1

まえがき

る．ここで、 i = 1, 2,.., ∞, ω = 2π / T = 2π f (角周波数)とする．
すなわち、21世紀の今日まで、人類の科学技術は、主としてニ
ュートン力学と古典解析学を武器として、極めて多岐∞に渡る多くの文明の利器を提供した．特に、20世紀末
f ( t ) = a0 + ∑  ai cos ( iωt ) + bi sin ( iωt ) に開発され爆発的な普及を遂げたコンピュータは、IT
i =1
関連産業を喚起し、その結果、あらゆる業種の在り方∞ b 
へ広汎な影響を与えた．コンピュータは従来型プロダ= a0 + ∑ ai2 + bi2 cos iωt − tan −1  i  
i =1
クトのインテリジェント化や多機能化を可能としただ ai  
けで無く、新しい科学技術の方法論を与えんとしてい∞ b 
る．すなわち、人間の物理的機能の強化のみならず脳= a0 + ∑ ci cos iωt − tan −1  i  
i =1
 ai  
機能の補完や情緒・精神面を勘案したプロダクト開発
のキーとなる非線形な複雑系の解析を可能とする．
（１）
筆者等は、家庭電化に伴う直火を使わない生活空間
が成り立つ．
が人間に与える影響を調べるため、燃焼現象に伴う炎(1) 式の高調波次数 i と高調波の振幅 ci の関係を両の1/f揺らぎ周波数解析を行い、これを基点として、化対数でFig．1に示すように描く．
学反応などの非生物系における1/fゆらぎ1)、人間の情
10
緒・精神活動に伴う1/fゆらぎ2)、さらに磁性体の磁区3)ホワ イトノイズ
挙動の1/fゆらぎ などを解明した．1/fゆらぎは非線形10な複雑系で観察される周波数特性であり、非線形現象の代表的な特徴である．非線形系は従来の線形系と比1 /f ゆ ら ぎ
べて現実の物理系を忠実に表現可能とするのみならず、
10
単調
人間の感性、例えば癒し(healing)効果などがあるとされている．すなわち、非線形な複雑系は、従来の単純10なヒューマン・インターフェイスと一線を画する人間の感性を前提とした機械と人間のインターフェイス構10
築の基幹となる一方法を提供する．
10
10
10
10
10
周波数
複雑系をシミュレーションするには系のマクロ的振る舞いを現す支配方程式を解くことによるトップダウ
Fig.1 Definition of 1/f fluctuation frequency
ン型と系を構成する局所部分の振る舞いを積み上げて
Fig．1で、周波数に無関係に振幅（パワースペクト
系全体の振る舞いを表現するボトムアップ型がある．
ル）が一定値をとる周波数特性はパルスやホワイトノ
前者は古典的な解析手法の拡張であり、後者はコンピイズに見られ、人間の感性ではランダム性が大きく不ュータを前提とする離散値系の手法であり、その代表快感を与える．他方、高周波数になると急激に振幅がとしてセルラー・オートマトン(cellar automaton)がある
減衰する波形は単調な信号、例えば単純な正弦波など
4)．多くの非線形な複雑系はこのセルラー・オートマ
であり、人間の感性に対して単調すぎて飽きられる感
トンで記述される5)．
覚を与える．
周波数に反比例して振幅が減衰する特性、以上のことを鑑み、本稿では筆者等が提案する1/fゆ
すなわち、振幅が1/fに比例して減衰する周波数特性をらぎ解析法と複雑系に関するいくつかの考察を行う．
1/fゆらぎ周波数特性と言い、人間の感性に対して心地よい感覚を与える．心地よさの測定は脳波でα波の発
2 1/fゆらぎ周波数の抽出
生度合いを測定することで検証される．
2.1 1/f ゆらぎ周波数とは
2.2 1/f ゆらぎ周波数を持つ波形の性質
0

パワースペクトル

-1

-2

-3

-4

0

一般に任意の周期波形 f ( t ) は平均値 a0 、余弦波
ai cos ( iωt ) および正弦波 bi sin ( iωt ) の和で表され

1

2

3

4

乱数を用いて生成した 1/f ゆらぎ周波数の波形を Fig．

08PR0003/08/0000-0064 ¥400 ©2008 SICE
64

2 に示す．Fig．2 のフーリエスペクトラム対高調波次数の関係をそれぞれの対数で描くと Fig．
3 となる．Fig．
3 で、直線は最小自乗法で得られた直線近似であり、その勾配は-1.002 であり、1/f ゆらぎ周波数特性が抽出された．

カラー動画像の任意の1フレームは赤、緑および青成分からなる．赤、緑、青画像はそれぞれ独立なモノク
ロ画像として表現できる．
Fig．5は時間 t1 , t2 , t3 におけるそれぞれのモノクロフ
レーム画像を示す．Fig．5で、全体のフレーム画像に共通なx,yスクリーン座標上の任意の画素位置における画素値を式(1)の周期関数 f ( t ) に対応させて時間(フレーム)軸方向の1/fゆらぎ周波数特性を求める．すなわち、画素(pixel)毎にフレーム方向に生ずる1/fゆらぎ周波数特性を求める．このため、動画像の1/fゆらぎ
周波数特性は1枚の静止画として抽出される．

Fig. 2 An example of waveform having
1/f frequency characteristic.
Fig. 4はFig. 2に示す波形のサンプリング個数を半分、すなわち、全データが65536点からなる波形の前半部分
32768個からなる波形のフーリエ・パワースペクトラム
対高調波次数の関係である．Fig. 4で、直線は最小自乗法で得られた近似であり、その勾配は-0.978であり、
1/fゆらぎ周波数特性は、やや大きな誤差であるが、近似的に抽出された．同様にFig．2の波形で、第16384
から第57344点まで40960個の波形から得られた周波数
Fig.5 Sequential 3 frame monochrome images
特性を直線近似して得られる勾配は-0.998となった．
拠って、1/fゆらぎ周波数特性を呈する波形はサンプ
2.4 フィルター
リング個数に拠らず、1/fゆらぎ周波数特性を近似的にいま、横軸の画素位置 j、縦軸の画素位置 k におけ与える．その近似精度はサンプリング個数に比例する
る第i 次高調波のフーリエ・パワースペクトラムをことが判明した．

Fj ,k ( iω ) = ci j ,k , i = 0,1,.., l − 1

(2)

とすれば、
y j ,k =  log F j ,k (ω )log Fj ,k ( 2ω )T

⋅ log Fj ,k ( lω )  ,x j ,k = [log ω log 2ω ⋅ log lω ] ,i = 1, 2,.., m, j = 1, 2,.., nT(3)に対して、直線近似
Fig.3

y = a + bx

log ci versus log i characteristic

Straight line denotes the least squares approximation
having -1.002 gradients.

(4)を適用するとすれば、以下の線形システム方程式が成り立つ．
Y j ,k = C X j ,k ,Y j , k =  log F j , k (ω )X j ,k = [ a

Fig.4

1
1
C =
⋅

1

log ci versus log i characteristic

Straight line denotes the least squares approximation
having -0.978 gradientslog F j , k ( 2ω )

⋅

b] ,
T

log ω log 2ω ⋅ log lω , j = 1, .., m ,k = 1, .., n(5)(5)式の線形システム方程式の係数行列 C は l（フレーム数）行 2 列であるから不適切(ill posed)である．こ2.3 動画像から 1/f ゆらぎ周波数特性の抽出

65

Tlog F j , k ( lω )  ,のため、(5)式の近似解ベクトル Xij,k*は誤差ノルム

r j , k = Y j , k − CX j , k *(6)を最小にする最小自乗法による解(Xk , j* = C T C

) CY
−1

T

j ,k

(7)

を採用する．
尚、実際の計算ではフーリエ係数の精度を勘案して全フレーム数の 1/4 項までに対応するフーリエ・パワースペクトラムから(4)式の直線近似の係数 a, b を決定した．Figs．3,4 が得られた直線近似の例である．Fig．
6 の直線近似は比較的よく成り立っている．実際は、近似がどの程度正しいかは(6)式の残差ノルム|rj,k|から評価できる．換言すれば、直線近似の精度を勘案して有意な 1/f ゆらぎ周波数部分のみを抽出することが可能である．例えば、(6)式から全体の平均誤差を

Fig．7 An Example of Music Frequency Characteristicが交互に配置されていることが判る．
3.2 動画像で計算し、Fig．6 に示すようにこの平均誤差を閾値とするフィルターを作成し最小自乗法が一定の精度以上で成り立つ部分のみを抽出することが可能である．

大部分の自然界の 1/f ゆらぎ周波数特性は動画像から抽出される．化学反応などの非生物系における 1/f
ゆらぎ 1)、人間の情緒活動に伴う 1/f ゆらぎ 2)、さらに磁性体の磁区挙動の 1/f ゆらぎ 3)などが報告されている．
Fig．8 は楽器（三味線）演奏に伴う 1/f ゆらぎ周波数の抽出例である．白のドットが 1/f ゆらぎを呈する
画素を示す．尚、この結果は最小自乗法の誤差ノルムを使った空間フィルター(2.4 節を参照)を使っている．

Fig．6 An example of space filter

Fig．8 1/f Fluctuation Frequency Extraction from a Music
Instrument (Shyamisen) playing
White dots denote the 1/f fluctuating pixel position．

mε mean =

3

n

∑∑ r

j ,k

j =1 k =1

m× n

(8)

4

例題

セルラー・オートマトン

3.1 音楽4,1 歴史的背景

音楽データはいわゆる音であり一次元データである
ため、比較的に容易に解析可能である．音楽データの周波数特性は、1)音楽全体として 1/f ゆらぎを呈するもの、2)一定の周波数帯域で 1/f ゆらぎを呈するもの、さらに 3)全体としても周波数帯域別にも 1/f ゆらぎを呈さないものに大別できる．
Fig．7 はいわゆる癒し効果を与える音楽として市販されている曲のゆらぎ解析結果である．Fig．7 で、直線は最小自乗近似に拠る周波数特性の傾きを表し、緩やかな傾きはほぼ-1 であり、急峻な傾きは-1 よりも小さい値である．直線の傾きを観察すれば、この音楽は周波数帯域毎に単調な旋律と 1/f ゆらぎを呈する旋律

セルラー・オートマトンは 1950 年代にノイマンとウラムが考え出したものである．1970 年代から生物の形
体生成のモデルに利用され、1986 年には格子ガスオー
トマトン法へと進化した．その後、格子ボルツマン法に発展し、複雑な流れの解析などに適用されてきた．
今日では物理･科学現象、材料･交通･電子回路特性、社会･経済現象などに幅広く適用されている 4)．
他方、スティーブン・ウルフラム（Stephen Wolfram、
1959 年 8 月 29 日 - ）はアメリカの Wolfram Research社の創業者で最高経営責任者であり、また、理論物理学者でもある．彼は 15 歳にして素粒子論の学術論文を執筆し、オックスフォード大学を 17 歳で卒業．その後

66

カリフォルニア工科大学に進み、高エネルギー物理学、場の理論、宇宙論の研究を行った．20 歳で理論物理学の研究によりカリフォルニア工科大学において Ph．D．
の学位を取得した．1982 年より現在では『複雑系』に分類される自然界の複雑さについて研究し、セル・オ
ートマトンに関する革新的研究を行った．さらに彼はそれらの成果を Review of Modern Physics 誌に掲載した 5)．
4.2 磁性体における複雑系a)磁区 磁性体は一般に原子レベルの大きさで電子スピンに起因する磁石を持っている．これらの微少磁石
は外部磁界に応じて系に蓄えられる磁気エネルギーを最小にするように運動する．原子レベルの磁石を個々に観察するのは不可能であるが、通常、磁性体中では微小磁石が複数個の凝集した塊の集合を形成する．これらの塊の集合（これを磁区と言う）、すなわち、磁区は電子顕微鏡や金属顕微鏡で観察される．

Fig．11 Preisach Modelー・オートマトンモデルである．これは、コンピュー
タが出現する以前の 1930 年代に於いても複雑な物理現象を表現する方法はセルラー・オートマトン型にならざるを得ないことを意味し、興味深い．
自己駆動粒子モデルで表現される道路の渋滞推移は
1/f ゆらぎ周波数特性を呈する．そしてそれらのシミュ
レーションモデルはセルラー・オートマトンで構築される 7)．さらに非線形波動の代表であるソリトン波も箱玉モデルと呼ばれるセルラー・オートマトンで表される 8)．
拠って、セルラー・オートマトンは複雑系解析の有力な一方法である.

Fig．9 A Magnetic Domain Image by Bitter Method
Sample: soft iron5. まとめ本稿では信号から1/fゆらぎ周波数特性を抽出する
方法とそれらの具体的例、および、セルラー・オート
マトンモデルに関して概述した．

参考文献
1) 寺西正晃，丸山和夫，早野誠治，齋藤兆古，堀井清之，自然界の画像が持つ 1/f 周波数成分の可視化，可視化情報シンポジウ
ム,B108,2005．

Fig．10 1/f Fluctuating parts accompanying with magnetic
domain movements．
White dots denote the 1/f fluctuating pixel position

2) 宮坂総，齋藤兆古，加藤千恵子,” 動画像の色彩情報可視化とそ

Fig．9 は表面磁区のみであるがビッター法で得られた磁区画像の一例である 3)．これらの磁区状態は外部磁界の強度と極性に対して運動し、Fig．10 に示すように 1/f ゆらぎ周波数特性を呈する．
b)Preisach Model F.Preisach は磁性体の磁化特性を表現するため磁化方向を軸とする 2 次元平面上で、外部磁界強度に応じて単位磁化の極性が反転するモデル、
すなわち、Fig．11 に示す Preisach モデルを考えた 6)．
Fig．11 で、∓ がゼロ、+1 が正方向の単位磁化、-1 が負方向の単位磁化であるから、明らかに Preisach の磁化モデルは外部磁界に応じて 3 値を取る一種のセルラ

3) 須永高志，寺西正晃，齋藤兆古, “ビッター法による可視化磁区の応用,” 第 35 回可視化情報シンポジウム 工学院大学 2007 年
７月 25 日, C213, Vol．
27, Suppl．No．
1 (2007 年 7 月)pp．227-22867

画像から周波数特性の抽出,” 日本 AEM 学会誌 Vol．15,No．
2(2007) pp．195-200
4) フリー百科事典『ウィキペディア（Wikipedia）』セル・オート
マトン、ja.wikipedia.org/wiki/セル・オートマトン
5) フリー百科事典『ウィキペディア（Wikipedia）』スティーブン・
ウルフラム、ja.wikipedia.org/wiki/スティーブン・ウルフラム
6) F.Preisach, Zeitscrift fur physic, 94, No.5 (1935)7) http://soliton.t.u-tokyo.ac.jp/nishilab/http://soliton.t.u-tokyo.ac.jp/
8) http://takahashi.math.sci.waseda.ac.jp/works/public/03-Surikagaku-478 -35-preprint.pdf

第 51 回自動制御連合講演会
2008 年 11 月 22 日，23 日 山形大学工学部

214

セルラー・オートマトンによる磁化特性の表現○宮坂総齋藤兆古（法政大学）

Representation of Magnetization Characteristics by Cellar Automaton* S.Miyasaka and Y.Saito (Hosei University)
Abstract－ After we classified the magnetization curve of a ferromagnetic material into three magnetization
regions, we applied the Preisach type cellar automaton model to each of three magnetizing regions. As a result,
it is shown that the magnetization curves can be represented by the cellar automaton model.
Keywords：1/f Fluctuation, Preisach Model, Cellar Automaton

1

まえがき

0,(1 0 0) →1,(0 1 1) →1,(0 1 0) →1,(0 0 1) →1,(0 0 0) →0)を適用した場合をFig.1に示す．
t=0
t=1
t=2
t=3
t=4

セルラー・オートマトンは、1940年代にノイマンによって提唱された．その後、1984年にウルフラムによって複雑系のシミュレーション行う決定的な方法として紹介され、近年、多くの分野で新世代シミュレーシ
ョン技術として研究・開発が行われている．特に、従来の方法ではモデル化が困難であった交通渋滞、
雪崩、
さらに火事などの自然災害のシミュレーションにも用いられており、今後さらなる研究開発と実用化が期待されている1-3)．
一方、磁性体は極めて複雑な磁化特性を呈するが、1930年代にこの磁化特性を表現するモデルがプライザ
ッハによって提唱された．このプライザッハモデルは極めて具体的な磁区理論に基づくため、他に代替えとなるモデルが無く、現代でも磁気記録理論や磁性体を含む磁界計算に広範に使われている4)．
我々は既に磁化特性を表現するプライザッハモデル
が、未だセルラー・オートマトンの概念が提唱されて無いにも拘わらず、一種のセルラー・オートマトンモ
デルであることを示した．
本論文では、まず古典磁区理論に基づき磁性体の磁化過程を3領域、すなわち、可逆的磁壁移動、非可逆的
磁壁移動、さらに磁化ベクトルの回転磁化領域へ分割した後、これらの3領域へプライザッハ型セルラー・オ
ートマトンモデルを適用した．その結果、磁性体の磁化特性がセルラー・オートマトン型モデルで表現可能
であることを報告する．

2

Fig.1 Example of Cellar Automaton based on Rule 30.

3磁化特性のセルラー・オートマトンモデル

magnetic flux density

3.1 プライザッハモデル
セルラー・オートマトンを用いたプライザッハモデ
ルのシミュレーションを行った結果をFig.2に示す．

magnetic

field

Fig.2 Magnetization Curve computed by Cellar Automaton,
Corresponding to the Preisach Model.

Fig.2 のシミュレーションに使ったセルラー・オート
マトンを Fig.3 に示す．Fig.3 で、黒色部分は＋１の値、灰色部分は－１の値を持つ．Fig.3 で､黒に転移する場
合はプラスに磁化,灰色に転移する場合はマイナスに磁化されると考える．Fig.2 は､Fig.3 の個々の画像の総和を縦軸の値、横軸は磁界とする．Fig.3 はプライザッ
ハが提案したプライザッハ線図そのものに他ならない．

セルラー・オートマトン

2.1 セルラー・オートマトンとはセルラー・オートマトンとは、格子サイトがいろい
ろな初期値を持つ不連続な系からなり、これらのサイ
トはそれぞれのサイトがいくつかの局所的な隣接サイトの値に基づいた新しい値と有限な数の過去の時間ス
テップを仮定すると、離散的な時間ステップで状態変化する系を与える1-3)．また、セルラー・オートマトンには、空間の広がりに応じて 1 次元、2 次元、さらに 3次元のセルラー・オートマトンがある．
2.2 セルラー・オートマトンの例
Fig.1は初期値(0 1 0 0 1 1 0 1 0 1 1 0 0 0 1 0)
でセルラー・オートマトンのルール30((左近傍 tの値
右近傍)→t+1の値：(1 1 1)→0,(1 1 0) →0,(1 0 1) →Frame1

Frame10

Frame21

Frame30

Frame41

Frame50

Fig.3 Cellar Automaton Images representing the Magnetization
Processes

08PR0003/08/0000-0068 ¥400 ©2008 SICE
68

区挙動画像の画素値からプライザッハ型セルラー・オ
ートマトンの磁化値を決定する．その結果得られたプライザッハ(型セルラー・オートマトン)線図を Fig.6 に示す．Fig.7 は Fig.6 から再現された磁化曲線である．

3.2 磁区挙動に基づく解析
磁性体表面磁区挙動をビッター法によって可視化した．得られた磁区挙動画像を、磁化過程によって 3 過程に分類した．磁性体の磁化過程は、可逆的磁壁移動
範囲、非可逆的磁壁移動範囲さらに回転磁化領域分けられる．本論文では、非可逆的磁壁移動範囲について考える．非可逆的磁壁移動範囲で、正に単位磁化された部分を白色（数値１）、磁化されてない部分を黒色（数値ゼロ）として 2 値化した．その結果、得られた磁区画像 Fig.4 は Fig.3 の白色部が空間的にランダムに位置するセルラー・オートマトン画像に対応する．
Fig.4 から Fig.3 と同様にして Fig.5 の磁化特性が得られる．

Frame1

Frame5

Frame10

Frame15

Frame20

Frame25

Fig.6 Cellar Automaton Images representing the Magnetization
Processes taking into account the practical magnetization
B@TD
1

0.00[s]0.167[s]0.5-10000

-5000

5000

10000

H@AêmD-0.5-10.333[s]0.50[s]Fig.7 Magnetization Curve by Cellar Automaton Images obtained
from the Preisach type cellar automaton modelFig.4 Binary Images of Magnetic Domains.

4 まとめ
B@TD
1

本論文では、最初に、2 次元セルラー・オートマト
ンから古典プライザッハ型磁化特性モデルを導き、レイリーループを描いた．次に、ビッター法によって可視化された磁区画像を、磁化過程の相違で 3 領域へ分類し、分類された領域の中で可逆的磁壁移動範囲に対する磁区画像を画素値に応じて１かゼロへ 2 値化してセルオートマトンモデルを作成した．最後に、実際の磁区挙動画像から得られたセルラー・オートマトンモ
デルの状態遷移ルール、および画素値からプライザッ
ハ型モデルのルールと磁化値を採用し、実際に観測した磁化特性に対応する磁化特性の再現に成功した。

0.5-10000

-5000

5000

10000

H@AêmD-0.5-1

Fig.5 Magnetization Curve by obtained Binary Images ofMagnetic Domains

しかし、Fig.4 の磁区挙動画像から得られた Fig.5 の磁化特性と 3.1 節で述べたシミュレーション結果には多少の相違がある．これは、Fig.4 に示す磁区画像は本来、中間的なグレイレベルを持つにも拘わらず、単純な 2 値化を行ったことに起因する誤差であり、磁化曲線の反転部で顕著に誤差が反映することに拠る．Fig.4の 2 値化した磁区画像から、Fig.3 のプライザッハ型セルラー・オートマトンモデルの状態遷移ルールを抽出する． Fig.4 から以下の状態遷移ルールが導かれる．
徐々に磁界を増加すると、しばらくはゆっくりと磁化され、その後、磁界が閾値を越えると、急激に磁化される．さらに、磁界の方向が反転する直前では、磁界が弱くなるためゆっくりと磁化される．これらのセル
ラー・オートマトンの状態遷移ルールをプライザッハ
型セルラー・オートマトンモデルへ適用し、Fig.4 の磁

参考文献
1)
2)

3)

4)

69

R.J.ゲイロード/P.R.ウェリン、荒井隆(訳)：Mathematica
複雑系のシミュレーション物理学と生物学の探求
R.J.ゲイロード/西舘数芽(著)、西舘数芽/西川清(訳)：Mathematica 自然現象の計算モデル化セルラーオートマ
タ・シミュレーション
西成活裕：セルラー・オートマトンによる複雑現象のモ
デル化、東京大学ホームページ
http://soliton.t.u-tokyo.ac.jp/nishilab/mypapers/rikouJ.pdf
須永高志、齋藤兆古、堀井清之：ビッター法による磁性材料の鉄損分布の可視化、可視化情報学会シンポジウム(2007) C207

第 51 回自動制御連合講演会
2008 年 11 月 22 日，23 日 山形大学工学部

215

サウンドデータにおける 1/ｆゆらぎ周波数の抽出○菅井桂子

齊籐兆古（法政大学大学院）

1/f fluctuation frequency extraction from sound data
* K. Sugai and Y. Saito (Graduate School of Hosei University)
Abstract－Most of the modern conveniences are designed to operate in a highly efficient manner with sophisticated human interface. According to the spreading use of micro-computers, modern conveniences are able toequip the human intelligent functions. Further, modern ergonomics have to take the mental environmental impression into account. Sound is one of the most effective signals appealing to human sensitivity. Therefore, it ispossible to say that sound is an extremely important factor when we design the intelligent human interface ofmodern conveniences as well as ergonomics considering the human mental impressive effects into account.
Since 1/f fluctuation frequency characteristic may give a healing effect to the human, this article makes a studyof the 1/f fluctuation frequency characteristics extracted from various types of music.
Key Words: 1/f Fluctuation, human interface, sound1.

2． 1/ｆゆらぎとは

緒論

現代の多くの機械は、単に高性能が要求されるだけ「１/ｆゆらぎ」は自然界に多く存在し 2)，例えば小で無く、より洗練されたヒューマンインターフェイス
川のせせらぎ，小鳥の囀り，爽やかなそよ風などの心を備えなければならない．この意味で、現代の人間工安らぐリズムが相当する．同様に，心地良い音楽を聴学は取り扱い易さに加えて人間へ与える心理的・精神
いたり，快い感じを抱いたり，安静にしているときの的な影響も考慮しなければならない．
脳波にも「１/ｆゆらぎ」が存在する．
音響は人間の感性へ訴える最も効果的な信号である．
「１/ｆゆらぎ」解析法として，信号へ離散フーリエこのことは近年、音楽療法が多くの医療施設で取り入変換を適用し，各周波数に対するパワースペクトラム
れられるなど、音楽の効果的な作用が医学的に認めらを計算する．周波数の低下とともにパワースペクトラ
れていることからもわかる．音楽には様々な病気の改ムが増加するような信号の中で，パワースペクトラム
善効果や癒し効果があるとされている．
の振幅が周波数に対して反比例する信号が「１/ｆゆら
音楽療法とは病気の改善や治療を目的に実施されるぎ」である．
テラピーである．慢性分裂病や自閉症などの改善、高視覚的に判りやすくするために，よく行われる方法
血圧や脳血管障害の克服、
老人性痴呆症の予防や治療、 は，フーリエ・パワースペクトラム対周波数の両対数
がん患者の延命効果、心身の緊張や痛みの緩和効果なグラフを描き，描かれる線図の傾きによってゆらぎのど、種々の医学的作用を持つことが証明され、すでに種類を大別する方法である．Fig.1 にフーリエ・パワー
多くの医療分野で導入されている．
この音楽療法には、 スペクトラムの例を示す．
効果的な音楽を聴覚情報として耳から聞き入れる受動
Fig.1 において，直線の傾きがゼロの場合は主にホワ
的音楽療法と、楽器演奏や歌を唄うといった能動的音
イトノイズである．また，直線の傾きが急になる程単楽療法がある．
調な信号である．そしてホワイトノイズと単調な信号
人間が耳で感じ取れる周波数は、医学的にみて15ヘの中間的な信号で傾きが約-1 の場合を「１/ｆゆらぎ」
ルツから20,000ヘルツといわれている．耳から入力さと呼び，人間が心地よいと感じる信号であるとされてれる音の周波数は人間の脳から脊椎にある各骨格部位
いる．
と対応しており、例えば延髄より上の脳神経系は4,000[Power]10ヘルツ以上の高周波音に対応し、頸椎は2,000から3,000
ホワイトノイズ
ヘルツ、また胸椎は800から2,000ヘルツ、腰椎から仙10
椎は125から800ヘルツという周波数に呼応している．
そのため音楽に含まれる高周波音は、副交感神経の分
1/f ゆらぎ
10布する延髄から大脳にかけての神経系を刺激し、その単調結果、脳神経系、ホルモン系、循環系、免疫系といっ
10た人間の健康を支えている生体機能に効果的に働く1)．
拠って、音はものづくりにおけるヒューマンインタ
10
ーフェイスを考える上で極めて重要な因子であるとい
10
10
10
10
10 [Frequency]
周波数
える．本論文では音楽中に含まれる1/f ゆらぎ周波数
Fig.1 1/f fluctuation frequencyに関して幾つかの考察を述べている．
0

パワースペクトル

-1

-2

-3

-4

0

1

2

3

4

08PR0003/08/0000-0070 ¥400 ©2008 SICE70

3. 音楽の 1/ｆゆらぎ周波数特性
3.1 実験方法
まず、音楽データの信号を時間領域にフーリエ変換
し、
フーリエ係数の絶対値を計算してフーリエ･パワー
スペクトラムを得る．得られたフーリエ・パワースペ
クトラムから、フーリエ・パワースペクトラム対周波数の両対数グラフを描く．
フーリエ変換は全サンプル数に等しい実部と虚部を与えるので、独立なパワースペクトラムは全サンプ
ル数の半分までであり、さらにフーリエ変換の精度を勘案し全サンプル数の 1/4 項までで高調波解析を行う．
フーリエ･パワースペクトラム対周波数の両対数図か
ら周波数特性を直線近似する一次関数を最小自乗法によって求め、その比例項（傾き）から揺らぎの種類を大別し特徴を吟味する．本稿では全体の傾きに加え、部分的な周波数帯域における傾きも求めた．

Fig.2 Japanese music(Syamisen)

3.2 実験結果
和楽やクラシック音楽、民族音楽、現代音楽など、幅広いジャンルで実験を行った。計算機のメモリー制
約から、曲中で主題の展開部にあたる部分の約1分から1分30秒を抜粋して計算した．
その結果、全体を通してのフーリエ･パワースペクト
ラムの傾きが-1で1/ｆゆらぎを呈する音楽、全体で見ると傾きが急で一見単調だが、部分的に見ると1/ｆゆらぎが存在する音楽、逆に全体で見ると傾きはゼロであるが部分的に1/ｆゆらぎが現れている音楽など、ほとんどの曲で1/ｆゆらぎ周波数特性の存在が確認できた．しかし、一方で、1/ｆゆらぎが現れても微細であり、ほとんど確認できない曲も存在した．
1/ｆゆらぎが確認できた例をFig.2とFig.3に示す．
Fig.2は三味線であり、全体で見ても部分で見ても1/ｆ
ゆらぎが現れていた．Fig.3はクラシックの管弦楽で、全体の傾きは-1.4と多少急であったが、部分的に見ると2箇所で1/ｆゆらぎが現れていた．
次に、ほとんど1/fゆらぎが存在しなかった例をFig.4とFig.5に示す．Fig.4は南米の民族音楽で、今回抜粋した部分ではほとんど傾きが無かった．Fig.5はピアノ協奏曲の一部だが、傾きが-2の部分の他にわずかな範囲
で1/ｆゆらぎが確認できた．

Fig.3 Classic musicFig.4 Folk music of South America

4. まとめフーリエ変換と最小自乗法を用いて音楽データから
1/f ゆらぎ周波数成分を抽出した．その結果、様々なジャンルの音楽で1/f ゆらぎ周波数成分の存在が確認され、さらに、それらは個々の音楽で固有の様式で周波数帯域に組み込まれていることが判明した．

Fig.5 Piano Concerto
2)

5. 参考文献1)

和合晴久,新･健康モーツァルト音楽療法 PART2:血液循
環系疾患の予防,P1/P13 (2004 年)

71

寺西正晃，丸山和夫，早野誠治，齋藤兆古，堀井清之，
自然界の画像が持つ 1/f 周波数成分の可視化，可視化情報シンポジウム,B108,2005.