超音波技術（R言語）
超音波伝搬状態の測定解析

超音波システム研究所

１． 準備
解析用データの確認測定データ（PSDATA ファイル）から解析用のMicrosoft Office Excel CSV ファイルがホルダーにあることを確認してください
ファイルがない場合は、CSV ファイルに保存してください

例D:¥20120321S¥ 20120321-0101.csv

２． 解析ソフトの立ち上げ

ダブルクリックして立ち上げる＞｜

左記のようなプロンプト表示が行われます

エラー表示が行われた場合には、その他の作業ファイル・・・を終了してからもう一度立ち上げてください３． 解析ソフトの読み込み
３－１：パッケージ－＞
ローカルにあるｚｉｐファイルからの・・・・

３－２：デスクトップのｚｉｐファイル
ＴＩＭＳＡＣ1.2.1 を選択する

３－３：読み込み画面の確認

３－４：パッケージの読み込みパッケージ－＞パッケージ読み込み・・・

選択する３－５：ｔｉｍｓａｃの選択

３－６：「ＯＫ」選択により読み込み以上で解析準備完了です参考

バイスペクトルバイスペクトルは以下のように
周波数 f1、f 2、f1 + f 2 のスペクトルの積で表すことができる。
B( f1 , f 2 ) = X ( f1 )Y( f 2 )Z ( f1 + f 2 )主要周波数が f1 であるとき、f1 + f1 = f 2、f1 + f 2 = f3 で表される f 2、f3 という周波数成分が存在すればバイスペクトルは値をもつ。
これは主要周波数 f1 の整数倍の周波数成分を持つことと同等であるので、バイスペクトルを評価することにより、高調波の存在を評価できる。
詳しい説明は専門書・・・を読んで確認してください

エクセルファイルのデータ列

解析コマンド
dev.off()par(mfrow=c(2,2))：

2 行 2 列のグラフ表示設定

data11 <- read.table("C:/20191220/20191220-0030/20191220-0030_12.csv",
skip=6, sep=",", nrows=6000)

plot(data11$V2)
data11 <- read.table("C:/20191220/20191220-0030/20191220-0030_12.csv",
skip=6, sep=",", nrows=6000)

spectrum(data11$V2,method="ar")data11 <- read.table("C:/20191220/20191220-0030/20191220-0030_12.csv",
skip=6, sep=",", nrows=6000)

bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0030/20191220-0030_12.csv",
skip=6, sep=",", nrows=6000)

autcor(data11$V2)

dev.off()
解説 終了
plot(data11$V2)
解説 data11 の 2 番目のデータ列（１ｃｈの測定データ）に対してプロット（音圧測定データのグラフ作成）を行う
spectrum(data11$V2,method="ar")解説 data11 の 2 番目のデータ列（１ｃｈの測定データ）に対してＡＲ（自己回帰）モデルによるスペクトル解析を行う
bispec(data11$V2)
解説 data11 の 2 番目のデータ列（１ｃｈの測定データ）に対してバイスペクトル解析を行う
autcor(data11$V2)
解説 data11 の 2 番目のデータ列（１ｃｈの測定データ）に対して自己相関の解析を行う
dev.off()

参考１）以下のようにテキストデータをコピーしてＲの画面にペーストすると１ｃｈと２ｃｈのデータ比較ができます
dev.off()par(mfrow=c(4,2))
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
plot(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
spectrum(data11$V2,method="ar")data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
autcor(data11$V2)

skip=6,

data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
plot(data11$V3)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
spectrum(data11$V3,method="ar")data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
bispec(data11$V3)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
autcor(data11$V3)

skip=6,skip=6,skip=6,
skip=6,

skip=6,
skip=6,
skip=6,

２）以下のようにテキストデータをコピーしてＲの画面にペーストするとすべての解析を連続的に行います
dev.off()par(mfrow=c(2,2))
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
plot(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
spectrum(data11$V2,method="ar")data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_12.csv",
sep=",", nrows=6000)
autcor(data11$V2)・・・・・・・・・・・・・・・・
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_24.csv",
sep=",", nrows=6000)
plot(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_24.csv",
sep=",", nrows=6000)
spectrum(data11$V2,method="ar")data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_24.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_24.csv",
sep=",", nrows=6000)
autcor(data11$V2)

skip=6,skip=6,skip=6,skip=6,skip=6,skip=6,skip=6,skip=6,３）以下のようにテキストデータをコピーしてＲの画面にペーストすると各種データの関係をグラフで確認できます
下記をコピー＆ペースト
setwd("D:/R-data")
df <- read.csv("20170905-0054/20170905-0054_03.csv",skip=2,header=T,row.names=1)
ndf <- na.omit(df)#欠損(NA)を除く
ndf
daBbb <- as.numeric(ndf[,1])#ndf[,1]を数字に変換
daCcc <- as.numeric(ndf[,2])#ndf[,2]を数字に変換

A1 <- ts(daBbb)
B1 <- ts(daCcc)
df <- read.csv("20170905-0054/20170905-0054_04.csv",skip=2,header=T,row.names=1)
ndf <- na.omit(df)#欠損(NA)を除く
ndf
daBbb <- as.numeric(ndf[,1])#ndf[,1]を数字に変換
daCcc <- as.numeric(ndf[,2])#ndf[,2]を数字に変換

