ENERGY SPECTRUM MONTHLY PLOT

Energy spectrum data - derived from a fast fourier transform of approximately 30 minutes of wave data - is plotted with frequency on the x axis and time on the y axis. The elevation of a "mountain" corresponds to the energy density, with the scale to the right of the plot. Every fourth data run is plotted with one of four colors.

This plot is useful for observing the evolution of wave energy over all the measured frequencies during the month.

Matlab code for generating the plot (no longer used)

%-------------------------------------------------------------------------------
% Name:
%   .dbase/web_programs/plotpro/matlab_scripts/mountain.m
% Description:
%   Creates the monthly spectral energy plot (mountains).
%   A few lines are pre-pended to this script.
% Called by:
%   .f90/plot_utils.make_plot
%-------------------------------------------------------------------------------

path(path,'/project/dbase/web_programs/plotpro/matlab_scripts');

wkg = wkg_path;


% LOAD DATA, GIVE SHORT NAMES FOR LOADED VARIABLES
% NOTE: x = frequency data, y = energy density data

eval(['freq=load(''',wkg,'freq.',pid,'01'');']);
eval(['ener=load(''',wkg,'ener.',pid,'01'');']);

x = freq(1,:);
y = ener;

dt = read_date_file([pid,'01'],wkg);
tim = (dt.day-1)*24 + (dt.hour + dt.min/60); % in hours

yearchar = num2str(dt.year(1));
monthnum = dt.month(1);

clear freq ener date


% CHECK TO SEE IF THERE IS ANY DATA

if ~ ( any(x>0) & any(y>0) )
  clear all;
  quit;
end


% SET FIGURE WINDOW

figure('units','normal','paperunits','normal','position',...
        [.1,.1,.5,.8],'paperposition',[.05,.05,.9,.93],'visible','off');

sz = size(y);
map = [ 0 0 0; 1 0 0; 0 1 0; 0 0 1; .5 .5 .5];
colormap(map);
cd = colormap;


% TAG FIRST RUN OF EACH DAY

pday = 0;
for i=1:sz(1)
  day = floor(tim(i)/24) + 1;
  if day ~= pday
    nday(i) = 1;
    pday = day;
  else
    nday(i) = 0;
  end
end


% SCALE Y VALUES

yscale = 7;
ym = max(max(y));
% Fixed scale: set ym=50.;
%ym=50.;
y = yscale*24*y./ym;

logx = log(x);


% PLOT DATA

eps = 5;
prvtime = 790;
j = 0;
for i=sz(1):-1:1
  if ( (prvtime - tim(i) > eps) | nday(i) == 1) 
    patch('xdata',logx,'ydata',tim(i)+y(i,:),'edgecolor','none','facecolor',cd(1+rem(j,5),:));
    j = j + 1;
    prvtime = tim(i);
  end
end

yaxisheight = .76;
set(gca,'position',[0.1 0.1 0.7 yaxisheight]);


% GET STATION INFO

info = read_info_file([pid,'01'],wkg);


% IDENTIFY IF SURGE BASIN OR STANDARD PROCESSING

surge = 0;
basin = 0;
standard = 0;

if ~isempty( findstr(info.stream_name,'SURGE') )
  surge = 1;
elseif ~isempty( findstr(info.stream_name,'BASIN') )
  basin = 1;
else
  standard = 1;
end

% WRITE STATION INFO ON PLOT

if ( surge == 1 )
  xtext=log(1/128);
elseif ( basin == 1 )
  xtext=log(1/22);
else
  xtext=log(.10);
end

ytext= 860;
txtspace = 24;

text(xtext,ytext+3*txtspace,[char(info.stnid),' ',char(info.stname)],...
  'horizontalalignment','center','fontsize',14);

stream_name = ['(',char(info.stream_name),')'];
if ( length(stream_name) > 5 )
  text(xtext, ytext+2*txtspace, stream_name, 'horizontalalignment',...
  'center','fontsize',12);
end

mname = get_monthname(monthnum,'long');
text(xtext,ytext + -0.3*txtspace,[char(mname),' ',yearchar],'horizontalalignment',...
'center','fontsize',14);

text(xtext,ytext+0.7*txtspace,'ENERGY SPECTRUM','horizontalalignment','center',...
     'fontsize',14);

% MAKE Y AXIS LABELS

ylim = 824;
set(gca,'ylim',[0,ylim],'ytick',[0:24:720],'yticklabel',num2str([1:31]'));
ylabel('Day of Month (UTC)','fontsize',12);


% MAKE X AXIS LABELS

numticks = 5;

if ( surge == 1 )
  lxlim = 1/1024; % (Hz)
  uxlim = 1/16;   % (Hz)
  f_sig_digs = 2;
  p_sig_digs = 4;
  f_labels = 0.20;
elseif ( basin == 1 )
  lxlim = 1/128; % (Hz)
  uxlim = 1/4;   % (Hz)
  f_sig_digs = 2;
  p_sig_digs = 3;
  f_labels = 0.20;
else
  lxlim = 1/25; % (Hz) 
  uxlim = 1/4; % (Hz)
  f_sig_digs = 2;
  p_sig_digs = 2;
  f_labels = 0.07;
end

% note: xtl = "x tick label"

diff = (log(uxlim) - log(lxlim))/(numticks-1);
xt = log(lxlim) + (0:(numticks-1))*diff;
xt(numticks) = log(uxlim);
fxtl = num2str(exp(xt'),f_sig_digs);
pxtl = num2str(1./exp(xt'),p_sig_digs);
set(gca,'xlim',[log(lxlim),log(uxlim)],'xtick',xt,'xticklabel',fxtl);
hold on;
plot([log(lxlim) log(lxlim)],[0 840],'k-');
box on; 
text(xt,-35*ones(1,numticks),pxtl,'horizontalalignment','center');
text(log(uxlim)+f_labels,-15,'Frequency (Hz)','horizontalalignment','left',...
     'fontsize',10);
text(log(uxlim)+f_labels,-35,'Period (s)','horizontalalignment','left',...
     'fontsize',10);


% PUT SCALE ON PLOT

if ym < .1
  rndval = .01;
elseif ym < 1
  rndval = .1;
elseif ym < 10
  rndval = 1;
else
  rndval = 10;
end

rnd = rndval*round(ym/rndval)
yscaleheight = (yaxisheight/ylim) * (rnd/ym)*yscale*24;
sh=subplot('position',[.87,.60,.01,yscaleheight]);

hold on;

%draw the scale line (3 units long)

ytick = [0,rnd/4,rnd/2,3*rnd/4,rnd];
yticklabs = num2str(ytick',4);
plot([0,0],[0,rnd]);
hold on;
set(gca,'ylim',[0,rnd],'ytick',ytick,'yticklabel',yticklabs,'xticklabel',[]);
text(+9,rnd/2,'Energy Density (m^2/Hz)','rotation',90,'horizontalalignment','center','fontsize',12);


% PUT BLUE BOX AROUND PLOT AND ADD WEB ADDRESS

blue_box(0.01);

[pr_axes,ax_orig]=page_relative_axes;

axes(pr_axes);
text(.1,.035,'(Colors delineate data records)','FontSize',11);
axes(ax_orig);

% PRINT THE PLOT

set(gcf, 'visible','off');

png_name = [wkg,'plot_',pid,'.png'];
print('-dpng','-r75',png_name);


quit