Fitting Curve with scattered dots

2 views (last 30 days)
Hi,
I have to plot the beam propagation of a Laser. I can plot the scattered dots, which we have measured in our experiments. Now I wanna Fit those dots. But it doesn't work. It says "Error using fit>iFit Insufficient data. You need at least 3 data points to fit this model." But I already have 22 data points. Can anybody help me?
close all;
clear variables;
%Files and Pics
data_list = dir('*.bmp');
%Zugreifen auf Ordner und BIlder und Namen zuweisen
name = {data_list.name};
%Laden der schwarzen Filterdatei
%black = uigetfile("E:\Master_Technologiemanagement\Spezialisierung\Technik\Laser\Praktikum\SB_V2\Pictures\dunkel\460_dunkel.bmp");
%Konvertieren in Matrix
imBlack = imread("E:\Master_Technologiemanagement\Spezialisierung\Technik\Laser\Praktikum\SB_V2\Pictures\dunkel\460_dunkel.bmp");
lambda = 0.633; %Wellenlänge Laser in um
x_array = [];
y_array = [];
z_array = [];
for i = 1:length(name)
picture = imread(name{i}); %Einlesen des Bildes in Matrix
current_name = name{i};
zaxe = str2double(name{i}(1:3)); %Auslesen des z Wertes aus Dateinamen
imBlack = double(imBlack);
%picture(:,:,i)=picture; %speichert das aktuelle (i-te)
%Bild als 2D Matrix in der 3D Matrix (:,:,i)
picture = double(picture);
filtered = minus(picture,imBlack);
%Ermittlung ROI
ROI = im2double(filtered);
ROIdiff = [1,1];
ROIoffset =[0, 0];
%1.Moment
[row,column] = size(filtered);
Summe = 0;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*x;
end
end
first_Moment_X = Summe;
Summe = 0;
Integral = sum(sum(filtered));
Schwerp_X = first_Moment_X./Integral;
Schwerp_X = double(Schwerp_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*y;
end
end
first_Moment_Y=Summe;
Schwerp_Y = first_Moment_Y/Integral;
%2.Moment
[row,column] = size(filtered);
Summe = 0.;
Integral = sum(sum(filtered));
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(x-Schwerp_X).^2;
end
end
sec_Moment_X = Summe;
Summe = 0;
Sigma_X = sec_Moment_X./Integral;
Sigma_X = double(Sigma_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).^2;
end
end
sec_Moment_Y=Summe;
Summe = 0;
Sigma_Y = sec_Moment_Y./Integral;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).*(x-Schwerp_X);
end
end
sec_Moment_XY=Summe;
Summe = 0;
Sigma_XY = sec_Moment_XY./Integral;
%Berechnung des Strahldurchmessers nach Formel 18 und 19 in Norm
dxPixel = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y+(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4*(Sigma_XY.^2))));
dyPixel = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y-(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4.*(Sigma_XY.^2))));
%Gesamtbild mit ROI, Schwerpunkt und Durchmesser
% imshow(picture);
axis on
hold on;
% rectangle('position',[ROIoffset(2), ROIoffset(1),2*dxPixel,2*dyPixel],'EdgeColor','b');
% plot(ROIoffset(2)+SchwerpX, ROIoffset(1)+Schwerp_Y,'r+', 'MarkerSize', 5, 'LineWidth', 1);
t=-pi:0.01:pi;
xOff=ROIoffset(2)+Schwerp_X+dxPixel.*0.5.*cos(t);
yOff=ROIoffset(1)+Schwerp_Y+dyPixel.*0.5.*sin(t);
% plot(xOff,yOff,'g')
hold on;
scatter(zaxe,dxPixel,50,'x','MarkerFacecolor','blue');
end
%Plotten der einzelnen Bilder
hold on;
[fitx,gofx]= fit(zaxe,dxPixel,"poly2");
fithandle = plot(fitx);
set(fithandle,'color','[0 .7 .7]');
[fitxw,gofxw] = fit(zaxe,dxPixel,'poly2','weights',1./dxPixel);
  5 Comments
Matt J
Matt J on 21 Apr 2022
@Berkin Bilgic which 3 lines and why can't you attach the results that those lines produce?

Sign in to comment.

Accepted Answer

