close all
clear all

V = (5000/2.2)/1025; 

filename = 'sp_156_eb93001';

g = 9.81;
ff = 0:.002:.58;
ff = ff';
alpha = 8.1e-3;
beta = .74;
w = 2*pi*ff;
df = ff(2)-ff(1);
 Hs = 1:2;
 clear S_Hs;
 for i = 1:length(Hs)
   U19 = sqrt(Hs(i)*g/0.21);
   w0 = g/U19;
   S_Hs(i,:) = 2*pi*(alpha*g^2)*exp(-beta*(w0./w).^4)./(w.^5);
   S_Hs(i,1) = S_Hs(i,2);
 
 
   Pb_PM(i) = BudalLimit(S_Hs(i,:)',ff,V);
 end

 


cmd = (['wc ' filename]);

[status result] = system(cmd);

if(status > 0) 
  error([cmd ' returns error']);
end

n_lines = str2num(strtok(result));
n_LinesPerRecord = 74;
n_Records = n_lines/n_LinesPerRecord;

fid = fopen(filename);

%tline = cell{n_LinesPerRecord,1};

for i = 1:n_Records
    
  for j = 1:n_LinesPerRecord
    tline{j,1} = fgetl(fid); 
  end
  disp(['Processing ' tline{1}]);
  [serialdate(i) freq{i} bw{i} energy{i} Hs(i) Tp(i)] = ParseDataWellSpectrum(tline);
  
  xi = freq{i}(1):.005:freq{i}(end);
  yi_energy = interp1(freq{i},energy{i},xi);
  
  Pb_meas(i) = BudalLimit(yi_energy,xi,V);
  
end

fclose(fid);

load('../History.mat');

for i = 1:4:n_Records-4
  print_filenm = ['WaveConditions_' datestr(serialdate(i),'mmmm_dd_yyyy_HH-MM')]
  if(~exist([print_filenm '.pdf']))  
    
    figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
    for j = 0:3
       if ((i+j) <= n_Records)
         data_filenm = ['WaveData_'  datestr(serialdate(i+j),'mmmm_dd_yyyy_HH-MM') '.mat'];
         subplot(2,2,j+1)
         plot(ff,S_Hs,freq{i+j},energy{i+j});
         xlabel('freq (Hz)');
         ylabel('Energy (m^2/Hz)');
         title(['Wave Conditions: ' datestr(serialdate(i+j))]);
         legend(['Hs = 1m  (P_b = ' num2str(Pb_PM(1),4) ' W)'],['Hs = 2m  (P_b = ' num2str(Pb_PM(2),4) ' W)'],['Measured (P_b  = ' num2str(Pb_meas(i+j),4) ' W']);   
         PxxHeave = energy{i+j};
         f = freq{i+j};
         sdate = serialdate(i+j);
         PeakPeriod = Tp(i+j);
         BudalLimit = Pb_meas(i+j);
         SigWaveHeight = Hs(i+j);
         save(data_filenm,'f','PxxHeave','sdate','PeakPeriod','BudalLimit','SigWaveHeight');
       end  
    end
  print(print_filenm,'-dpdf');
  close all;  
  end
end

first_date_idx = 19;

%last_date_idx = find(MeanSerialDate<datenum(2013,5,6,12,0,0));
%last_date_idx = last_date_idx(end);
last_date_idx = 800; %length(MeanSerialDate);

xi = MeanSerialDate(first_date_idx):(MeanSerialDate(last_date_idx)-MeanSerialDate(first_date_idx))/1000:MeanSerialDate(last_date_idx);
yi_max = interp1(serialdate,Pb_meas,xi);
yi_meas = interp1(MeanSerialDate(first_date_idx:last_date_idx),MeanRectifiedPower(first_date_idx:last_date_idx),xi);
yi_Hs = interp1(serialdate,Hs,xi);
yi_Tp = interp1(serialdate,Tp,xi);


save('../../Sept2013DeploymentSummaryData','xi','yi_max','yi_meas','yi_Hs','yi_Tp');





figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
plot(yi_Hs,yi_meas,'x');
title([datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
xlabel('Significant Wave Height (m)');
ylabel('Rectified Power (W)');
print('RectifiedPowerVsWaveHeight','-dpdf');


figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
plot(yi_Tp,yi_meas,'x');
title([datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
xlabel('Peak Period (Tp)');
ylabel('Rectified Power (W)');
print('RectifiedPowerVsPeakPeriod','-dpdf');



figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');

[AX,H1,H2] = plotyy(xi,yi_meas,xi,yi_Hs);
set(AX(1),'FontSize',20);
set(AX(2),'FontSize',20);
title(['Electrical Power and Sig. Wave Height for ' datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);

legend('Rectified Power','Sig. Wave Height');

xlabel('Date');
datetick('x','keepticks');


set(AX(1),'YLim',[0 1000]);
set(AX(1),'YTick',0:200:1000);

set(AX(2),'YLim',[0 2]);
set(AX(2),'YTick',0:.2:2);

set(get(AX(1),'Ylabel'),'String','Rectified Power (W)'); 
set(get(AX(2),'Ylabel'),'String','Sig. Wave Height (m)');
set(get(AX(1),'Ylabel'),'FontSize',20);
set(get(AX(2),'Ylabel'),'FontSize',20);
set(AX(2),'XTick',[]);

print('ElectricalPowerVsWaveHeight','-dpdf');


figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
plot(xi,yi_meas,xi,yi_max);
title(['Electrical Power for ' datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
legend('Rectified Power','Max Absorbable Power');
xlabel('Date');
ylabel('Watts');
datetick('x','keepticks');
print('ElectricalPowerComparison','-dpdf');


figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
[AX,H1,H2] = plotyy(xi,100*yi_meas./yi_max,xi,yi_Hs);
set(H1,'LineWidth',2)
set(AX(1),'FontSize',20);
set(AX(2),'FontSize',20);
title(['Performance for ' datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);

xlabel('Date');
datetick('x','keepticks');

set(AX(1),'YLim',[0 20]);
set(AX(1),'YTick',0:2:20);

set(AX(2),'YLim',[0 2]);
set(AX(2),'YTick',0:.2:2);

set(get(AX(1),'Ylabel'),'String','Percentage of Theoretical Maximum'); 
set(get(AX(2),'Ylabel'),'String','Sig. Wave Height (m)');
set(get(AX(1),'Ylabel'),'FontSize',20);
set(get(AX(2),'Ylabel'),'FontSize',20);
set(AX(2),'XTick',[]);

print('PercentageOfTheoreticalMax','-dpdf');



figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
plot(yi_Hs,100*yi_meas./yi_max,'x');
title([datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
xlabel('Significant Wave Height (m)');
ylabel('Percentage of Theoretical Maximum');
print('PercentageOfMaxVsWaveHeight','-dpdf');


figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
plot(yi_Tp,100*yi_meas./yi_max,'x');
title([datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
xlabel('Peak Period (s)');
ylabel('Percentage of Theoretical Maximum');
print('PercentageOfMaxVsPeakPeriod','-dpdf');


figure('Position',[19 45 1538 1093],'PaperPositionMode','manual','PaperType','usletter','PaperPosition',[0 0 11 8.5],'PaperOrientation','landscape');
set(gca,'FontSize',20);
scatter3(yi_Hs(10:920),yi_Tp(10:920),100*yi_meas(10:920)./yi_max(10:920));
title([datestr(MeanSerialDate(first_date_idx)) ' through ' datestr(MeanSerialDate(last_date_idx))]);
xlabel('Significant Wave Height (m)');
ylabel('Peak Period (s)');
zlabel('Percentage of Theoretical Maximum');

print('PercentageOfMaxVsPeakPeriodAndWaveHeight','-dpdf');

