clear
close all
% Sprawozdanie nr 3 fizyka - obliczenia
%%

%
%% Pomiary

% światłość żarówki
Iz=19; % [cd]

% błędy pomiarowe przyrządów
err_U = 0.1; %[V]
err_I = 5;    %[uA]
err_r = 0.1;  %[cm]

% dla stałej odległości
r1=20; %[cm]
Ur1 = [1;2;4;6;8;10;12;14;15]; %[V]
Ir1 = [4;10;22;33;45;56;68;82;85]; %[uA]

r2=40;
Ur2 = [1;2;4;6;8;10;12;14;15]; %[V]
Ir2 = [0;1;3;5;7.5;10;12;14;15.5]; % [uA]

r3=60;
Ur3 = [1;2;4;6;8;10;12;14;15]; %[V]
Ir3 = [0;0;1;2;2.5;3.5;4.5;5.5;6]; % [uA]
% liczba próbek w seriach pomiarowych
n1 = size(Ur1);

% dla stałego napięcia
U1 = 5; %[V]
Ru1 = [20;25;30;35;40;45;50]; %[cm]
Iu1 = [28;15.5;9.5;6;4.3;3;2]; %[uA]

U2 = 10; %[V]
Ru2 = [20;25;30;35;40;45;50]; %[cm]
Iu2 = [56;32;20;13.5;9.5;7;5]; %[uA]

U3 = 15; %[V]
Ru3 = [20;25;30;35;40;45;50]; %[cm]
Iu3 = [85;48;30;20.5;14.5;11;8]; %[uA]
% liczba próbek w seriach pomiarowych
n2 = size(Ru1);

% zamiana jednostek do SI

% długość z [cm] na [m]
r1 = r1/100;
r2 = r2/100;
r3 = r3/100;

Ru1 = Ru1/100;
Ru2 = Ru2/100;
Ru3 = Ru3/100;

err_r = err_r/100;

% prąd z [uA] na [A]
Ir1 = Ir1*1e-6;
Ir2 = Ir2*1e-6;
Ir3 = Ir3*1e-6;
Iu1 = Iu1*1e-6;
Iu2 = Iu2*1e-6;
Iu3 = Iu3*1e-6;

err_I = err_I*1e-6;

%% 1)

E1= Iz*(1./Ru1).^2;
E2= Iz*(1./Ru2).^2;
E3= Iz*(1./Ru3).^2;

%% 3)

err_u1 = ones(n1)*err_U/sqrt(3);
err_u2 = ones(n2)*err_U/sqrt(3);

err_i = ones(n1)*err_I/sqrt(3);

err_e1 = sqrt(((-2*Iz)*((1./Ru1).^3)*(err_r/sqrt(3))).^2);
err_e2 = sqrt(((-2*Iz)*((1./Ru2).^3)*(err_r/sqrt(3))).^2);
err_e3 = sqrt(((-2*Iz)*((1./Ru3).^3)*(err_r/sqrt(3))).^2);


%% 2)

figure();

subplot(2,1,1);
plot(E1,Iu1*1e6,'--',E2,Iu2*1e6,'--',E3,Iu3*1e6,'--');
hold on
plot(E1,Iu1*1e6,'+',E2,Iu2*1e6,'+',E3,Iu3*1e6,'+');
hold off
grid on
ylabel('I[\muA]');
xlabel('E[lx]');
legend('U=5V','U=10V','U=15V');
title('I(E)');

subplot(2,1,2);
plot(Ur1,Ir1*1e6,'--',Ur2,Ir2*1e6,'--',Ur3,Ir3*1e6,'--');
hold on
plot(Ur1,Ir1*1e6,'+',Ur2,Ir2*1e6,'+',Ur3,Ir3*1e6,'+');
grid on
ylabel('I[\muA]');
xlabel('U[V]');
xlim([0,16]);
legend('r=20cm','r=40cm','r=60cm');
title('I(U)');

%% 4)

y1 = log(Iu1);
y2 = log(Iu2);
y3 = log(Iu3);

