MATLABのデータ抽出についておたずねしたいことがあります。

質問

MATLABのデータ抽出についておたずねしたいことがあります。
MATLAB初心者ですが、以下のようなプログラムを現在書いています。

f = input('file name? ---> ','s');
m = csvread(f,66,2,[66,2,1089,3]);
time=(0:0.002:2.046);
t=reshape(time,1024,1);
current=m(:,1)*1000;
voltage=m(:,2)*10;
m2=[t current voltage];
a=m2(1,2);
b=m2(2:40,2);
if a>max(b);
g=a;
end
for I=2:1024;
for J=1:1024;
for K=39+I;
if K>1024;
K=1024;
end
c=m2(I,2);
d=m2(J:I-1,2);
e=m2(I+1:K,2);
if c>max(e)&&c>max(d);
g=c;
else
end
end
end
end
plotyy(t,current,t,voltage)
ylim([0,700])

ifの条件を満たしたときに、m2の2列目の値だけでなく条件を満たしたその行の3列全てのデータを抽出し、条件を満たした行だけの新たな行列を作りたいのですが、どのような命令を書けばよいでしょうか。

お手数をおかけしますが、ご教授よろしくお願いいたします。

回答No.5

Kules · ベストアンサー

Kules part5です。

>この点についてですが、私の説明不足で申し訳ありません。
例えば、k=60で60番目と番目が同じだった場合、75番目の値を60番目に入れるのではなく、60番目の値は抽出せずに、75番目の値だけを抽出したいということだったのです。

なるほど～ということは、
10,30,50,90番目が同じ値で全体的にも最大値であるとすると、
k=10…p=[10,30]になるので何もしない
k=30…p=[10,30,50]になるので何もしない
k=50…p=[30,50]なので、50番目を抽出
k=90…p=90なので、90番目を抽出
という感じで動作させればいいんですかね？
これならわざわざIdxを使わなくてもlogicalで対処できそうですね。
ということで書いてみます。A No.2とNo.3を混ぜたものになりそうです。
f = input('file name? ---> ','s');
m = csvread(f,66,2,[66,2,1089,3]);
t=(0:0.002:2.046)';
m2=[t,m];
m2(:,2)=m2(:,2)*1000;
m2(:,3)=m2(:,3)*10;
N=length(t);
F=false(N,1);
for k=1:N;
Start=max([1,k-30]);
End=min([k+31,N]);
p=find(m2(:,2)==max(tmp));
if k==max(p(p<End));
F(k)=true;
end;
end;
m3=m2(F,:);

そんなに解説が必要な感じではありませんね。
先ほどはany(p==k)でpのどれか１つでもkと一致していれば
pの最大値をIdxに代入していましたが、
今回はpのEnd以下の最大のものとkが一致した時のみ
Fをtrueにしています。

こういうことでいいのですかね？
違ってたらor何か思った通りに動かなかったらor何かエラーが出たら
またその旨補足願います。

参考になれば幸いです。

お礼

回答ありがとうございました。今回も何の問題もなく動作し、所望の結果を得ることができました。
この度は何度も回答いただき本当にありがとうございます。

このプログラムはまだ先を作らなければならないので、また何か分からない点がありましたら質問させていただくことがあるかもしれません。
そのときは、また知識をお貸しいただけると大変助かります。

ありがとうございました。

回答No.4

Kules

さらにさらにKulesです。
やっぱりミスが出ましたね(笑)
ちゃんと動作確認せずに書くとこういうことが起こるので注意が必要ですね^^;

原因はわかりました。
p=find(tmp==max(tmp));→p=find(m2(:,2)==max(tmp));

Idx(k)=max(p);→Idx(k)=max(p(p<=End));
にそれぞれ変えて下さい。

どういうことかと言えば、
元々30番目と45番目が同じ値でその周辺における最大値の時、
k=30で45、k=45でも45をIdxに代入して欲しかったんですが、
pをtmp内でのインデックスにしてしまったため、
本来はkの値に関わらずp=[30,45]という値が出て欲しいところなのに
kの値を変えるとpの値も変わるようになってしまっていたということです。
(k=45の時、15～76番目でtmpを作っているため、30番目はtmpの中では16番目,45番目はtmpの中では31番目)
ということはp=[16,31]と出るのでk=40と一致するわけがありません。
そのため、pはあくまでも元の配列のインデックスから取りたいので
tmp==　と比較している箇所をm2(:,2)==　に変更しました。
ところがこれの場合、今度はpの探索範囲が広すぎて、本来出してはいけない値を出すことがあります。
(例えば30，45，90が同じ値の場合、k=45では15～76番目から探しているので答えは45だが、
単純にmax(p)とすると90も入ってしまう)
そのため、pはEnd以下の値から探すように
max(p)ではなくmax(p(p<End))としています。

