% New description:
% This script is for loading a MUCE simulation file and populating with: 1.run num.,
% 2.The surface float distance from the anchor, 3. Lagrangian position of maximum cable bending above cable floats.
% 4. bending factor above cable floats, 5. mean tension below buoy 6. Mean tension below anchor float.
% 6. Wear factor at the release below the anchor float.
% Finally, loading the data in a matrix to allow import into Excel. The array is
% named SimDataSummary
% This is a slightly modified version for testing. It adds maximum tension
% and std. dev. to the array and only loads scope 1.25 runs.

close all
clear all

% Define run parameters so proper file is loaded

           H_S = [ 2.5 2.5 3.5 3.5 4 4 4.5 4.5 ];  %meters
           PER = [ 7.9 7.9 9.3 9.3 10 10 10.6 10.6 ];  %Seconds
           SCOPES = [1.25 1.35 1.25 1.35 1.25 1.35 1.25 1.35];
          
    num = 0; %Create index for array
    
for dir_num = [1 3 5 7]
   
    suffix = ['MUCE-600lbs-200d-24mSnubber-AsBuilt',num2str(dir_num)];  
 
  path = ['./',suffix,'/']; 

    scope = SCOPES(dir_num);
    wave_per = PER(dir_num);
    wave_amp = H_S(dir_num)/2 ;
    
      
for run_num = [1 2 3 4]
    
    
    
% Load appropriate file
datafile = ['Res-',num2str(run_num),'c-',num2str(scope),'s-200d-',num2str(wave_amp),'A-',num2str(wave_per),'T-7W-8SL'];
load([path,datafile,'.mat']);
load([path,datafile,'-spec.mat']);

num = num+1; % Increment index with each pass

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 2. Calculate distance from anchor for surface float
xDis = x(end,1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 3. Calculate bend factor at nodes above cable floats
%Bend_n is the nodes where we want to perform the rf analysis, where the bending appears most severe in cable, approx. 1.9 meters past
%the end of the floated cable section.

bend_n = [first_UC_node-10:first_UC_node+40];

% First, create array x_b of the actual node numbers of whole bending section nodes of interest
% This includes whole floated section + more.
 x_b = [first_FC_node-10:first_UC_node+150]; % For whole floated section

% Create subset array of the actual Lagrangian position for bending section 
s_b = s(x_b)';

% Create array of Lagrangian distance between bending section nodes, s_diff
s_diff = diff(s_b')';

% Create array of the nodes_of_interest array numbers that corrospond to
% the actual node numbers. nodes is the cable index of nodes-of-interest.
[~, nodes_bending, ~] = intersect(nodes,x_b);

% Create array of full time resolution node to node angle differences in section
% of interest using the nodes_of_interest index
rad_diff_t = diff(phi_t(:,nodes_bending)')';

% Calculate degrees per meter by dividing deg_diff by s_diff (=node
% spacing). s_diff array is expanded to match size of deg_diff_t for
% function to work OK.
rad_m = rad_diff_t./(ones(length(rad_diff_t),1)*s_diff);

% Convert to degrees
deg_m = rad_m(:,:)*(180/pi);

% Find the full resolution (nodes-of-interest) deg_m array number that corrosponds to
% the actual node numbers we are interested in for the RF analysis, bend_n.

[~, node_b, ~] = intersect(x_b,bend_n);

% Find the Lagrangian positions that corrosponds to
% the actual node numbers bend_n.

node_s = s(bend_n);

% Create bending time history of specific nodes
bending_t = (deg_m(:,node_b));

% *********** start rainflow analysis ********************
% Find out size of the array holding the the time series of angle
% differences, bending_t
size_c = size(bending_t,2);
% Allocate space fo an array of the bending factor calculation results
bf_val = ones(1,size_c);
% Build 1 row array
for i = 1:size_c
    % For each column of bending_t, calculate the extrema over time as the
    % first step in the rainflow analysis
bending_e = sig2ext(bending_t(:,i));
    % Then do a rainflow analysis
rf = rainflow(bending_e);
    % Calculate "total bending travel" at specific node based on 2 X
    % the sum of the rainflow bending cycle amplitudes in resulting rf matrix row 1
bending_f = 2 * sum(rf(1,:));
    % Fill out the array with the calculation results for each column.
bf_val(1,i)= bending_f;
end

% Find the most extreme value in the array ABOVE THE FLOATED SECTION ONLY, and its index I
% First, find out where in bf_val the first upper cable (UC) node in bf_val
% is.
[c, fstUC, ib] = intersect(bend_n,first_UC_node);
% Actual firstUC in bf_val actually -1 because diff operation
% above shifts all index's left one place.
firstUC = fstUC-1;
% Find the location of the max value starting with the first UC node
[bfmax,I] = max(bf_val(firstUC:end));
% Calculate node # in RF bending section of maximum position
max_node = I+firstUC-1;
% Find lagrangian value for max node
bfmax_pos = node_s(max_node);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 4. Calculate mean tension below surface float
% Use the end node -1

node_f = size(T_t); % find size of tension matrix
tension_f = (T_t(:,node_f(1,2))); % use size of matrix -1 to identify node. NO! Changed to last node
% because nodes of interest was broken and short last 136 nodes.

% Calculate mean value of surface float tension
mean_tf = mean(tension_f);
% Calculate max value of surface float tension
max_tf = max(tension_f);


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 5. Calculate mean tension and wear factor at node below anchor float
% Chain node is 6 and anchor node is 7
Bend_c = 6 ;% for chain node
Bend_f = 7 ;% for anchor float

% Create array of angle differences between chain and anchor nodes
rad_diff_t = diff(phi_t(:,Bend_c:Bend_f)')' ;

% Convert Radians to degree angle between nodes
deg_t = rad_diff_t(:,:)*(180/pi);

% Calculate local extrema in the signal time history
bending_e = sig2ext(deg_t);

% Rainflow analysis
rf = rainflow(bending_e);

% For WEAR FACTOR, create tension time history of specified node
tension_t = (T_t(:,Bend_c));

% Calculate mean value of tension
mean_t = mean(tension_t);

% Calculate "wear factor" at specific node based on 2 X the sum of the rainflow bending cycle
% amplitudes in rf row 1 X mean tension

wear_f = mean_t * 2 * sum(rf(1,:));

%create index name for current data set that matches the directory and run
%number.

run_str = [num2str(dir_num),num2str(run_num)];
dir_val = str2num(run_str);

% load array
SimData(num,:) = [dir_val,xDis,bfmax_pos,bfmax,mean_tf,mean_t,wear_f];

   end
end