% -----------------------------------------------------------
%  sw06_filespcgrm_all_sw53_shru
%
%     PROGRAM SPECIFIC TO SHRU3_SW53
%
%  This program loops through data files in the shru directory 
%    and plots the spectrogram every 1 record (~64sec). It can 
%    also write audio files if applied. 
%
%  The directory structure the shru creates is as the following:
%    Root directory (the number of shru system)
%       |
%       --- Sub-directory named by date
%              |
%              --- data file named by data and time
%
%  The frequency band of spectrogram plots are 0-4900Hz by default
%    or any range set by users.
%
%  This program reads a certain record into a structure containing
%    all headers and data.
%
%    Using Kaiser window with fft of order 6 (which is close to hanning)
%
%   Lots of bookkeeping done here!
%
%
%             
%
%    subroutines needed for this program
%    1) sw06_readFile_shru
%    2) sw06_getAllDRH_shru
%    3) sw06_getDRH_shru
%    4) sw06_gettime_shru
%    5) specgram_yt_t0
%
%  Modifications:
%
%     Version converts .wav files to mp3
%     Also uses 'new' naming convention that includes timestamp.
%
% --------------------------------------------------------------

% start fresh
clear all; pack;

% -------------------------------------------------------------
%  some variables to set

sensitivity    = 170;                % dB muP/V, it needs double check!!
fft_size       = 1024;               % fft size (1024 ~.11 sec)
window_length  = 1024;               % Kaiser of order 6 window 
overlap_length = 0;                  % overlap (will slow down)
% crange         = [40 110];           % color ranges for freq bands
crange         = [50 120];           % color ranges for freq bands
norm_fac       = 10 ;                % normalization for wav file

iswritewav     = 1;                  % save audio in a wav file
do_spec        = 1;                  % not implemented yet, do spectrogram?

% where shru subroutines live  - Newhall specific
program_dir = '/home/gringo/arthur/sw06/shru/programs'; 

 %  for shru3 - Newhall specific
rootdir = '/home/gringo/arthur/sw06/shru/data/shru3_sw53';        % Data files for SHRU3_SW53 ONLY

% Prefix string for plot and filename
title_string0 = ['SHRU3\_sw53, '];   pixname0 =['SHRU3_sw53'];

% where to store the images
pixdir = [rootdir '/spectrograms'];

% where to store the audio files
if iswritewav,
    wavdir = [rootdir '/audio'];
end
 
% -------------------------------------------------------------

% add the program directory into the matlab path
addpath(program_dir)

% get data subdirectories (by date) from root
%    each day has a subdir, only get data dirs
subdir = dir([rootdir '/D*']);

f_range_setting = ['f_range = [ 0 max(SPL_f)];'];     % default freq band

% check to see if the pixdir and wavdir exist
if ~exist(pixdir),
    mkdir(pixdir);
end

if ~exist(wavdir,'dir'),
    mkdir(wavdir);
end

figure(999);                     % open a figure to plot spectrograms

% go through all subdirs

tic;