お手数掛けました。もう少し別の書き方もできますが、変更箇所が
一番少なそうな方法で修正しています。

参考になれば幸いです。

補足

回答ありがとうございます。今回もエラーなく動作しました。
毎回本当にありがとうございます。

大変申し上げにくいのですが、もう一点お願いできますでしょうか。

＞tが遅い方を取るということは、例えばk=60で60番目と75番目が同じだった場合は75番目の値Current,Voltageの値を60番目に入れるということでしょうか？

この点についてですが、私の説明不足で申し訳ありません。
例えば、k=60で60番目と番目が同じだった場合、75番目の値を60番目に入れるのではなく、60番目の値は抽出せずに、75番目の値だけを抽出したいということだったのです。

この場合、今書かれているプログラムから、m3の行列からtが同じ行を削除するような命令文を書けばよいのでしょうか。

回答No.3

Kules

三度Kulesです。

>所望の動作を確認できました。

ちゃんと動いたとのことで何よりです。

>この点についてですが、実はExcelではfalseになるようにしていますが、その後私が手を加えてtが遅い方のデータを手動で記入していました。

ｔが遅い方を取るということは、例えばk=60で60番目と75番目が同じだった場合は75番目の値Current,Voltageの値を60番目に入れるということでしょうか？
そうなるとかなりめんどくさいですね…単純にfalse/trueでは分岐できなさそうです。
（配列の並び上もっといい方法はあると思いますが、とりあえずプログラムでしたいことをそのままMatlabでのプログラムに起こしたいと思います）
ということで書いてみます。

f = input('file name? ---> ','s');
m = csvread(f,66,2,[66,2,1089,3]);
t=(0:0.002:2.046)';
m2=[t,m];
m2(:,2)=m2(:,2)*1000;
m2(:,3)=m2(:,3)*10;
N=length(t);
Idx=NaN(N,1);
for k=1:N;
Start=max([1,k-30]);
End=min([k+31,N]);
tmp=m2(Start:End,2);
p=find(tmp==max(tmp));
if any(p==k);
Idx(k)=max(p);
end;
end;
Idx(isnan(Idx))=[];
m3=m2(Idx,:);

以下は解説です。
・Idx=NaN(N,1);
確かdoubleになりますが、全ての値をNaN(Not a Number,Excelにもありますよね？)にします。

・tmp=m2(Start:End,2);
p=find(tmp==max(tmp));
tmpでStartからEndに相当するCurrentの配列を取りだし、
pにはtmpの値がtmpの最大値と一致するインデックスを取りだしています。

・any(p==k);
Idx(k)=max(p);
kはインデックスに相当しますが、先ほど取りだしたpの中に１つでもkと一致するものがあった場合、
Idx(k)をNaNからpの最大値(つまりtが遅いもの)に書き変えています。pが１つしかなかった場合は
その値がそのまま入ります。

・Idx(isnan(Idx))=[];
値が入らずNaNのままになっているものを取り除きます。

こんな感じでどうでしょう？
前回もそうなんですが、実際Matlab上で動作させずに書いているので動く保証はありません(笑)
書き方もおそらくKules独特のものですし…
何かエラーが出て動かない部分があればその旨補足願います。

参考になれば幸いです。

補足

回答ありがとうございます。お手数をおかけしまして大変申し訳ありません。

いただいたプログラムを実行したところ、エラーなく動作しましたが、
一つ目のピーク値しか検出されず、それ以降のピークがm3に入っていません。
前回のプログラムでは全てのピーク値が検出されたのですが…。

