%data = dlmread('physicalData/series3/21BT3','\t');                                  %USER INPUT REQUIRED - Specify the file containing the data to be analysed.
data = dlmread('Series2_EWTEC/release_new/release_Series2/1BT2.txt','\t',1,0);

%%%% THE FOLLOWING LINE MAY CAUSE PROBLEMS %%%% comment out if not required
%data = data(1:12800,:); %Due to non-unique values the physical data needs to be truncated in advance

columns = 5; %Columns of input data. If using the submission files for the blind test this should be 5.

%% USER SHOULD NOT NEED TO EDIT CODE AFTER THIS POINT %%
startWindowTime = 35.3; %Start time of the data we're interested in, relative to the experiment.
endWindowTime = 50.3;   %End time of the data we're interested in, relative to the experiment.

interpolatedFs = 4;     %Sample frequency upon which the data will be interpolated (only interested in frequency content up to 2Hz so 4Hz sampling is all that's required).

ti = (0:(1/interpolatedFs):(endWindowTime - startWindowTime)) + startWindowTime;    %Creating a new time vector for the desired time window and sampling frequency.
windowData = zeros(length(ti),columns);                                             %Predefining the size of a new array for the windowed data.
windowData(:,1) = ti';                                                              %Assigning the new time vector to the first column of the new windowed data array.

zeroOffset = zeros(1,columns);                                                      %Predefining a vector for the zero offsets of the data.

padding = 10;                                               %To achieve a smoother spectrum in frequency space the windowed data is padded with zeros to produce a longer time series.
paddedData = zeros(length(windowData)*padding,columns);     %The duration of the padded data is a multiple of the duration of the windowed data, given by the variable 'padding'.

for j = 2:columns
    windowData(:,j) = interp1(data(:,1),data(:,j),windowData(:,1));                                 %Populating the windowed data array for all 4 assessment varaibles.
    zeroOffset(1,j) = windowData(1,j);                                                              %Populating the zero offsets for all 4 assessment varaibles (really only necessary for the mooring load).
    paddedData(length(paddedData)-length(windowData(:,j))+1:end,j)= windowData(:,j)-zeroOffset(1,j);%Populating the padded data array for all 4 assessment varaibles.
end

N = length(paddedData(:,1));                    %Assigning the variables required for the FFT analysis
f = interpolatedFs*(0:(N/2))/N;                 %Specifying the frequency vector
Y = zeros(length(paddedData),columns);          %Predefining an array for the FFT output for each assessment variable.
varDensity = zeros(length(1:N/2+1),columns);    %Predefining an array for the variance density of each assessment variable.

for j = 2:columns
    Y(:,j) = fft(paddedData(:,j));                      %Performing an FFT of the windowed, zeroed and padded data
    P2 = (abs(Y(:,j))/interpolatedFs).^2;               %Calculating the variance density for double-sided spectrum
    varDensity(:,j) = P2(1:N/2+1);
    varDensity(2:end-1,j) = 2*varDensity(2:end-1,j);    %Calculating the variance density for single-sided spectrum
end

for j = 2:columns                                       %Ploting the results for all coloumns of the input data
    plot(f,varDensity(:,j)) 
    title('Single-side variance density spectrum of ...')
    xlabel('f (Hz)')
    ylabel('Psi(f)')
    hold on
end
%%