"parfor" skips or fail with filtfilt function

2 views (last 30 days)
Hi, I am trying to use 'parfor' to speedy a loop that performs 3 steps to a simple time series (SignCorr): 1 - Detrend 2 - Linear regression 3 - Low Pass Filtering Apparently, 'parfor' skip the step 3. I can run the same code just changing 'parfor' to 'for' and the results are different. Then I noticed that the difference is the LowPass Filtering Step. I am attaching the code only with the LPF and the 'parfor' results in a generic error: "Struct contents reference from a non-struct array object.".
Any idea? Thanks!
Fs = 1/2; % sampling frequency (Hz)
tp = 180; % Time Points
% LowPass Filter design
FpassLOW = 0.1; % passband frequency (Hz)
FstopLOW = 0.1 + (Fs/tp); % stopband frequency (Hz)
ApassLOW = 1;% passband ripple (dB)
AstopLOW = 50; % stopband attenuation (dB)
hLOW = fdesign.lowpass('Fp,Fst,Ap,Ast',FpassLOW,FstopLOW,ApassLOW,AstopLOW,Fs);
HdLOW = design(hLOW, 'equiripple');
EPIreshTMP = zeros(size(finalEPIresh));
tic
parfor tt3 = 1:size(EPIresh,2) % EPIresh = 180x176004 Matrix (Time x Spatial Dimensions)
SignCorr = EPIresh(:,tt3);
MMM = mean(SignCorr);
% Low Pass Filtering
s_filt2 = repmat(s_corr,6,1); % we use the center copy, avoinding phase distortion
s_filtAB1LOW = filtfilt(HdLOW.Numerator,1,s_filt2);
s_filtAB1LOW = s_filtAB1LOW((1+2*180):(3*180));
EPIreshTMP(:,tt3) = s_filtAB1LOW(:) + MMM; %giving the mean back
end
toc
newEPIresh = EPIreshTMP;
  6 Comments
Jan
Jan on 13 Sep 2018
Edited: Jan on 13 Sep 2018

@Brunno: If you post the complete error message, we could see, which line is causing the problem. It is very useful for a discussion, if the readers do not have to guess, where the problem occurs. So please do not post just a part, which seems to be relevant, but everything printed in red to the command window.

Brunno Machado de Campos
Brunno Machado de Campos on 13 Sep 2018

@Jan: As I mentioned: attached a print of the error.

And to be crystal clear, the Line 17 is: "parfor tt3 = 1:size(EPIresh,2)"

Not sure of how this is more informative than what I posted before but anyway. Thank you.

Sign in to comment.

Accepted Answer

OCDER
OCDER on 12 Sep 2018
Looks like your object HdLOW is not being passed to the parfor correctly.
The output of design seems to be a System Object. https://www.mathworks.com/help/dsp/ref/design.html
parfor processes Objects differently than a pure structure. https://www.mathworks.com/help/distcomp/objects-and-handles-in-parfor-loops.html
Try the following:
EPIresh = rand(180,100000); % simulating my image data
Frgval = [rand(180,6),ones(180,1)]; % simulating my
% regression matrix (confounding series + column of ones
Fs = 0.5; % sampling frequency (Hz)
tp = 180; % Time Points
% LowPass Filter design
FpassLOW = 0.1; % passband frequency (Hz)
FstopLOW = 0.1 + (Fs/tp); % stopband frequency (Hz)
ApassLOW = 1; % passband ripple (dB)
AstopLOW = 50; % stopband attenuation (dB)
hLOW = fdesign.lowpass('Fp,Fst,Ap,Ast',FpassLOW,FstopLOW,ApassLOW,AstopLOW,Fs);
HdLOW = design(hLOW, 'equiripple');
EPIreshTMP = zeros(size(EPIresh));
tic
HdLOW_Numerator = HdLOW.Numerator; %<===Pull out what you need from the object
parfor tt3 = 1:size(EPIresh,2)
SignCorr = EPIresh(:,tt3);
MMM = mean(SignCorr);
% Detrend
SignCorr = detrend(SignCorr,'linear');% detrend
SignCorr = SignCorr + MMM; %giving the mean back
% Regression
b = regress(SignCorr,Frgval);
WReg = Frgval*b;
s_corr = SignCorr - WReg;
s_corr = s_corr + b(end); % giving the mean back
% Low Pass Filtering
s_filt2 = repmat(s_corr,6,1); % I put the real data in a
% sanduiche of itself, avoinding phase distortion
s_filtAB1LOW = filtfilt(HdLOW_Numerator,1,s_filt2); %<===Change this to not call the object
s_filtAB1LOW = s_filtAB1LOW((1+2*180):(3*180)); %¨taking the
% center copy of the series
EPIreshTMP(:,tt3) = s_filtAB1LOW(:)+ b(end);%giving the mean back
end
toc
newEPIresh = EPIreshTMP;
  2 Comments
Brunno Machado de Campos
Brunno Machado de Campos on 12 Sep 2018
That's worked properly! I missed that "differences" that you mentioned on for and parfor. Thank you!
Brunno Machado de Campos
Brunno Machado de Campos on 12 Sep 2018
Just to add another observation: In a 'parfor' loop, we can't even check if a structure variable is empty or not (isempty(XXX)).
Thank you all!

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!