Fit a 3D Array allong 3 rd dimension

My question is ralated to this one: https://de.mathworks.com/matlabcentral/answers/45680-faster-method-for-polyfit-along-3rd-dimension-of-a-large-3d-matrix#answer_55902

I have a series of images of the size 256,256,100.

I need to do a linear fit along the 3rd dimension (100 images in this series) in order to fit for the relative change at each pixel along the 100 images.

I need to get the fitted series and the fit parameters (p1, p2) for a linear fit. Each pixel is supposed to decrease in along these 100 images in a linear way, but the slope can be different between the pixels.

the above mentioned answer will give me the complete corrected series but not the fit parameters.

Thanks!

Answer 1

[Thomas Rosin]

Try this function

temp=your3Dmatrix;
for i=1:256
    for j=1:256
        datafit=temp(i,j,:);
        datafit=datafit(:)';
        p=polyfit(1:100,datafit,1);
        intercept(i,j)=p(2);
        slope(i,j)=p(1);
    end
end

You should get a 2d matrix of p1 and p2.

Answer 2

[Martin Offterdinger]

Thanks, Thomas. Your solution works (I have tried it already previously), but the problem is, that is very slow. Therefore the solution, which i mentioned is very nice, albeit it does not tell me the fit coefficients.

[Thomas Rosin]

the fit coefficient are in the variable 'a' so doing this should work

dataCube = rand(1000,1000,300);
x = rand(300,1);
sizeCube = size(dataCube);
p1 = zeros(sizeCube(1),sizeCube(2));
p2 = zeros(sizeCube(1),sizeCube(2));
for ii = 1:sizeCube(1);
   for iii = 1:sizeCube(2);
        signal = squeeze(dataCube(ii,iii,:));
        a = polyfit(x,signal,z)
        p1(ii,iii) = a(1);
        p2(ii,iii) = a(2);
     end
end

The reason why this code is faster is because the matrix size is preset.

Answer 3

[Nikolai Khokhlov]

Hi! I have the same task, but with quadratic fit parameters. It looks the math trick from your link is not so helpful to get the parameters of the fits, unfortunately. But, Matlab has the parallel tollbox. For me, it is twice faster with parfor loop. May be there is 2D solution, I am just at start.

%just activate parallel tollbox
parfor ii=1:2 mm = ii; end

% matrix3D is 3D matrix/cube
angle_0 = zeros(s1,s2); %positions of parabolas' centers
for xx = 1:s1
    parfor yy = 1:s2
        pixelValue = squeeze( matrix3D(xx,yy,:));
        p = polyfit(angles, pixelValue, 2);
        angle_0(xx,yy) = -p(2)/2/p(1); %get axis of parabola
    end
end

But I think, final step will be GPU. Particularly, I'm working on the experiment like in Article. They say, GPU is 100 times faster, then CPU...

Answer 4

[yanqi liu]

sir，may be use the isosurface to make 3D data

Answer 5

[Nikolai Khokhlov]

I made a next iteration of my previoius answer with parfor parallel calculations. Now it works 5-6 times faster compare to double for loop. The idea is to work with 1D array as parfor works with 1D. I'm working wIth the frames from the files, thus the loading 2D matrixes and rearrange them to 1D array is the first step. The second one is to use parallel computing. The last step is to rearrange the array back to 2D

dir_name = 'd:\DataDir\';

file_name = '*deg.dat'; % my data files have difefrent angles in degrees at the end, like this: Test_1.50deg.dat / Test_-2.00deg.dat

f_n = dir([dir_name file_name]);

[~,idx] = sort([f_n.datenum]); % sort file by time as I get the measurements of angles and save them consistently in time

f_n = f_n(idx);

fileNames = {f_n.name}; %get array of files' names

angles = zeros(length(fileNames),1);

temp = load([dir_name fileNames{1}]);

s1 = size(temp,1);% get the resolution

s2 = size(temp,2);

matrix2D = zeros(s1*s2,length(fileNames)); %zeros to preset the matrix and, thus, to make code faster

% activate parallel tolbox (may be there is another 'right' method)

parfor ii=1:2 mm = ii; end

parfor dd = 1:length(fileNames)

%extracting array of angles from files' names

fileN = fileNames{dd};

k1 = strfind(fileN,'_');

k2 = strfind(fileN,'deg');

angles(dd) = str2double(fileN(k1+1:k2(1)-1)); %get array of angles

%get 2D image

TwoDArray = load([dir_name fileN]);

%convert image to 1D array

matrix2D(:,dd) = reshape(TwoDArray.',1,[]); %reshape 2D matrix to 1D array

end;

angle_0 = zeros(s1*s2,1);

parfor xx = 1:s1*s2

pixelValue = matrix2D(xx,:);

p = polyfit(angles, pixelValue', 2); %make a fit

angle_0(xx) = -p(2)/2/p(1); %get axis of parabola from parameters of the fit to 1D array

end;

angle_0 = flipud(rot90(reshape(angle_0,s2,s1))); %reshape 1D array to 2D matrix back