A2 <- ts(daBbb)
B2 <- ts(daCcc)
df <- read.csv("20170905-0054/20170905-0054_05.csv",skip=2,header=T,row.names=1)
ndf <- na.omit(df)#欠損(NA)を除く
ndf
daBbb <- as.numeric(ndf[,1])#ndf[,1]を数字に変換
daCcc <- as.numeric(ndf[,2])#ndf[,2]を数字に変換

A3 <- ts(daBbb)
B3 <- ts(daCcc)
df <- read.csv("20170905-0054/20170905-0054_06.csv",skip=2,header=T,row.names=1)
ndf <- na.omit(df)#欠損(NA)を除く
ndf
daBbb <- as.numeric(ndf[,1])#ndf[,1]を数字に変換
daCcc <- as.numeric(ndf[,2])#ndf[,2]を数字に変換

A4 <- ts(daBbb)
B4 <- ts(daCcc)

par(mfrow=c(4,2))
ts.plot(A1,main="A1")ts.plot(B1,main="B1")
ts.plot(A2,main="A2")ts.plot(B2,main="B2")
ts.plot(A3,main="A3")ts.plot(B3,main="B3")ts.plot(A4,main="A4")ts.plot(B4,main="B4")
下記をコピー＆ペースト
par(mfrow=c(2,2))plot(A1,B1)
plot(A2,B2)
plot(A3,B3)
plot(A4,B4)

実施例：バイスペクトルの変化を確認する
dev.off()par(mfrow=c(3,2))
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_13.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_14.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_15.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_15.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_17.csv",
sep=",", nrows=6000)
bispec(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_18.csv",
sep=",", nrows=6000)
bispec(data11$V2)

skip=6,skip=6,
skip=6,skip=6,skip=6,skip=6,

実施例：自己相関の変化を確認する
dev.off()par(mfrow=c(3,2))
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_13.csv",
sep=",", nrows=6000)
autcor(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_14.csv",
sep=",", nrows=6000)
autcor(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_15.csv",
sep=",", nrows=6000)
autcor(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_15.csv",sep=",", nrows=6000)
autcor(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_17.csv",
sep=",", nrows=6000)
autcor(data11$V2)
data11 <- read.table("C:/20191220/20191220-0022/20191220-0022_18.csv",
sep=",", nrows=6000)
autcor(data11$V2)

skip=6,skip=6,
skip=6,skip=6,
skip=6,
skip=6,

解析の詳細・解析結果の解釈・・・については以下の参考書籍・・・の専門書を読んでください
参考書籍
１：統計数理
１）叩いて超音波で見る―非線形効果を利用した計測佐藤 拓宋 (著) 出版社: コロナ社 (1995/06)２）電気系の確率と統計
佐藤 拓宋 (著)
出版社: 森北出版 (1971/01)
３）不規則信号論と動特性推定
宮川 洋 (著), 佐藤拓宋 (著), 茅 陽一 (著)
出版社: コロナ社 (1969)
４）赤池情報量規準 AIC―モデリング・予測・知識発見
赤池 弘次 (著), 室田 一雄 (編さん), 土谷 隆 (編さん)
出版社: 共立出版 (2007/07)
５）ダイナミックシステムの統計的解析と制御
赤池 弘次 (著), 中川 東一郎 (著)
出版社: サイエンス社(1972)
参考資料
統計的な考え方を利用した超音波
http://ultrasonic-labo.com/?p=12202
超音波の洗浄・攪拌・加工に関する「論理モデル」
http://ultrasonic-labo.com/?p=3963
物の動きを読むhttp://ultrasonic-labo.com/?p=1074

参考＜＜超音波の音圧測定・解析＞＞
１）多変量自己回帰モデルによるフィードバック解析により超音波伝搬状態の安定性・変化について解析評価します
２）インパルス応答特性・自己相関の解析により対象物の表面状態・・に関する解析評価を行います
３）パワー寄与率の解析により超音波（周波数・出力）、形状、材質、測定条件・・データの最適化に関する解析評価を行います
４）その他（表面弾性波の伝搬）の非線形（バイスペクトル）解析により対象物の振動モードに関するダイナミック特性の解析評価を行います
この解析方法は、複雑な超音波振動のダイナミック特性を時系列データの解析手法により、超音波の測定データに適応させることで実現しています。
表面弾性波を利用した超音波制御技術
http://ultrasonic-labo.com/?p=14311
表面弾性波の利用技術
http://ultrasonic-labo.com/?p=7665
精密測定プローブ
http://ultrasonic-labo.com/?p=11267
超音波と表面弾性波（オリジナル超音波システムの開発技術）
http://ultrasonic-labo.com/?p=14264

解析ソフトについて

TIMSAC for R package
統計数理研究所

November 1, 2006TIMSAC(TIMe Series Analysis and Control program package) は，統計数理研究所で開発された時系列データの解析，予測，制御のための総合的プログラムパッケージです．・・・

TIMSAC は FORTRAN で 書 か れ た プ ロ グ ラ ム で す が ， ユ ー ザ ー が 作 成 し た
FORTRAN，C，Java のプログラムにこのライブラリをリンクすることにより，より扱い易い環境が実現されました．
バイスペクトルの解析関数

bispec() : バイスペクトルの計算
自己相関の解析関数

autcor() : 直接法による自己共分散関数の計算

３）TIMSAC for R package
http://jasp.ism.ac.jp/ism/timsac/

以上