% Class 10 - Root Finding with Bracketed Methods
% Coded by Nigel Reuel on 2.14.19
%
clc
close all
%
% Histogram example 
% Change # bins here
nbins = 8;
D = xlsread('Dist');
figure
H = histogram(D,nbins);
xlabel('Measured Rate');
ylabel('Frequency');
% Get normal distribution parameters
P = fitdist(D,'Normal')
% Where 
mu =P.mu;
sigma = P.sigma;
% Or you can solve this one directly
mu = mean(D);
sigma = std(D);
% Simulate 5 random measurements from this distribution
for i = 1:5
    Measure = normrnd(mu,sigma)
end
% BONUS: plot PDF on top of a normalized histogram 
% This is easily done using the 'normalized' option with histogram function
figure
h = histogram(D,nbins,'Normalization','pdf');
% Funtion for a normal distribution (to plot on top!)
X = linspace(0.01,0.02,100);
f = exp(-(X-mu).^2./(2*sigma^2))./(sigma*sqrt(2*pi));
hold on;
plot(X,f,'LineWidth',1.5)
xlabel('Measured Rate')
ylabel('Probability Density')
hold off
%
% First method for solving a root problem
% Use graphical method
cd = 0.25;
t = 4;
g = 9.81;
VFunc = @(m) sqrt(g*m/cd)*tanh(sqrt(g*cd/m)*t)-36;
fplot(VFunc,[120 200])
% Examples of ROOT finding 
% Book problem 5.5
F1 = @(x) -12-21.*x+18.*x.^2-2.75.*x.^3;
% Graphically
fplot(F1,[-1 7])
% From inspection first root is -0.415
% Using the incremental search
ns = 1000; % Number of steps (more = smaller step size)
xb = incsearch(F1,-1,7,ns)
% Bisect method
[root,fx,ea,iter]=bisect(F1,-1,0,1)
% False position method
[root,fx,ea,iter]=falsepos(F1,-1,0,0.01)
% Example 5.20 from Book
u = 1800; %m/s
mo = 160000; %kg
% Gravity is the same as above
q = 2600; %kg/s
% At what time is the upward vel = 750m/s
Vup_func = @(t) u*log(mo/(mo-q*t))-g*t
Root_Func = @(t) u*log(mo/(mo-q*t))-g*t-750
% Use bisect method
[time_sol,~,~,~]= bisect(Root_Func,10,50,1)
% Check your answer by using the time in the
% velocity function
Vup_func(time_sol)