Parforループ内でのインデックス指定

Parforループ内である行列データをFFTする処理を行おうと考えています。

そこで以下のコードを作成したのですが、

• Fの使い方が原因でParforループを実行できません。
• Fに対する有効なインデックスは、PARFORループ内に制限されます。
• Parforループにおいて異なる方法でインデックス化されている変数Fは、潜在的に反復間の依存を引き起こします。

という警告が表示されてしまいます。初歩的なミスをしているのだと思いますが、原因がわかっておりません。

l=5;        %Parforループ動作回数
J = (1:1:l); % 合計回数

    parfor b = 1:length(J)

     wind(b,:) = B{1,b}(1:Sp).*hann(Sp)';  %窓関数で処理する

     F(b,:)=fft(wind(b,:));

     fft_out=20*log10(abs(F(b,:)));

     fbin = Sp/(4*OSR);

     Inband = (0:Sp/(2*OSR));

     signal_bins = fbin+(-(nb-1)/2:(nb-1)/2);

     noise_bins = setdiff(Inband,signal_bins);

     snr(b,:) = 10*log10(sum(abs(F(b,signal_bins+1)).^2)/sum(abs(F(b,noise_bins+1)).^2));

    end

Answer 1

[mizuki]

以下のメッセージが原因について述べています。

parforループにおいて異なる方法でインデックス化されている変数Fは、潜在的に反復間の依存を引き起こします。

変数 snr を定義している行に

 F(b, signal_bins + 1)

というコードがあります。この中の b が問題となります。for ループの場合、1番目から順番に1,2,3,...とループ内の計算が実行されますが、parfor ループの場合は1番目から順番には実行されません。

そのため、b が 3 から実行されてしまった場合、parfor を実行している MATLAB ワーカは F(3,:) の情報を保持して呼び出すことができますが、別のワーカからは呼び出すことができません。

問題を簡単にした例のコードが以下になります。

   parfor b = 1:5
       F(b,:) = b*ones(1,2);
       F(b) 　　% ここでエラー
       F(b,:)   % これはok
   end

参考: parfor ループ内の変数の分類(スライス化された変数) https://jp.mathworks.com/help/distcomp/classification-of-variables-in-parfor-loops.html

対策としては、F(b, signal_bins+1) を含むコードを別の parfor で呼び出す、ループごとに代入するようなスライス化の変数定義を避けるなどがあります。

以下だとうまくいきそうです。

parfor b = 1:5
    wind(b,:) = B{1,b}(1:Sp).*hann(Sp)';  %窓関数で処理する
    F(b,:)=fft(wind(b,:));
end

parfor b=1:5
    fft_out=20*log10(abs(F(b,:)));
    fbin = Sp/(4*OSR);
    Inband = (0:Sp/(2*OSR));
    signal_bins = fbin+(-(nb-1)/2:(nb-1)/2);
    noise_bins = setdiff(Inband,signal_bins);
    snr(b,:) = 10*log10(sum(abs(F(b,signal_bins+1)).^2)/sum(abs(F(b,noise_bins+1)).^2));
end

[Jiro Doke]

Mizukiさんが仰る通り、parfor内ではループ変数（b）をインデックスとして使う変数はすべて同じ形式のインデックスでなくてはなりません。つまり、

F(b,:) や F(b,signal_bins+1) や F(b,noise_bins+1) など異なったインデックスは使えません。

Mizukiさんの提案でも良いと思いますし、もう一つは一時変数を使う方法もあるかと思います。

l=5;
J = (1:1:l);

parfor b = 1:length(J)

      wind(b,:) = B{1,b}(1:Sp).*hann(Sp)';

      tmp=fft(wind(b,:));                 % <--- ここでは tmp に代入

      fft_out=20*log10(abs(tmp));
      fbin = Sp/(4*OSR);
      Inband = (0:Sp/(2*OSR));
      signal_bins = fbin+(-(nb-1)/2:(nb-1)/2);
      noise_bins = setdiff(Inband,signal_bins);
      snr(b,:) = 10*log10(sum(abs(tmp(signal_bins+1)).^2) ./ ...
                  sum(abs(tmp(noise_bins+1)).^2));

      F(b,:) = tmp;                       % <--- ここで F(b,:) に代入
  end

[O.E]

Mizuki さんJiro Doke さん丁寧な解説と解答ありがとうございました。

疑問が解決いたしました。