x1 = log(E1);
x2 = log(E2);
x3 = log(E3);

coeff1 = polyfit(x1,y1,1);
coeff2 = polyfit(x2,y2,1);
coeff3 = polyfit(x3,y3,1);

a1= exp(coeff1(2)); b1= coeff1(1);
a2= exp(coeff2(2)); b2= coeff2(1);
a3= exp(coeff3(2)); b3= coeff3(1);



%% 5)


range = (x1(n2(1))- x1(1))/8;
step = range/100;
x1 = (x1(1)-range):step:(x1(n2(1))+range);
y1 = coeff1(1)*x1 + coeff1(2);

range = (x2(n2(1))- x2(1))/8;
step = range/100;
x2 = (x2(1)-range):step:(x2(n2(1))+range);
y2 = coeff2(1)*x2 + coeff2(2);

range = (x3(n2(1))- x3(1))/8;
step = range/100;
x3 = (x3(1)-range):step:(x3(n2(1))+range);
y3 = coeff3(1)*x3 + coeff3(2);

figure();

loglog(exp(x1),exp(y1)*1e6,'--',exp(x2),exp(y2)*1e6,'--',exp(x3),exp(y3)*1e6,'--');
hold on
plot(E1,Iu1*1e6,'+',E2,Iu2*1e6,'+',E3,Iu3*1e6,'+');
hold off
grid on
ylabel('I[\muA]');
xlabel('E[lx]');
xlim([50,800]);
ylim([1,200]);
legend('U=5V','U=10V','U=15V');
title('I(E)');

%% 6)


y1 = log(Iu1);
y2 = log(Iu2);
y3 = log(Iu3);

x1 = log(E1);
x2 = log(E2);
x3 = log(E3);

eps1 = y1 - (x1*coeff1(1) +coeff1(2));
eps2 = y2 - (x2*coeff2(1) +coeff2(2));
eps3 = y3 - (x3*coeff3(1) +coeff3(2));

s1=sqrt(sum(eps1.^2)/(n2(1)-2));
s2=sqrt(sum(eps2.^2)/(n2(1)-2));
s3=sqrt(sum(eps3.^2)/(n2(1)-2));

u_B1 = s1*sqrt(sum(x1.^2)/(n2(1)*sum(x1.^2)-(sum(x1))^2));
u_B2 = s2*sqrt(sum(x2.^2)/(n2(1)*sum(x2.^2)-(sum(x2))^2));
u_B3 = s3*sqrt(sum(x3.^2)/(n2(1)*sum(x3.^2)-(sum(x3))^2));

u_A1 = s1*sqrt(n2(1)/(n2(1)*sum(x1.^2)-(sum(x1))^2));
u_A2 = s2*sqrt(n2(1)/(n2(1)*sum(x2.^2)-(sum(x2))^2));
u_A3 = s3*sqrt(n2(1)/(n2(1)*sum(x3.^2)-(sum(x3))^2));


disp('wspolczynniki regresji liniowej  A i B');
disp([coeff1(1),coeff1(2)]);
disp([coeff2(1),coeff2(2)]);
disp([coeff3(1),coeff3(2)]);

disp('bledy wspolczynnikow regresji liniowej  A i B');
disp([u_A1,u_B1]);
disp([u_A2,u_B2]);
disp([u_A3,u_B3]);


%% 7)

u_b1 = u_A1;
u_b2 = u_A2;
u_b3 = u_A3;


u_a1 = sqrt((exp(coeff1(2))*u_B1)^2);
u_a2 = sqrt((exp(coeff2(2))*u_B2)^2);
u_a3 = sqrt((exp(coeff3(2))*u_B3)^2);

format long

disp('wspolczynniki a i b we wzorze I(E)');
disp([a1,b1]);
disp([a2,b2]);
disp([a3,b3]);

disp('bledy wspolczynnikow a i b');
disp([u_a1,u_b1]);
disp([u_a2,u_b2]);
disp([u_a3,u_b3]);