%% Liegt eine Breit-Wigner Verteilung vor?
%  Visualisierung von Verteilungsfunktionen

% Vers. 1.1 03.2016 
% Copyright W.Schweizer
% Simulation physikalischer Systeme - Computational Physics mit MATLAB 
% 2016

clear,close all,clc
%% Berechnung der Wahrscheinlichkeit - Vergleich mit Normalverteilung
tic
load BreitWignerdata
testwerte = E(1:100);
a=0.5;%1;%0.25;1.63  % empirischer Wert
normh = @(x) exp(-a*x.^2)/sqrt(pi)*sqrt(a);   %Normalverteilung
bwh = @(x) 0.5 +1/pi*atan(x);                 %Breit-Wigner Verteilung
xt = sort(testwerte);
xintn = [];
xintbw = [];
xintbw = bwh(xt);
for xmaxi=xt
    xintn = [xintn,integral(normh,-inf,xmaxi)];
end
toc
%% Plotten der beiden Wahrscheinlicheiten
xmin = floor(min(xt));
xmax = ceil(max(xt));
fh1 = subplot(2,1,1);
fa1 = gca;
plot(xt,xintn,'-.k',xt,xintbw,'-k'),shg
xlim([xmin,xmax])
ylim([0,1.05])
ylabel('Wahrscheinlichkeit')
xlabel('Zufallsvariable')
%% Auslesen der Achsendaten
fh2 = subplot(2,1,2)
fa2 = gca
% Umskalieren der Daten
XTick = get(fa1,'XTick')
YTick = get(fa1,'Ytick')
hold on
%% 
b= (xintbw(end)-xintbw(1))/(xt(end)-xt(1));      %Berechnen der Gerade
a= xintbw(1)-b*xt(1);
gerade = [b a];
ygerade = polyval(gerade,XTick);
plot(XTick,ygerade,'k')                          %Plotten der Geraden
ygtick = bwh(polyval([1/b -a/b],YTick))     
set(fa2,'YTick',ygtick)
skalierungsfaktor = polyval(gerade,xt)./xintbw;
plot(xt,xintbw.*skalierungsfaktor,'kp'),shg      %Skalieren Breit-Wigner  
plot(xt,xintn.*skalierungsfaktor,'ko'),shg       %Skalieren normalverteilt
ylabel('Wahrscheinlichkeit')
xlabel('Zufallsvariable')
title('Y-Achse umskaliert')


