clear; clc;
% volny monoliticky priznaticky nosnik
% data z dokumentu Alumina standard a MKP
load dlzky.txt;
load sirky.txt;
load tloustky.txt;
load hmotnosti.txt;
load emoduly.txt;
load gmoduly.txt;
load micka.txt;
load frekvencie.txt; 
frekvencie_analytika=zeros(12,1);
frekvencie_analytika_bez=zeros(12,1);
frekvencie_alumina=frekvencie(:,1); % data z Alumina
frekvencie_MKP=frekvencie(:,2); % frekvencie z MKP
teleso=[1 2 3 4 5 6 7 8 9 10 11 12];

%% ohyb

for n=1:12
    a=dlzky(n,1);        % dlzka [m]
    b=sirky(n,1);        % sirka [m]
    t=tloustky(n,1);        % tloustka [m]
    hmotnost=hmotnosti(n,1);
    rho=hmotnost/(a*b*t);      % hustota [kgm-3]
    E=emoduly(n,1);        % E [Pa]
    D=E*(b*(t^3))/12; % ohybova tuhost
    rhoA=rho*b*t; % menovatel
    frekv=omega(a,D,rhoA); % frekvencia [Hz
    frekvencie_analytika_bez(n,1)=frekv; % frekvencia bez korekcie
    korekcia=1/(sqrt(1+6.585*((t/a)^2))); % frekvencia s korekciou
    frekv=frekv*korekcia;
    frekvencie_analytika(n,1)=frekv;
end
% rozdiel medzi korigovanou a frekvenciou z clanku
rozdiel=(frekvencie_analytika-frekvencie_alumina);
% rozdiel medzi frekvenciou z clankou a frekvenciou bez korekcie
rozdiel_bez=(frekvencie_analytika_bez-frekvencie_alumina);
% relativne rozdiely [%]
rel_rozdiel=100*rozdiel./frekvencie_alumina;
rel_rozdiel_bez=100*rozdiel_bez./frekvencie_alumina;
rozdiel_MKP_al=(frekvencie_MKP-frekvencie_alumina);
rel_rozdiel_MKP_al=100*rozdiel_MKP_al./frekvencie_alumina;

figure(1)
subplot(2,1,1)
plot(teleso,frekvencie_analytika,'b*',teleso,frekvencie_analytika_bez,'r*',...
    teleso,frekvencie_alumina,'g',teleso,frekvencie_MKP,'ko')
legend('analytika s korekciou','analytika','alumina','MKP')
grid on
set(gca,'xtick',[1 2 3 4 5 6 7 8 9 10 11 12])
ylabel('Hz')
subplot(2,1,2)
plot(teleso,rel_rozdiel,'b*',teleso,rel_rozdiel_bez,'r*',...
    teleso,rel_rozdiel_MKP_al,'ko')
grid on
ylabel('%')
set(gca,'xtick',[1 2 3 4 5 6 7 8 9 10 11 12])


%% torz

alumina_torz=frekvencie(:,5);
MKP_torz=frekvencie(:,6);
rel_rozdiel_torz=frekvencie(:,8);

figure(2)
subplot(2,1,1)
plot(teleso,alumina_torz,'g',teleso,MKP_torz,'ko')
grid on
ylabel('Hz')
set(gca,'xtick',[1 2 3 4 5 6 7 8 9 10 11 12])
legend('alumina','MKP')
subplot(2,1,2)
plot(teleso,rel_rozdiel_torz,'ko')
grid on
ylabel('%')
set(gca,'xtick',[1 2 3 4 5 6 7 8 9 10 11 12])

format long

frekvencie_analytika_bez

frekvencie_analytika

frekvencie_MKP