% Original loop, use this one!
for idir = 1:length(subdir),   

   % had a crash, start with D0811 (index 13+2)
   %%  for idir = 15:length(subdir),  
 
    
    shru_dir  = [rootdir '/' subdir(idir).name];
    shru_file = dir([shru_dir '/*.DAT']);          %  get all data files
    spec_dir  = [pixdir '/' subdir(idir).name];
      if ~exist(spec_dir),   mkdir(spec_dir);   end
    audio_dir = [wavdir '/' subdir(idir).name];
      if ~exist(audio_dir),   mkdir(audio_dir);   end  
 
    % go through all data files in certain directory
    for ifile = 1:length(shru_file),

 
        fprintf('Working on file: %s',shru_file(ifile).name);
        
        filename = [shru_dir '/' shru_file(ifile).name];
        fid = fopen(filename,'rb','ieee-be');       % big endian for shru
        fprintf('\n\tReading file...  please wait...\n');

        [SHRU,status] = sw06_readFile_shru(fid);    % read all records for ifileth file, SHRU is a structure
        fclose(fid);
        
        title_string = [title_string0 shru_file(ifile).name ', Starting Time: '];    % spectrogrm title
        

        % go through every record
        for irec = 1:length(SHRU),
            
            fprintf('Processing #%d/%d record...\n',irec,length(SHRU))

            [TM] = sw06_gettime_shru(SHRU(irec).hdr);
            start_time = datenum(TM.yr,1,1)-1+TM.yday;             %in days
            start_time_sec = start_time*24*60*60;                  %in seconds

            STM = SHRU(irec).hdr.atime([1 2 4 5]);                       % get hr min

            title_string_final = [title_string datestr(start_time)];      % the title for the spectrogram
            
            Info = strvcat('         sampling freq',sprintf('         [%7.2fHz]',SHRU(irec).hdr.rhfs),...
                '        FFT pts',sprintf('       [%d]',fft_size),...
                '        window pts',sprintf('       [%d]',window_length),...
                '        overlap pts',sprintf('      [%d]',overlap_length),...
                '          conversion factor',sprintf('         [%ddB %sPa/V]',sensitivity,'\mu'));
            
           % Calculate Power Spectrogram -----------------------------------------------------

            if (do_spec == 1),   
      
              [SPL,SPL_f,SPL_t_sec] = specgram_yt_t0( SHRU(irec).data*10^(sensitivity/20), ...
                 fft_size, SHRU(irec).hdr.rhfs, kaiser(window_length,6), overlap_length, start_time_sec );
  
   
              eval(f_range_setting)             % execute the string of code above... set the freq band
              SPL_t = SPL_t_sec/60/60/24;       % the time axis of the power spectrogram [in days]
                
                for ifrq = 1:size(f_range,1);   % dummy loop, but keep it for multiple freq bands                  
   
                    clf;
                    itemp = find(SPL_f>=f_range(ifrq,1)&SPL_f<=f_range(ifrq,2));       % select freq band
                    imagesc(SPL_t,SPL_f(itemp),SPL(itemp,:)); axis xy;        % plot the spectrogram

                        % write title, labels and ticks 

                        % settings of fontsize
                        titlefontsize = 16;
                        axislabelfontsize = 14;
                        ticklabelfontsize = 10;
                        infolabelfontsize = 10;
                        colorbarfontsize = 10;

                        xtick = SPL_t(find( mod(SPL_t,15/60/60/24) < SPL_t(2)-SPL_t(1) ));
                        dytick = fix((f_range(ifrq,2)-f_range(ifrq,1))/12/100)*100;    % auto adjust ytick
                        ytick = [f_range(ifrq,1):400:f_range(ifrq,2)];
                        xlabelstr = ['Time (hr:min:sec)'];
                        ylabelstr = ['Frequency (Hz)'];

                        set(gca,'position',[.13 .11 .7 .815],...
                            'xtick',xtick,'xticklabel',datestr(xtick,'HH:MM:SS'),...
                            'ylim',[f_range(ifrq,1) f_range(ifrq,2)],...
                            'ytick',ytick,...
                            'FontSize',ticklabelfontsize);
                        set(get(gca,'title'),'string',title_string_final,'FontSize',titlefontsize)
                        set(get(gca,'xlabel'),'string',xlabelstr,'FontSize',axislabelfontsize)
                        set(get(gca,'ylabel'),'string',ylabelstr,'FontSize',axislabelfontsize)

                    % plot the colorbar
                    caxis(crange(ifrq,:))       % set the color range
                    hcb = colorbar;           % plot the colorbar
                    
                        % write info and label
                        ytick = [crange(ifrq,1):10:crange(ifrq,2)];

                        set(hcb,'position',[.865 .11 .03 .4],...
                            'ytick',ytick,...
                            'FontSize',colorbarfontsize);
                        set(get(hcb,'title'),'string',strvcat(Info,' ',' '),...
                            'fontsize',infolabelfontsize,...
                            'HorizontalAlignment','center');
                        set(get(hcb,'ylabel'),'string',strvcat('Power Spectral Density', ...
                            '(dB re 1\muPa per 1Hz)'),...
                            'fontsize',colorbarfontsize);
                    
                    % save the spectrogram 
                    pixname = sprintf('%s_%s_%s_r%3.3d.png',pixname0,shru_file(ifile).name(1:end-4),STM,irec);
                           
                        
                        % switch ifrq                         %% keep it for multiple freq bands
                        %     case 1
                        %         pixname = [pixdir '/1stfreq_' pixname];
                        %     case 2
                        %         pixname = [pixdir '/2ndfreq_' pixname];
                        % end
                        
                    pixname_lg = [spec_dir '/' pixname];
                        
                    % write spec files if applied

                    print('-dpng','-r150', pixname_lg)
                   

                  end   % if do_spec
  
                    % write audio files if applied  -------------------------------------------------
                    if iswritewav, 
                        wavname = sprintf('%s_%s_%s_r%3.3d.wav',pixname0,shru_file(ifile).name(1:end-4),STM,irec);
                        wavname = wavname(1:end-4);                          % remove extension
                        wavname_lg = [audio_dir '/' wavname];
                        % normalize
                        % wavwrite(SHRU(irec).data*norm_fac,SHRU(irec).hdr.rhfs,8,wavname_lg);
                    
                        wavwrite(SHRU(irec).data*norm_fac,SHRU(irec).hdr.rhfs,16,wavname_lg); 
   
                  
                        % convert to mp3
                        [s,w] = system(['lame --silent ', wavname_lg, '.wav   ', wavname_lg, '.mp3 ' ]);
                        [s,w] = system(['rm ',wavname_lg, '.wav   ']);

                    end    
                    % ------------------------------------------------------------------------
                
                end   % frequency band loop
 
        end % records loop
        
        tm_chk = toc;
        fprintf('Duration since starting: %f minutes (%f hrs)...\n',tm_chk/60,tm_chk/60/60);

    end % files loop
    
end % directories loop

rmpath(program_dir)
close(999)