原因はどこだと考えられますか(汗

回答No.2

Kules

再びKulesです。

なるほど…Excelでの書き方を見ると、
H68の値はそれより前30点より大きく、かつ後ろ31点より大きい時H68そのもので、そうでない時は空白
ということでしょうか？
ということは、雰囲気としてはcurrentが振動していて、前後62点の中での極大値になるところを探しにいっているという感じですかね？

なるほど…ではそれを踏まえて書いてみます。
これを決めるのに関係なさそうなものは省略します。
f = input('file name? ---> ','s');
m = csvread(f,66,2,[66,2,1089,3]);
t=(0:0.002:2.046)';
m2=[t,m];
m2(:,2)=m2(:,2)*1000;
m2(:,3)=m2(:,3)*10;
N=length(t);
F=false(N,1);
for k=1:N;
Start=max([1,k-30]);
End=min([k+31,N]);
if m2(k,2)==max(m2(Start:End,2));
F(k)=true;
end;
end;
m3=m2(F,:);

m3が抽出した行列です。
以下は解説です。
・t=(0:0.002:2.046)';
書かれたプログラムでは1行1024列の行列を作ってからreshapeで1024行1列の行列にしていますが、
これは転置すれば済む話なので'１つにしています。
・F=false(N,1);
logical型の変数を使っています。1024行1列false(偽)が並んだ状態です。
・Start=max([1,k-30]);
配列を探す初めのインデックスを決めています。ループ変数のkが小さいと30個前を探すことができませんので(k-30が0以下になってしまうので)、その時は1になるようにしています(1と0以下の数なら大きいのは1、kが31以上ならk-30の方が大きくなる)
・End=min([k+31,N]);
先ほどの逆です。minを使っているのはインデックスが要素数(=N)を超えないようにするための細工です。
・if m2(k,2)==max(m2(Start:End,2));
F(k)=true;
end;
currentのk番目の要素がStart～Endの中での最大値と等しい時、Fのk番目の要素をtrueにする処理です。
・m3=m2(F,:);
m2の中で、Fがtrueになっている行だけ全て抜き出すような形です。

これで大体は所望の動作をします。ただ、これだと本当にしたいのと少し違う動作になるかも知れません。
というのは、例えばk=60の時30～91を探すことになると思いますが、仮に60番目と35番目が同じ値で、しかもそれが30～91の中での最大値だった場合、質問者様の書かれたExcelの関数ではfalseですが、私が書いたものではtrueになります。
これを避けたいというのであれば、少しメンドクサイ書き方になりますが、forループの中身を
Start1=max([1,k-30]);
End1=max([1,k-1]);
Start2=min([k+1,N]);
End2=min([k+31,N]);
if m2(k,2)>max(m2([Start1:End1,Start2:End2],2));
F(k)=true;
end;

とすれば同じ動作になります。
Matlab独特の(というかKuled独特の？)書き方になっているのでかなりわかりにくいかも知れません。
わからない部分orエラーが出て動かない部分があればその旨補足願います。

参考になれば幸いです。

補足

何度も回答ありがとうございます。

Kulesさんのプログラムを実行したところ、所望の動作を確認できました。
この度は本当に助かりました。プログラムの内容もわかりやすい説明をいただき、把握することができました。

また、最後に補足をいただいた、

＞例えばk=60の時30～91を探すことになると思いますが、仮に60番目と35番目が同じ値で、しかもそれが30～91の中での最大値だった場合、質問者様の書かれたExcelの関数ではfalseですが、私が書いたものではtrueになります

この点についてですが、実はExcelではfalseになるようにしていますが、その後私が手を加えてtが遅い方のデータを手動で記入していました。ですのでMATLABのプログラムでもそのようにしたいのですが、プログラムが思いつかず困っています。

何度も申し訳ありませんが、なにかいい方法はないでしょうか…。

回答No.1

Kules

ええと…プログラム自体にエラーの元となることが書かれていますし、
最終目的が
plotyy(t,current,t,voltage)
ylim([0,700])
だとすると
この3重のforループは結局何もしていないことになるのですが。
（最終的に使っているt,current,voltageは何も変化していない）
具体的にどのようなものを出そうとしているのか、
もう少し補足願います。
特に
・プログラム中に出てくるgは何か？
について。m2からいくらかの成分を抽出して、その値の大小によってgにcを代入したりしなかったりしていますが、このgがどこにも使われていないため意味のない構文になっています。
gにどういう役割を持たせたいのでしょうか？

以上よろしくお願い致します。

補足

回答ありがとうございます。質問が分かりづらくてすいません。

作成中のプログラムをそのまま載せてしまいましたので分かりづらくなってしまいました。

plotyy以下は今は無視していただいて構いません。
また、if文の条件を満たした場合の実行文g=c;も、どのように命令文を書けばいいの分からず、仮に書いておいたものです…。説明をしませんで申し訳ありません。

今、最終的にやりたいことは、m2の2列目(current)の中でのピーク値の抽出です。
Excelでは以下のような数式を用いて判別していました。

=IF(H68>MAX(H38:H67),(IF(H68>MAX(H69:H99),H68,"")),"")

上式のH列がm2の2列目(current)と考えてください。
この条件に当てはまった複数の行を抽出し、新たな行列を作成したいと考えています。
(2列目(current)だけでなく、1列目(t)、3列目(voltage )全て抽出したいのです)
この新たな行列の作成の仕方が分からず、g=cと、仮に書いていたのです。

と、このような説明で大丈夫でしょうか…。

ちなみに、MATLAB2010ではfindpeaksという関数もあるようですが、私の使っているものは2009bなので、使えないみたいでした。