%% Small Deflection Blade Moments
% Function to calculate blade bending moment distributions with small
% deflection assumptions (not accounting for blade deflection).
% Add option to include saddle induced bending moments

function [moments] = SmallDefMoments(start,finish)

% Calculating initial actuator angles
inact_alpha = atand((start.act_in(1,3) - start.act_in(2,3)) / ...
    (start.act_in(1,1) - start.act_in(2,1))) *-1;

outact_alpha = atand((start.act_out(1,3) - start.act_out(2,3)) / ...
    (start.act_out(1,1) - start.act_out(2,1))) *-1;

% Calculating vertical forces
F_in = finish.F_in;
F_out = finish.F_out;
Fz_in = F_in * cosd(inact_alpha);
Fz_out = F_out * cosd(outact_alpha);

% Calculating spanwise bending moments
momentBend = [0, (F_in * start.act_in(2,3) * 0.0044482 +...
    F_out * start.act_out(2,3) * 0.004482)];   % Test stand BM
count = 1; %starting counter for adding inboard actuator moment point
for i = 1:length(start.blade)
    if start.blade(i,3) < start.act_in(2,3)
        momentBend(i+1,:) = [start.blade(i,3), (F_in * (start.act_in(2,3)...
            - start.blade(i,3)) * 0.0044482 + F_out * (start.act_out(2,3)...
            - start.blade(i,3)) * 0.004482)];
        count = count + 1;
    else
        momentBend(i+1,:) = [start.blade(i,3), (F_out * (start.act_out(2,3)...
            - start.blade(i,3)) * 0.004482)];
    end
end

% Adding in bending moment at actuator points
in_mom = [start.act_in(2,3), (F_out * (start.act_out(2,3) - start.act_in(2,3))...
    * 0.004482)];
momentBend = [momentBend(1:count,:); in_mom; momentBend(count+1:end,:)];
momentBend(end+1,:) = [start.act_out(2,3), 0];
momentBend(:,2) = momentBend(:,2) / 1000;

% User option to include saddle induced bending moments at maximum force
x = input('Include saddle induced moments? (Yes = 1, No = 0) ');
momentT = momentBend;
indI = momentBend(:,1);
indO = momentBend(:,2);
if x == 1
    % Calculating final inboard actuator and saddle angles
    inact_tau = atand((finish.act_in(1,3) - finish.act_in(2,3)) / ...
        (finish.act_in(1,1) - finish.act_in(2,1))) *-1;
    outact_tau = atand((finish.act_out(1,3) - finish.act_out(2,3)) / ...
        (finish.act_out(1,1) - finish.act_out(2,1))) *-1;
    
    insad_tau1 = atand((finish.saddle_in(2,1) - finish.saddle_in(1,1)) / ...
        (finish.saddle_in(2,3) - finish.saddle_in(1,3)));
    insad_tau2 = atand((finish.saddle_in(3,1) - finish.saddle_in(4,1)) / ...
        (finish.saddle_in(3,3) - finish.saddle_in(4,3)));
    insad_tau = (insad_tau1 + insad_tau2) / 2;
    
    outsad_tau1 = atand((finish.saddle_out(2,1) - finish.saddle_out(1,1)) / ...
        (finish.saddle_out(2,3) - finish.saddle_out(1,3)));
    outsad_tau2 = atand((finish.saddle_out(3,1) - finish.saddle_out(4,1)) / ...
        (finish.saddle_out(3,3) - finish.saddle_out(4,3)));
    outsad_tau = (outsad_tau1 + outsad_tau2) / 2;
    
    % Calculating induced moment contributions
    MindI = finish.saddleZ_in * finish.F_in * sind(insad_tau - inact_tau) * 0.004482 * 0.001;
    MindO = finish.saddleZ_out * finish.F_out * sind(outsad_tau - outact_tau) * 0.004482 * 0.001;
    
    % Adding to calculated bending moments
    % Outer first
    for i = 1:length(momentBend)
        momentT(i,2) = momentBend(i,2) + MindO;
        indO(i,2) = MindO;
    end
    % Inner induced moment added to blade inboard of actuator and
    % subtracted at outboard portions
    for i = 1:length(momentBend)
        if momentBend(i,1) < start.act_in(2,3)
            momentT(i,2) = momentBend(i,2) + MindI;
            indI(i,2) = MindI;
        else
            momentT(i,2) - momentBend(i,2) - MindI;
            indI(i,2) = -MindI;
        end
    end    
else
end

hold on
plot(momentT(:,1),momentT(:,2))

moments.Total = momentT;
moments.Bending = momentBend;
moments.InducedOuter = indO;
moments.InducedInner = indI;
        




end