Cris LaPierre
Cris LaPierre on 21 Apr 2022
There are a couple issues to fix.
The error you are currently getting is because zaxe and dxPixel are scalars. In order to use them as inputs to the fit function, they need to be column vectors. That is fixed by using indexing to capture the result of each loop and assign it to a specific row of the corresponding vector.
zaxe(i,1) = str2double(name{i}(1:3)); %Auslesen des z Wertes aus Dateinamen
...
dxPixel(i,1) = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y+(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4*(Sigma_XY.^2))));
dyPixel(i,1) = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y-(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4.*(Sigma_XY.^2))));
You will then get a new error that one of your values (in dxPixel) is NaN, which is not allowed. You can fix that with a call to prepareCurveData.
[zaxe,dxPixel] = prepareCurveData(zaxe,dxPixel);
The final code might look like this.
% Load images
unzip('bmps.zip')
%Files and Pics
data_list = dir('*.bmp');
%Zugreifen auf Ordner und BIlder und Namen zuweisen
name = {data_list.name};
%Konvertieren in Matrix
imBlack = imread("460_dunkel.bmp");
lambda = 0.633; %Wellenlänge Laser in um
x_array = [];
y_array = [];
z_array = [];
for i = 1:length(name)
picture = imread(name{i}); %Einlesen des Bildes in Matrix
current_name = name{i};
zaxe(i,1) = str2double(name{i}(1:3)); %Auslesen des z Wertes aus Dateinamen
imBlack = double(imBlack);
%picture(:,:,i)=picture; %speichert das aktuelle (i-te)
%Bild als 2D Matrix in der 3D Matrix (:,:,i)
picture = double(picture);
filtered = minus(picture,imBlack);
%Ermittlung ROI
ROI = im2double(filtered);
ROIdiff = [1,1];
ROIoffset =[0, 0];
%1.Moment
[row,column] = size(filtered);
Summe = 0;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*x;
end
end
first_Moment_X = Summe;
Summe = 0;
Integral = sum(sum(filtered));
Schwerp_X = first_Moment_X./Integral;
Schwerp_X = double(Schwerp_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*y;
end
end
first_Moment_Y=Summe;
Schwerp_Y = first_Moment_Y/Integral;
%2.Moment
[row,column] = size(filtered);
Summe = 0.;
Integral = sum(sum(filtered));
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(x-Schwerp_X).^2;
end
end
sec_Moment_X = Summe;
Summe = 0;
Sigma_X = sec_Moment_X./Integral;
Sigma_X = double(Sigma_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).^2;
end
end
sec_Moment_Y=Summe;
Summe = 0;
Sigma_Y = sec_Moment_Y./Integral;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).*(x-Schwerp_X);
end
end
sec_Moment_XY=Summe;
Summe = 0;
Sigma_XY = sec_Moment_XY./Integral;
%Berechnung des Strahldurchmessers nach Formel 18 und 19 in Norm
dxPixel(i,1) = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y+(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4*(Sigma_XY.^2))));
dyPixel(i,1) = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y-(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4.*(Sigma_XY.^2))));
%Gesamtbild mit ROI, Schwerpunkt und Durchmesser
t=-pi:0.01:pi;
xOff=ROIoffset(2)+Schwerp_X+dxPixel(i).*0.5.*cos(t);
yOff=ROIoffset(1)+Schwerp_Y+dyPixel(i).*0.5.*sin(t);
end
scatter(zaxe,dxPixel,50,'x','MarkerFacecolor','blue');
[zaxe,dxPixel] = prepareCurveData(zaxe,dxPixel);
Warning: Removing NaN and Inf from data
%Plotten der einzelnen Bilder
hold on;
[fitx,gofx]= fit(zaxe,dxPixel,"poly2");
fithandle = plot(fitx);
set(fithandle,'color','[0 .7 .7]');
[fitxw,gofxw] = fit(zaxe,dxPixel,'poly2','weights',1./dxPixel);
hold off
  1 Comment
Bersin Tekmen
Bersin Tekmen on 21 Apr 2022
Damn I am so happy right now. Thank you so much. This one worked for me. Thank you.

Sign in to comment.

More Answers (1)

Voss
Voss on 21 Apr 2022
The problem is that zaxe and dxPixel are overwritten on each iteration of the for loop, so after the loop finishes, they have the value they had on the last iteration. That's a problem for fit because they are both scalars.
If you want to use all zaxe and dxPixel values used in each iteration of the for loop you can make them into vectors with one element per for-loop iteration, and use those vectors in fit after the loop. Here I've done that and the vectors are called zaxe_all and dxPixel_all.
(Also, note that the loop iterates over all .bmp files, including 460_dunkel.bmp and 560_fokus.bmp. I got a NaN somewhere that fit didn't like, corresponding to 460_dunkel.bmp, so I removed that from name. It is unclear whether 560_fokus.bmp should also be removed.)
(Another thing I want to point out is you can use zaxe_all and dxPixel_all to create your scatter plot after the loop, rather than doing a scatter plot for each iteration inside the loop. However, you're specifying a constant MarkerFaceColor but using a Marker that's not affected by MarkerFaceColor, so I wasn't sure what the intent was, so I didn't change that.)
close all;
clear variables;
%Files and Pics
data_list = dir('*.bmp');
%Zugreifen auf Ordner und BIlder und Namen zuweisen
name = {data_list.name}
name = 1×23 cell array
{'460.bmp'} {'460_dunkel.bmp'} {'470.bmp'} {'480.bmp'} {'490.bmp'} {'500.bmp'} {'540.bmp'} {'544.bmp'} {'548.bmp'} {'552.bmp'} {'556.bmp'} {'560_fokus.bmp'} {'564.bmp'} {'568.bmp'} {'572.bmp'} {'576.bmp'} {'580.bmp'} {'582.bmp'} {'620.bmp'} {'630.bmp'} {'640.bmp'} {'650.bmp'} {'660.bmp'}
%Laden der schwarzen Filterdatei
%black = uigetfile("E:\Master_Technologiemanagement\Spezialisierung\Technik\Laser\Praktikum\SB_V2\Pictures\dunkel\460_dunkel.bmp");
%Konvertieren in Matrix
% imBlack = imread("E:\Master_Technologiemanagement\Spezialisierung\Technik\Laser\Praktikum\SB_V2\Pictures\dunkel\460_dunkel.bmp");
imBlack = imread("460_dunkel.bmp");
name(strcmp(name,'460_dunkel.bmp')) = []; % remove 460_dunkel.bmp
zaxe_all = zeros(size(name));
dxPixel_all = zeros(size(name));
lambda = 0.633; %Wellenlänge Laser in um
x_array = [];
y_array = [];
z_array = [];
for i = 1:length(name)
picture = imread(name{i}); %Einlesen des Bildes in Matrix
current_name = name{i};
zaxe = str2double(name{i}(1:3)); %Auslesen des z Wertes aus Dateinamen
imBlack = double(imBlack);
%picture(:,:,i)=picture; %speichert das aktuelle (i-te)
%Bild als 2D Matrix in der 3D Matrix (:,:,i)
picture = double(picture);
filtered = minus(picture,imBlack);
%Ermittlung ROI
ROI = im2double(filtered);
ROIdiff = [1,1];
ROIoffset =[0, 0];
%1.Moment
[row,column] = size(filtered);
Summe = 0;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*x;
end
end
first_Moment_X = Summe;
Summe = 0;
Integral = sum(sum(filtered));
Schwerp_X = first_Moment_X./Integral;
Schwerp_X = double(Schwerp_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*y;
end
end
first_Moment_Y=Summe;
Schwerp_Y = first_Moment_Y/Integral;
%2.Moment
[row,column] = size(filtered);
Summe = 0.;
Integral = sum(sum(filtered));
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(x-Schwerp_X).^2;
end
end
sec_Moment_X = Summe;
Summe = 0;
Sigma_X = sec_Moment_X./Integral;
Sigma_X = double(Sigma_X);
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).^2;
end
end
sec_Moment_Y=Summe;
Summe = 0;
Sigma_Y = sec_Moment_Y./Integral;
for y=1:row
for x=1:column
Summe = Summe + double(filtered(y,x))*(y-Schwerp_Y).*(x-Schwerp_X);
end
end
sec_Moment_XY=Summe;
Summe = 0;
Sigma_XY = sec_Moment_XY./Integral;
%Berechnung des Strahldurchmessers nach Formel 18 und 19 in Norm
dxPixel = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y+(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4*(Sigma_XY.^2))));
dyPixel = 2*sqrt(2)*sqrt(Sigma_X + Sigma_Y-(sign(Sigma_X-Sigma_Y).*sqrt((Sigma_X-Sigma_Y).^2+4.*(Sigma_XY.^2))));
%Gesamtbild mit ROI, Schwerpunkt und Durchmesser
% imshow(picture);
axis on
hold on;
% rectangle('position',[ROIoffset(2), ROIoffset(1),2*dxPixel,2*dyPixel],'EdgeColor','b');
% plot(ROIoffset(2)+SchwerpX, ROIoffset(1)+Schwerp_Y,'r+', 'MarkerSize', 5, 'LineWidth', 1);
t=-pi:0.01:pi;
xOff=ROIoffset(2)+Schwerp_X+dxPixel.*0.5.*cos(t);
yOff=ROIoffset(1)+Schwerp_Y+dyPixel.*0.5.*sin(t);
% plot(xOff,yOff,'g')
hold on;
scatter(zaxe,dxPixel,50,'x','MarkerFacecolor','blue');
zaxe_all(i) = zaxe;
dxPixel_all(i) = dxPixel;
end
%Plotten der einzelnen Bilder
hold on;
[fitx,gofx]= fit(zaxe_all(:),dxPixel_all(:),"poly2");
fithandle = plot(fitx);
set(fithandle,'color',[0 .7 .7]);
[fitxw,gofxw] = fit(zaxe_all(:),dxPixel_all(:),'poly2','weights',1./dxPixel_all);

Community Treasure Hunt

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

Start Hunting!

Translated by