text:
% Set-up script for ISCE processing to TRAIN.
%
%     Copyright (C) 2016  Bekaert David
%     davidbekaert.com
%
%     This program is free software; you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation; either version 2 of the License, or
%     (at your option) any later version.
%
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%     GNU General Public License for more details.
%
%     You should have received a copy of the GNU General Public License along
%     with this program; if not, write to the Free Software Foundation, Inc.,
%     51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
%
% By David Bekaert
% modifications:
% DB    20/12/2012      Include loading and storing of the IFG's too for
%                       phase-based methods

close all
clear all
clc

%% SETUP SCRIPT compatible with InsarApp.py
% inputs from user
% IFG folders are "date1-date2" or "date1_date2".
% How are date1 and date2 defined?
date_format = 'yyyymmdd';           % this needs to be recognized by matlab.
                                    % mmm for jan format, MM is min two digit,
                                    % while hh is hour 2 digit. You can
                                    % also include - or / if needed. e.g.
                                    % 'yyyy-mmm-dd'

% MASK file. In case there is a file called "MASK" with the same
% size of the data then it will be applied.
% in addtion, nodata in the lonlat data can be masked as well.
nodata_lonlat =0;                   % either a number of [] to ignore

% also save the phase data in a stamps-like format which allows to use plotting
% functionality from stamps.
save_data = 'y';                    % string either 'y' or 'n'

% radar or geo-coordinates
process_type = 'geo';               % a string either 'radar' or 'geo'

% satellite UTC time
UTC_time = '14:30';                 % a string defined as 'HH:MM'

% full path to your DEM file. Easiest is if this is a grd.
DEMfile = '/Volumes/Bekaert_disk2/Projects2/DEM/Taiwan/dem.grd';

% TRAIN directory (best not to change such it remains consistent between programs)
TRAIN_dir = 'aps_train';            % the location where the TRAIN dir will be set-up

% plot debug figures
plot_debug = 'y';                   % string either 'y' or 'n'

% filename of the unwrapped data
unw_filename = 'filt_topophase.unw';    % note the '.geo' string is added automatically if needed

%% CHECK if the DEM file exist
if exist(DEMfile,'file')~=2
    error([DEMfile ' does not exist']);
end

%% check if the type is correct
geo_str = '';
if strcmpi(process_type,'geo')
    geo_str = '.geo';
end
if ~strcmpi(process_type,'geo') && ~strcmpi(process_type,'radar')
    error(['process_type can only be ''geo'' or ''radar'''])
end

%% check on the filenames if needed
[temp1,temp2,temp3] =fileparts(unw_filename);
if strcmpi(geo_str,'.geo') && ~strcmpi(temp3,'.geo')
    fprintf(['Looks like you want to run train in .geo mode\nYour unw_filename does not have a ".geo" extension\n'])
    fprintf(['***You can ignore this if you do not plan to run the TRAIN phase-based correction methods.***\n'])
    % Check if action is needed
    repeat=1;
    while repeat==1
        action_flag= input('Do you want to update it to .geo? [y/n] ','s');
        if strcmpi(action_flag,'y')
            repeat=0;
        elseif strcmpi(action_flag,'n')
            repeat=0;
        end
    end
    if strcmpi(action_flag,'y')
        unw_filename = [unw_filename '.geo'];
    else
        repeat=1;
        while repeat==1
            action_flag= input('Do you want to continue? [y/n] ','s');
            if strcmpi(action_flag,'y')
                repeat=0;
            elseif strcmpi(action_flag,'n')
                repeat=0;
            end
        end
        if strcmpi(action_flag,'n')
           error('Abord by user. Update the unw_filename')
        end
    end
end

%% generate output directory to save along the way data
TRAIN_dir = [pwd filesep TRAIN_dir filesep 'SMALL_BASELINES' filesep];
if exist(TRAIN_dir,'dir')~=7
    mkdir(TRAIN_dir)
else
    error(['Directory exist: remove and rerun: ' TRAIN_dir])
end

%% GENERATE a list of IFG COMBINATIONS
list1 = dir('1*');      % get all the dir pre 2000
list2 = dir('2*');      % get all the dir post 2000
list = [list1 ;list2];  % combine the lists
clear list1 list2

% location of the lon and lat file
latfile = [];
latxmlfile = [];
lonfile = [];
lonxmlfile = [];
losfile = [];
losxmlfile = [];
hgtfile = [];
hgtxmlfile = [];
unwfile = [];
% information on the sensor
insarappxmlfile = [];
% extracting the ifg combinations
for k = 1:length(list)
    ifgday(k,:) = list(k).name;

    % search for a lon, lat and xml file
    if exist([ ifgday(k,:) filesep 'lon.rdr' geo_str],'file')==2 && exist([ ifgday(k,:) filesep 'lat.rdr' geo_str],'file')==2 && exist([ ifgday(k,:) filesep 'lon.rdr' geo_str '.xml'],'file')==2 && exist([ ifgday(k,:) filesep 'lat.rdr' geo_str '.xml'],'file')==2 && isempty(lonfile)
        lonfile = [ ifgday(k,:) filesep 'lon.rdr' geo_str];
        lonxmlfile = [ ifgday(k,:) filesep 'lon.rdr' geo_str '.xml'];
        latfile = [ ifgday(k,:) filesep 'lat.rdr' geo_str];
        latxmlfile = [ ifgday(k,:) filesep 'lat.rdr' geo_str '.xml'];
    end
    % search for the los file
    if exist([ ifgday(k,:) filesep 'los.rdr' geo_str],'file')==2 && exist([ ifgday(k,:) filesep 'los.rdr' geo_str '.xml'],'file')==2 && isempty(losfile)
        losfile = [ ifgday(k,:) filesep 'los.rdr' geo_str];
        losxmlfile = [ ifgday(k,:) filesep 'los.rdr' geo_str '.xml'];
    end
    % search for the hgt file
    if exist([ ifgday(k,:) filesep 'z.rdr' geo_str],'file')==2 && exist([ ifgday(k,:) filesep 'z.rdr' geo_str '.xml'],'file')==2 && isempty(hgtfile)
        hgtfile = [ ifgday(k,:) filesep 'z.rdr' geo_str];
        hgtxmlfile = [ ifgday(k,:) filesep 'z.rdr' geo_str '.xml'];
    end
    % search for the InSAR app xml file
    if exist([ ifgday(k,:) filesep 'insarApp.xml'],'file')==2 && isempty(insarappxmlfile)
        insarappxmlfile = [ ifgday(k,:) filesep 'insarApp.xml'];
    end
    % search for the unwrapped interferogram files
    if exist([ ifgday(k,:) filesep unw_filename],'file')==2
        unwfile{k} = [ifgday(k,:) filesep unw_filename];
    else
        unwfile{k} = [];
    end
end
clear list
% removing the separation between folders
if length(ifgday)<2*length(date_format)
    error('You folder names are not of format "DATE1[-/_/ /]DATE2')
end
ifgday1 = ifgday(:,1:length(date_format));
ifgday2 = ifgday(:,end-length(date_format)+1:end);
clear ifgday

% converting the dates in numerical dates
ifgday = [str2num(datestr(datenum(ifgday1,date_format),'yyyymmdd')) str2num(datestr(datenum(ifgday2,date_format),'yyyymmdd')) ];
clear ifgday1 ifgday2
ifgday = [datenum(num2str(ifgday(:,1)),'yyyymmdd') datenum(num2str(ifgday(:,2)),'yyyymmdd')];
day = unique(ifgday);
n_ifg = size(ifgday,1);
ifgday_ix = [];
for k=1:n_ifg
   ifgday_ix(k,1) = find(ifgday(k,1)==day);
   ifgday_ix(k,2) = find(ifgday(k,2)==day);
end
master_ix = 1;
master_day = day(master_ix);

%% GET LONLAT coordinates
if isempty(lonfile) || isempty(latfile) || isempty(latxmlfile)
    error('Could not file lon.rdr, lat.rdr, or lat.rdr.xml file')
end

% getting the file WIDTH's
[latwidth] = get_parm_xml(latxmlfile,'width');
[lonwidth] = get_parm_xml(lonxmlfile,'width');
[loswidth] = get_parm_xml(losxmlfile,'width');
if latwidth~=lonwidth || latwidth~=loswidth
    error(['Width of the lon, lat and los are not consistent'])
else
    WIDTH = latwidth;
    clear latwidth lonwidth lonwidth;
end
fprintf(['Width is ' num2str(WIDTH) '\n'])

% loading the lonlat information
% loading lon
LON = load_isce(lonfile);
if strcmpi(plot_debug,'y')
    figure('name','LON');
    imagesc(LON)
    axis equal; axis tight; colorbar; axis xy
end

% loading lat
LAT = load_isce(latfile);
if strcmpi(plot_debug,'y')
    figure('name','LAT');
    imagesc(LAT)
    axis equal; axis tight; colorbar; axis xy
end

% check if the lonlat needs to be masked
MASKED_flag = 'n';
if ~isempty(nodata_lonlat)
    MASK = ~((LON==nodata_lonlat) | (LAT==nodata_lonlat));        % mask of points to keep;
    MASKED_flag = 'y';
else
    MASK = logical(ones(size(LON)));
end

% check if a MASK file exist
if exist('MASK','file')==2
    MASK_temp = load_isce(MASK);
    fprintf('Need to write piece of code which combines these MASK files together \n ')
    MASKED_flag = 'y';

    keyboard
end

% combining
lonlat = [LON(MASK) LAT(MASK)];
%lonlat = [reshape(LON,[],1) reshape(LAT,[],1)];
n_ps = size(lonlat,1);
clear LON LAT

% getting a position file that relates a vector to matrix position
[LENGTH] = get_parm_xml(latxmlfile,'length');
[J,I] = meshgrid([1:LENGTH],[1:WIDTH]);
IND = sub2ind([WIDTH LENGTH], I(MASK), J(MASK));
%IND = sub2ind([WIDTH LENGTH], reshape(I,[],1), reshape(J,[],1));
clear I J
% generate template xml for aps file
templatexmlfile = lonxmlfile;
aps_save([TRAIN_dir 'isce2train.mat'],IND,geo_str,WIDTH,LENGTH,MASKED_flag,templatexmlfile)

%% getting the heights if found
phase_based_good = 0;
if ~isempty(hgtfile)
    fprintf('Found heights and will store...\n')
    HGT = load_isce(hgtfile);
    hgt = HGT(MASK);
    clear HGT;
    aps_save([TRAIN_dir 'hgt2.mat'],hgt)
    phase_based_good = 1;
end

%% getting the unwrapped interferograms if found
if ~isempty(unwfile)
    fprintf('Found unwrapped IFGs, will store...\n')
    % try and allocate the memory in advance
    ph_uw = zeros([n_ps n_ifg]);
    for k=1:length(unwfile)
        if exist(unwfile{k},'file')==2
            try
                [DATA] = load_isce(unwfile{k});
                DATA = DATA(:,:,2);                 % first band is amplitude, second band is unwrapped phase
                % apply the mask and store the data
                ph_uw(:,k) = DATA(MASK);
                clear DATA
            catch
                keyboard
            end
        end
    end
    phase_based_good = phase_based_good+1;
    aps_save([TRAIN_dir 'phuw_sb2.mat'],ph_uw)
end

if phase_based_good==2
   fprintf('Found all the information to run phase-based TRAIN corrections\n')
   phase_based_good = 'y';
else
    phase_based_good = 'n';
end
aps_save([TRAIN_dir 'isce2train.mat'],phase_based_good)

%% getting the incidence angle
[TEMP] = load_isce(losfile);
INC = TEMP(:,:,1);
HEADING = TEMP(:,:,2);
clear TEMP
if strcmpi(plot_debug,'y')
    figure('name','INC');
    imagesc(INC)
    axis equal; axis tight; colorbar; axis xy
    figure('name','HEADING');
    imagesc(HEADING)
    axis equal; axis tight; colorbar; axis xy
end
la = INC(MASK)./180*pi;
%la = reshape(INC,[],1)./180*pi;
clear INC
heading = nanmean(HEADING(MASK));
clear HEADING

%% getting the sensor information
% getting the file WIDTH
[SENSOR] = get_parm_xml(insarappxmlfile,'Sensor name');
fprintf(['Sensor is ' num2str(SENSOR) '\n'])

% getting the wavelength hard-coded
if strcmpi(SENSOR,'alos')
    lambda = 23.6;
elseif strcmpi(SENSOR,'radarsat2') || strcmpi(SENSOR,'radarsat') || strcmpi(SENSOR,'rsat') || strcmpi(SENSOR,'rsat2') ||   strcmpi(SENSOR,'envisat') || strcmpi(SENSOR,'ers') || strcmpi(SENSOR,'env') || strcmpi(SENSOR,'ers1')  || strcmpi(SENSOR,'ers2')
    lambda = 5.6;
else
   error('Sensor wavelength not yet hard-coded')
end
% convert to m
lambda =lambda./100;

%% saving the data
% stamps-like format
psver =2;
save([TRAIN_dir 'psver.mat'],'psver')
aps_save([TRAIN_dir 'ps2.mat'],lonlat,n_ifg,n_ps,day,ifgday_ix,ifgday,master_day,master_ix)
aps_save([TRAIN_dir 'la2.mat'],la)
dlmwrite([TRAIN_dir 'lambda.1.in'],lambda)
dlmwrite([TRAIN_dir 'heading.1.in'],heading)
dlmwrite([TRAIN_dir 'WIDTH.txt'],WIDTH)
fid = fopen([TRAIN_dir 'processor.txt'],'w');
fprintf(fid,'isce');
fclose(fid)
% will make few dummy directories
ERA_dir = [TRAIN_dir 'ERA'];
MERRA_dir = [TRAIN_dir 'MERRA'];
MODIS_dir = [TRAIN_dir 'MODIS'];
MERIS_dir = [TRAIN_dir 'MERIS'];
if exist([ERA_dir],'dir')~=7
    mkdir(ERA_dir)
end
if exist([MERRA_dir],'dir')~=7
    mkdir(MERRA_dir)
end
if exist([MERIS_dir],'dir')~=7
    mkdir(MERIS_dir)
end
if exist([MODIS_dir],'dir')~=7
    mkdir(MODIS_dir)
end
cd(TRAIN_dir)
% setting the stamps parameters
getparm                                 % this is for PS
setparm('small_baseline_flag','y')      % change to SB
setparm('percent_rand',[],-1)           % resetting all PS to SB variables
setparm('weed_standard_dev',[],-1)
setparm('scla_drop_index',[],-1)
setparm('unwrap_method',[],-1)
setparm('merge_resample_size',[],-1)
getparm                                 % now this will be for SB
% setting the TRAIN parameters
getparm_aps
setparm_aps('UTC',UTC_time)
setparm_aps('era_datap',ERA_dir)
setparm_aps('merra_datap',MERRA_dir)
setparm_aps('meris_datap',MERIS_dir)
setparm_aps('modis_datap',MODIS_dir)
setparm_aps('demfile',DEMfile)

author: dbekaert

repo: TRAIN

save_path: github-repos/MATLAB/dbekaert-TRAIN


text:
function outLink = hyperlink(url, urlText, altText1, altText2)
% Converts a url to a clickable link in order to improve usability
% when using the MATLAB desktop environment
%
% USAGE:
%   outLink = hyperlink(url, urlText, altText1, altText2)
%
% INPUTS:
%    url:         url address
%    urlText:     url for java
%    altText1:    alternative text before link
%    altText2:    alternative text after link
%
% OUTPUT:
%    outLink:     clickable link

    if nargin < 2
        urlText = url;
    end
    if nargin < 3
        altText1 = '';
    end
    if nargin < 4
        altText2 = '';
    end

    outLink = [altText1, url, altText2];
    if usejava('desktop')
        outLink = ['', urlText, ''];
    end
end

author: opencobra

repo: cobratoolbox


text:
function neqn = p14_neqn ( )

%*****************************************************************************80
%
%% P14_NEQN returns the number of equations for problem P14.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    17 February 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Output, integer NEQN, the number of equations.
%
  neqn = 51;

  return
end

author: johannesgerer

repo: jburkardt-m


function scores=test_svms(models, imnames,featdir)
W=cat(2, models.w);
b=cat(2, models.b);
for i=1:numel(imnames)
    feats=load_feats(featdir, imnames{i});
    feats=xform_feat(feats, models(1).meannrm);
    scores{i}=bsxfun(@plus, feats*W,b);
    fprintf('Done %d/%d\n',i, numel(imnames));
end

function feats=xform_feat(feats, nrm)
feats=feats*20./nrm;

function feats=load_feats(featdir, image_id)
x1=load(fullfile(featdir,'CPMC_segms_150_sp_approx_LBP_f_pca_2500_noncent',[image_id '.mat']));
x2=load(fullfile(featdir,'CPMC_segms_150_sp_approx_SIFT_GRAY_f_g_pca_5000_noncent', [image_id '.mat']));
x3=load(fullfile(featdir,'CPMC_segms_150_sp_approx_SIFT_GRAY_mask_pca_5000_noncent', [image_id '.mat']));
 feats= [x1.D;x2.D; x3.D]';


text:
function F=frame(ftype,varargin);
%FRAME  Construct a new frame
%   Usage: F=frame(ftype,...);
%
%   `F=frame(ftype,...)` constructs a new frame object *F* of type
%   *ftype*. Arguments following *ftype* are specific to the type of frame
%   chosen.
%
%   Time-frequency frames
%   ---------------------
%
%   `frame('dgt',g,a,M)` constructs a Gabor frame with window *g*,
%   time-shift *a* and *M* channels. See the help on |dgt| for more
%   information.
%
%   `frame('dgtreal',g,a,M)` constructs a Gabor frame for real-valued
%   signals with window *g*, time-shift *a* and *M* channels. See the help
%   on |dgtreal| for more information.
%
%   `frame('dwilt',g,M)` constructs a Wilson basis with window *g* and *M*
%   channels. See the help on |dwilt| for more information.
%
%   `frame('wmdct',g,M)` constructs a windowed MDCT basis with window *g*
%   and *M* channels. See the help on |wmdct| for more information.
%
%   `frame('filterbank',g,a,M)` constructs a filterbank with filters *g*,
%   time-shifts of *a* and *M* channels. For the ease of implementation, it
%   is necessary to specify *M*, even though it strictly speaking could be
%   deduced from the size of the windows. See the help on |filterbank| for
%   more information on the parameters. Similarly, you can construct a
%   uniform filterbank by selecting `'ufilterbank'`, a positive-frequency
%   filterbank by selecting `'filterbankreal'` or a uniform
%   positive-frequency filterbank by selecting `'ufilterbankreal'`.
%
%   `frame('nsdgt',g,a,M)` constructs a non-stationary Gabor frame with
%   filters *g*, time-shifts of *a* and *M* channels. See the help on
%   |nsdgt| for more information on the parameters. Similarly, you can
%   construct a uniform NSDGT by selecting `'unsdgt'`, an NSDGT for
%   real-valued signals only by selecting `'nsdgtreal'` or a
%   uniform NSDGT for real-valued signals by selecting `'unsdgtreal'`.
%
%   Wavelet frames
%   --------------
%
%   `frame('fwt', w, J)` constructs a wavelet frame with wavelet definition
%   *w* and *J* number of filterbank iterations. Similarly, a redundant time
%   invariant wavelet representation can be constructed by selecting `'ufwt'`.
%   See the help on |fwt| and |ufwt| for more information.
%
%   `frame('wfbt', wt)` constructs a wavelet filterbank tree defined by
%   the wavelet filterbank tree definition *wt*. Similarly, an undecimated
%   wavelet filterbank tree can be constructed by selecting `'uwfbt'`. See the
%   help on |wfbt| and |uwfbt| for more information.
%
%   `frame('wpfbt', wt)` constructs a wavelet packet filterbank tree
%   defined by the wavelet filterbank tree definition *wt*. Similarly, an
%   undecimated wavelet packet filterbank tree can be constructed by selecting
%   `'uwpfbt'`. See the help on |wpfbt| and |uwpfbt| for more information.
%
%   Pure frequency frames
%   ---------------------
%
%   `frame('dft')` constructs a basis where the analysis operator is the
%   |dft|, and the synthesis operator is its inverse, |idft|. Completely
%   similar to this, you can enter the name of any of the cosine or sine
%   transforms |dcti|, |dctii|, |dctiii|, |dctiv|, |dsti|, |dstii|,
%   |dstiii| or |dstiv|.
%
%   `frame('reddft',red)` constructs so called harmonic Parseval tight
%   frame or redundant dft with redundancy `red`. The frame accepts any
%   `red`, but |frana| will only work for signal lengths `Ls` for which
%   the number of coefficients `Ls*red` is an integer.
%
%   `frame('dftreal')` constructs a normalized |fftreal| basis for
%   real-valued signals of even length only. The basis is normalized
%   to ensure that is it orthonormal.
%
%   Special / general frames
%   ------------------------
%
%   `frame('gen',g)` constructs a general frame with a synthesis matrix *g*.
%   The frame atoms must be stored as column vectors in the matrix.
%
%   `frame('identity')` constructs the canonical orthonormal basis, meaning
%   that all operators return their input as output, so it is the dummy
%   operation.
%
%   Container frames
%   ----------------
%
%   `frame('fusion',w,F1,F2,...)` constructs a fusion frame, which is
%   the collection of the frames specified by *F1*, *F2*,... The vector
%   *w* contains a weight for each frame. If *w* is a scalar, this weight
%   will be applied to all the sub-frames.
%
%   `frame('tensor',F1,F2,...)` constructs a tensor product frame, where the
%   frames *F1, *F2*,... are applied along the 1st, 2nd etc. dimensions. If
%   you don't want any action along a specific dimension, use the `identity`
%   frame along that dimension. Any remaining dimensions in the input
%   signal are left alone.
%
%   Wrapper frames
%   --------------
%
%   Frames types in this section are "virtual". They serve as a wrapper for
%   a different type of frame.
%
%   `frame('erbletfb',fs,Ls,...)` constructs an Erb-let filterbank frame for
%   a given samp. frequency *fs* working with signals of length *Ls*. See
%   |erbfilters| for a description of additional parameters as all
%   parameters other than the frame type string 'erbletfb' are passed to it.
%   NOTE: The resulting frame is defined only for a single signal length
%   *Ls*. Shorter signals will be zero-padded, signals longer than *Ls*
%   cannot be processed.
%   The actual frame type is 'filterbank' or 'filterbankreal'.
%
%   `frame('cqtfb',fs,fmin,fmax,bins,Ls,...)` constructs a CQT filterbank
%   frame for a given samp. frequency *fs* working with signals of length
%   *Ls*. See |cqtfilters| for a description of other parameters.
%   NOTE: The resulting frame is defined only for a single signal length
%   *Ls*. Shorter signals will be zero-padded, signals longer than *Ls*
%   cannot be processed.
%   The actual frame type is 'filterbank' or 'filterbankreal'.
%
%   Examples
%   --------
%
%   The following example creates a Modified Discrete Cosine Transform frame,
%   analyses an input signal and plots the frame coefficients:::
%
%      F=frame('wmdct','gauss',40);
%      c=frana(F,greasy);
%      plotframe(F,c,'dynrange',60);
%
%   See also: frana, frsyn, plotframe

complainif_notenoughargs(nargin,1,'FRAME');

if ~ischar(ftype)
  error(['%s: First argument must be a string denoting the type of ' ...
         'frame.'],upper(mfilename));
end;

ftype=lower(ftype);

% True if the frame only works with real-valued input.
F.realinput=0;

% True if the frame only works with a fixed length.
F.fixedlength = 0;

% Handle the windowed transforms
switch(ftype)
 case {'dgt','dwilt','wmdct','filterbank','ufilterbank',...
       'nsdgt','unsdgt','wfbt','uwfbt','wpfbt'}
  F.g=varargin{1};

 case {'dgtreal','filterbankreal','ufilterbankreal',...
      'nsdgtreal','unsdgtreal'}
  F.g=varargin{1};
  F.realinput=1;

 case {'fwt','ufwt'}
  F.g=varargin{1};
  F.J=varargin{2};
  complainif_notposint(F.J,'J','FRAME');
end;

% Input param checking
switch(ftype)
  case 'fusion'
     wtmp = varargin{1};
     % Check w
     if ~isnumeric(varargin{1}) || ...
        ~(isscalar(wtmp) || numel(wtmp) == numel(varargin) -1)
       error('%s: Weights are not in a correct format.',upper(mfilename));
     end

     % Check frame objects
     for ii=2:numel(varargin)
        complainif_notvalidframeobj(varargin{ii},'FRAME');
     end
  case 'tensor'
     % Check frame objects
     for ii=1:numel(varargin)
        complainif_notvalidframeobj(varargin{ii},'FRAME');
     end
end

% For parsing optional parameters to the transforms.
vargs={};
definput=struct();

%% ---- Pre-optional parameters
% Common operations to deal with the input parameters.
switch(ftype)
  case {'dgt','dgtreal'}
    F.a=varargin{2};
    F.M=varargin{3};

    vargs=varargin(4:end);
    definput.keyvals.lt=[0 1];
    definput.flags.phase={'freqinv','timeinv'};

  case {'dwilt','wmdct'}
    F.M=varargin{2};
  case {'filterbank','ufilterbank','filterbankreal','ufilterbankreal'}
    F.a=varargin{2};
    F.M=varargin{3};

    [F.a,~]=comp_filterbank_a(F.a,F.M,struct());

  case {'nsdgt','unsdgt','nsdgtreal','unsdgtreal'}
    F.a=varargin{2};
    F.M=varargin{3};

    % Sanitize 'a' and 'M'. Make M a column vector of length N,
    % where N is determined from the length of 'a'
    F.a=F.a(:);
    F.N=numel(F.a);
    F.M=bsxfun(@times,F.M(:),ones(F.N,1));
  case {'ufwt'}
    vargs=varargin(3:end);
    definput.flags.scaling={'sqrt','noscale','scale'};
  case {'uwfbt'}
    vargs=varargin(2:end);
    definput.flags.scaling={'sqrt','noscale','scale'};
  case {'wpfbt'}
    vargs=varargin(2:end);
    definput.flags.interscaling={'intsqrt','intnoscale','intscale'};
  case {'uwpfbt'}
    vargs=varargin(2:end);
    definput.flags.interscaling={'intsqrt','intnoscale','intscale'};
    definput.flags.scaling={'sqrt','noscale','scale'};
end;

[F.flags,F.kv]=ltfatarghelper({},definput,vargs);

F.type=ftype;
F.origargs=varargin;
F.vargs=vargs;

%% ------ Post optional parameters

% Default value, works for all bases
F.red=1;

% Default value, frame works for all lengths
F.length=@(Ls) Ls;

switch(ftype)
  case 'gen'
    F.g=varargin{1};
    F.frana=@(insig) F.g'*insig;
    F.frsyn=@(insig) F.g*insig;
    F.length = @(Ls) size(F.g,1);
    F.red = size(F.g,2)/size(F.g,1);

  case 'identity'
    F.frana=@(insig) insig;
    F.frsyn=@(insig) insig;

  case 'dft'
    F.frana=@(insig) dft(insig,[],1);
    F.frsyn=@(insig) idft(insig,[],1);

  case 'reddft'
    F.red = 1;
    if nargin > 1, F.red = varargin{1}; end
    if ~isscalar(F.red) || F.red < 1
        error('%s: Redundancy must be greater or equal to 1.',upper(mfilename));
    end
    F.frana     = @(insig) dft(insig,F.red*size(insig,1),1);
    F.frsyn     = @(insig) postpad(idft(insig,[],1),size(insig,1)/F.red);
    F.lengthcoef= @(Ncoef) Ncoef/F.red;
    F.clength   = @(L) L*F.red;

  case 'dftreal'
    F.frana=@(insig) fftreal(insig,[],1)/sqrt(size(insig,1));
    F.frsyn=@(insig) ifftreal(insig,(size(insig,1)-1)*2,1)*sqrt((size(insig,1)-1)*2);
    F.length=@(Ls) ceil(Ls/2)*2;
    F.lengthcoef=@(Ncoef) (Ncoef-1)*2;
    F.realinput=1;
    F.clength = @(L) floor(L/2)+1;

  case 'dcti'
    F.frana=@(insig) dcti(insig,[],1);
    F.frsyn=@(insig) dcti(insig,[],1);

  case 'dctii'
    F.frana=@(insig) dctii(insig,[],1);
    F.frsyn=@(insig) dctiii(insig,[],1);

  case 'dctiii'
    F.frana=@(insig) dctiii(insig,[],1);
    F.frsyn=@(insig) dctii(insig,[],1);

  case 'dctiv'
    F.frana=@(insig) dctiv(insig,[],1);
    F.frsyn=@(insig) dctiv(insig,[],1);

  case 'dsti'
    F.frana=@(insig) dsti(insig,[],1);
    F.frsyn=@(insig) dsti(insig,[],1);

  case 'dstii'
    F.frana=@(insig) dstii(insig,[],1);
    F.frsyn=@(insig) dstiii(insig,[],1);

  case 'dstiii'
    F.frana=@(insig) dstiii(insig,[],1);
    F.frsyn=@(insig) dstii(insig,[],1);

  case 'dstiv'
    F.frana=@(insig) dstiv(insig,[],1);
    F.frsyn=@(insig) dstiv(insig,[],1);

  case 'dgt'
    F.coef2native=@(coef,s) reshape(coef,[F.M,s(1)/F.M,s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(comp_dgt(insig,F.g,F.a,F.M,F.kv.lt,F.flags.do_timeinv,0,0));
    F.frsyn=@(insig) comp_idgt(F.coef2native(insig,size(insig)),F.g,F.a,F.kv.lt,F.flags.do_timeinv,0);
    F.length=@(Ls) dgtlength(Ls,F.a,F.M,F.kv.lt);
    F.red=F.M/F.a;

  case 'dgtreal'
    F.coef2native=@(coef,s) reshape(coef,[floor(F.M/2)+1,s(1)/(floor(F.M/ ...
                                                      2)+1),s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(comp_dgtreal(insig,F.g,F.a,F.M,F.kv.lt,F.flags.do_timeinv));
    F.frsyn=@(insig) comp_idgtreal(F.coef2native(insig,size(insig)),F.g,F.a,F.M,F.kv.lt,F.flags.do_timeinv);
    F.length=@(Ls) dgtlength(Ls,F.a,F.M,F.kv.lt);
    F.red=F.M/F.a;
    F.lengthcoef=@(Ncoef) Ncoef/(floor(F.M/2)+1)*F.a;
    F.clength = @(L) L/F.a*(floor(F.M/2)+1);

  case 'dwilt'
    F.coef2native=@(coef,s) reshape(coef,[2*F.M,s(1)/F.M/2,s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(comp_dwilt(insig,F.g,F.M));
    F.frsyn=@(insig) comp_idwilt(F.coef2native(insig,size(insig)),F.g);
    F.length=@(Ls) dwiltlength(Ls,F.M);

  case 'wmdct'
    F.coef2native=@(coef,s) reshape(coef,[F.M,s(1)/F.M,s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(comp_dwiltiii(insig,F.g,F.M));
    F.frsyn=@(insig) comp_idwiltiii(F.coef2native(insig,size(insig)),F.g);
    F.length=@(Ls) dwiltlength(Ls,F.M);

  case 'filterbank'
    F.red=sum(F.a(:,2)./F.a(:,1));
    F.length=@(Ls) filterbanklength(Ls,F.a);
    F.lengthcoef=@(Ncoef) Ncoef/F.red;
    F.native2coef=@(coef) cell2mat(coef(:));
    F.coef2native=@(coef,s) vect2cell(coef,round(s(1)/F.red*F.a(:,2)./F.a(:,1)));
    F.frana=@(insig) F.native2coef(comp_filterbank(insig,F.g,F.a));
    F.frsyn=@(insig) comp_ifilterbank(F.coef2native(insig,size(insig)),...
                                      F.g,F.a,round(size(insig,1)/F.red));
    F.destructor=@() clear('comp_filterbank','comp_ifilterbank');

  case 'filterbankreal'
    F.red=2*sum(F.a(:,2)./F.a(:,1));
    F.length=@(Ls) filterbanklength(Ls,F.a);
    F.lengthcoef=@(Ncoef) 2*Ncoef/(F.red);
    F.native2coef=@(coef) cell2mat(coef(:));
    F.coef2native=@(coef,s) vect2cell(coef,round(2*s(1)/F.red*F.a(:,2)./F.a(:,1)));
    F.frana=@(insig) F.native2coef(comp_filterbank(insig,F.g,F.a));
    F.frsyn=@(insig) 2*real(comp_ifilterbank(F.coef2native(insig,size(insig)),F.g,F.a,...
                                             round(2*size(insig,1)/F.red)));
    F.destructor=@() clear('comp_filterbank','comp_ifilterbank');

  case 'ufilterbank'
    F.red=sum(F.a(:,2)./F.a(:,1));
    F.length=@(Ls) filterbanklength(Ls,F.a);
    F.lengthcoef=@(Ncoef) round(Ncoef/F.red);
    F.coef2native=@(coef,s) reshape(coef,[s(1)/F.M,F.M,s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(ufilterbank(insig,F.g,F.a));
    F.frsyn=@(insig) ifilterbank(F.coef2native(insig,size(insig)),F.g,F.a);

  case 'ufilterbankreal'
    F.red=2*sum(F.a(:,2)./F.a(:,1));
    F.length=@(Ls) filterbanklength(Ls,F.a);
    F.lengthcoef=@(Ncoef) round(Ncoef/F.red*2);
    F.coef2native=@(coef,s) reshape(coef,[s(1)/F.M,F.M,s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(ufilterbank(insig,F.g,F.a));
    F.frsyn=@(insig) 2*real(ifilterbank(F.coef2native(insig,size(insig)),F.g, ...
                                        F.a));

  case 'nsdgt'
    F.coef2native=@(coef,s) mat2cell(coef,F.M,s(2));
    F.native2coef=@(coef) cell2mat(coef(:));
    F.length=@(Ncoef) sum(F.a);
    F.lengthcoef=@(Ncoef) sum(F.a);
    F.red=sum(F.M)/sum(F.a);
    F.frana=@(insig) F.native2coef(nsdgt(insig,F.g,F.a,F.M));
    F.frsyn=@(insig) insdgt(F.coef2native(insig,size(insig)),F.g,F.a);

  case 'unsdgt'
    F.coef2native=@(coef,s) reshape(coef,[F.M(1),s(1)/F.M(1),s(2)]);
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(unsdgt(insig,F.g,F.a,F.M));
    F.frsyn=@(insig) insdgt(F.coef2native(insig,size(insig)),F.g,F.a);
    F.length=@(Ncoef) sum(F.a);
    F.lengthcoef=@(Ncoef) sum(F.a);
    F.red=sum(F.M)/sum(F.a);

  case 'nsdgtreal'
    F.coef2native=@(coef,s) mat2cell(coef,floor(F.M/2)+1,s(2));
    F.native2coef=@(coef) cell2mat(coef(:));
    F.frana=@(insig) F.native2coef(nsdgtreal(insig,F.g,F.a,F.M));
    F.frsyn=@(insig) insdgtreal(F.coef2native(insig,size(insig)),F.g,F.a,F.M);
    F.length=@(Ncoef) sum(F.a);
    F.lengthcoef=@(Ncoef) sum(F.a);
    F.red=sum(F.M)/sum(F.a);
    F.clength=@(L) sum(floor(F.M/2)+1);

  case 'unsdgtreal'
    F.coef2native=@(coef,s) reshape(coef,floor(F.M(1)/2)+1,s(1)/ ...
                                    (floor(F.M(1)/2)+1),s(2));
    F.native2coef=@(coef)   reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(unsdgtreal(insig,F.g,F.a,F.M));
    F.frsyn=@(insig) insdgtreal(F.coef2native(insig,size(insig)),F.g,F.a,F.M);
    F.length=@(Ncoef) sum(F.a);
    F.lengthcoef=@(Ncoef) sum(F.a);
    F.red=sum(F.M)/sum(F.a);
    F.clength=@(L) numel(F.M)*(floor(F.M(1)/2)+1);

  case 'fusion'
    F.w=varargin{1};
    F.frames=varargin(2:end);

    if any(cellfun(@(fEl) fEl.realinput,F.frames))
        error(['%s: Real-valued-input-only frames are not currently ',...
               'supported in the fusion frame.'],upper(mfilename));
    end

    F.Nframes=numel(F.frames);
    F.w=bsxfun(@times,F.w(:),ones(F.Nframes,1));
    F.length = @(Ls) comp_framelength_fusion(F,Ls);
    F.red=sum(cellfun(@framered,F.frames));

    % These definitions binds F itself, so they must execute last
    F.frana=@(insig) comp_frana_fusion(F,insig);
    F.frsyn=@(insig) comp_frsyn_fusion(F,insig);

  case 'tensor'
    % This frame type is currently broken. It must be reworked to reshape
    % to the standard layout in order not to break all the assumptions.
    F.frames=varargin;
    F.Nframes=numel(F.frames);
    for ii=1:F.Nframes
        if F.frames{ii}.realinput
            error(['It is not safe to embed a real-valued-input-only frame ' ...
                   'into the tensor frame.']);
        end;
    end;

    F.frana=@(insig) comp_frana_tensor(F,insig);
    F.frsyn=@(insig) comp_frsyn_tensor(F,insig);

    F.length=@(Ls) comp_framelength_tensor(F,Ls);

    F.red=prod(cellfun(@framered,F.frames));

  case {'fwt','dwt'}
    % We have to initialize F.g here already
    [F.g, F.info]=fwtinit({'strict',F.g});
    F.red= 1/(F.g.a(1)^(F.J)) + sum(1./(F.g.a(1).^(0:F.J-1))*sum(1./F.g.a(2:end)));
    F.frana=@(insig) wavcell2pack(comp_fwt(insig,F.g.h,F.g.a,F.J,'per'));
    F.frsyn=@(insig) comp_ifwt(...
                        wavpack2cell(insig,fwtclength(size(insig,1)/F.red,F.g,F.J)),...
                               F.g.g,F.g.a,F.J,size(insig,1)/F.red,'per');
    F.length=@(Ls) fwtlength(Ls,F.g,F.J);
  case {'wfbt'}
    [F.g,F.info]=wfbtinit({'strict',F.g});
    F.red = sum(1./treeSub(F.g));
    % comp_ specific
    [F.wtPath, F.rangeLoc, F.rangeOut] = treeBFranges(F.g);

    F.coef2native = @(coef,s) wavpack2cell(coef,wfbtclength(s(1)/F.red,F.g));
    F.native2coef = @(coef) wavcell2pack(coef);
    F.frana=@(insig) F.native2coef(comp_wfbt(insig,F.g.nodes(F.wtPath),...
                                             F.rangeLoc,F.rangeOut,'per'));
    F.frsyn=@(insig) comp_iwfbt(F.coef2native(insig,size(insig)),...
                                F.g.nodes(F.wtPath(end:-1:1)),...
                                [nodesInLen(F.wtPath(end:-1:1),size(insig,1)/F.red,1,F.g);size(insig,1)/F.red],...
                                F.rangeLoc(end:-1:1),F.rangeOut(end:-1:1),...
                                'per');
    F.length=@(Ls) wfbtlength(Ls,F.g);
  case {'wpfbt'}
    F.g=wfbtinit({'strict',F.g});
    F.red = sum(cellfun(@(aEl) sum(1./aEl),nodesSub(nodeBForder(0,F.g),F.g)));
    % comp_ specific
    F.wtPath = nodeBForder(0,F.g);
    F.rangeLoc = nodesLocOutRange(F.wtPath,F.g);
    [F.pOutIdxs,F.chOutIdxs] = treeWpBFrange(F.g);

    F.coef2native = @(coef,s) wavpack2cell(coef,...
                    s(1)./cell2mat(cellfun(@(aEl) aEl(:),...
                    reshape(nodesSub(nodeBForder(0,F.g),F.g),[],1),...
                    'UniformOutput',0))./F.red);
    F.native2coef = @(coef) wavcell2pack(coef);

    F.frana=@(insig) F.native2coef(...
                        comp_wpfbt(insig,F.g.nodes(F.wtPath),...
                                   F.rangeLoc,'per',F.flags.interscaling));
    F.frsyn=@(insig) comp_iwpfbt(F.coef2native(insig,size(insig)),...
                                 F.g.nodes(F.wtPath(end:-1:1)),...
                                 F.pOutIdxs,F.chOutIdxs,...
                                 size(insig,1)/F.red,...
                                 'per',F.flags.interscaling);
    F.length=@(Ls) wfbtlength(Ls,F.g);
  case {'ufwt'}
    F.g=fwtinit({'strict',F.g});
    F.coef2native = @(coef,s) reshape(coef,[s(1)/(F.J*(numel(F.g.a)-1)+1),F.J*(numel(F.g.a)-1)+1,s(2)]);
    F.native2coef = @(coef) reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(comp_ufwt(insig,F.g.h,F.g.a,F.J,F.flags.scaling));
    F.frsyn=@(insig) comp_iufwt(F.coef2native(insig,size(insig)),F.g.g,F.g.a,F.J,F.flags.scaling);
    F.length=@(Ls) Ls;
    F.red=(F.J*(numel(F.g.a)-1)+1);
  case {'uwfbt'}
    F.g=wfbtinit({'strict',F.g});

    % comp_ specific
    [F.wtPath, F.rangeLoc, F.rangeOut] = treeBFranges(F.g);
    F.nodesUps = nodesFiltUps(F.wtPath,F.g);
    F.red = sum(cellfun(@numel,F.rangeOut));

    F.coef2native = @(coef,s) reshape(coef,[s(1)/F.red,F.red,s(2)]);
    F.native2coef = @(coef) reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(...
                        comp_uwfbt(insig,F.g.nodes(F.wtPath),F.nodesUps,...
                                   F.rangeLoc,F.rangeOut,F.flags.scaling));
    F.frsyn=@(insig) comp_iuwfbt(F.coef2native(insig,size(insig)),...
                                 F.g.nodes(F.wtPath(end:-1:1)),...
                                 F.nodesUps(end:-1:1),F.rangeLoc(end:-1:1),...
                                 F.rangeOut(end:-1:1),F.flags.scaling);
    F.length=@(Ls) Ls;
  case {'uwpfbt'}
    F.g= wfbtinit({'strict',varargin{1}});
    F.red = sum(cellfun(@(fEl) numel(fEl.g),F.g.nodes));
    % comp_ specific
    F.wtPath = nodeBForder(0,F.g);
    F.nodesUps = nodesFiltUps(F.wtPath,F.g);
    F.rangeLoc = nodesLocOutRange(F.wtPath,F.g);
    [F.pOutIdxs,F.chOutIdxs] = treeWpBFrange(F.g);

    F.coef2native = @(coef,s) reshape(coef,[s(1)/F.red,F.red,s(2)]);
    F.native2coef = @(coef) reshape(coef,[size(coef,1)*size(coef,2),size(coef,3)]);
    F.frana=@(insig) F.native2coef(...
                        comp_uwpfbt(insig,F.g.nodes(F.wtPath),F.rangeLoc,...
                                    F.nodesUps,F.flags.scaling,...
                                    F.flags.interscaling));
    F.frsyn=@(insig) comp_iuwpfbt(F.coef2native(insig,size(insig)),...
                                  F.g.nodes(F.wtPath(end:-1:1)),...
                                  F.nodesUps(end:-1:1),F.pOutIdxs,F.chOutIdxs,...
                                  F.flags.scaling,F.flags.interscaling);
    F.length=@(Ls) Ls;

  %%%%%%%%%%%%%%%%%%%%
  %% WRAPPER FRAMES %%
  %%%%%%%%%%%%%%%%%%%%
  case {'erbletfb','cqtfb'}
    switch(ftype)
        case 'erbletfb'
            [g,a,~,L] = erbfilters(varargin{:});
        case 'cqtfb'
            [g,a,~,L] = cqtfilters(varargin{:});
    end
    % Search for the 'complex' flag
    do_complex = ~isempty(varargin(strcmp('complex',varargin)));
    if do_complex
       F = frameaccel(frame('filterbank',g,a,numel(g)),L);
    else
       F = frameaccel(frame('filterbankreal',g,a,numel(g)),L);
    end
    F.fixedlength = 1;

  otherwise
    error('%s: Unknown frame type: %s',upper(mfilename),ftype);

end;

% This one is placed at the end, to allow for F.red to be defined
% first.
if ~isfield(F,'lengthcoef')
    F.lengthcoef=@(Ncoef) Ncoef/framered(F);
end;

author: ltfat

repo: ltfat

save_path: github-repos/MATLAB/ltfat-ltfat


text:
function h = helicorder(varargin)

%HELICORDER: Helicorder generates a multi-line display of input waveform
%   data. Helicorder can display waveform data from multiple stations and
%   channels as long as the time range is the same for each. Helicorder
%   can also display detected events over the raw data.
%
%USAGE: helicorder()----------------------- Empty helicorder
%       helicorder(wave)------------------- Default properties
%       helicorder(wave,prop_1,val_1,...)-- User-defined properties
%
%VALID PROP/VAL:
%  'mpl'-->(Minutes Per Line)
%     single numeric value specifying number of minutes per helicorder
%     line. (This applies to all waveforms in the helicorder)
%     DEFAULT = 10
%
%  'ytick'-->(Y-Axis Ticks)
%     single numeric value specifying number of minutes between y-axis
%     tick marks and labels. 'ytick' depends on 'mpl', for example, if
%     'mpl' is 20, an attempt to set 'ytick' to 10 will be result in
%     'ytick' rounded to nearest multiple of 'mpl' (20 in this case). Note
%     also that manually setting 'mpl' will result in an automatic
%     adjustment of 'ytick'.
%     DEFAULT = 30
%
%  'e_sst'-->(Event Start/Stop Times)
%     If wave contains only one wavefrom object, 'e_sst' can be entered
%     as an Nx2 array of matlab times. If there is more than one waveform
%     in wave (M waveforms), 'e_sst' must be a 1xM cell array, each cell
%     containing Nx2 array of start/stop times. A single waveform can also
%     have multiple sets of start/stop times as in a waveform with multiple
%     event families. In this case 'e_sst' should be entered as a 1xL cell
%     array for a particular waveform (with L separate groups to be
%     highlighted).
%     EXAMPLE: h = helicorder(w,'e_sst',A) where w is a 1x2 waveform object
%        A is a 1x2 cell array
%        A(1,1) is a Nx2 numeric array
%        A(1,2) is a 1x3 cell
%        A(1,2){1} is a Nx2 numeric array
%        A(1,2){2} is a Nx2 numeric array
%        A(1,2){3} is a Nx2 numeric array
%     DEFAULT = [] (No events)
%
%  'trace_color'-->(Trace Color)
%     If wave contains only one wavefrom object, 'trace_color' can be
%     entered as a 1x3 array of RGB values (between 0 and 1). For wave
%     arguments longer than 1, 'trace_color' should be entered as a 1xN
%     cell array, each containing a 1x3 array of RGB values.
%
%  'event_color'-->(Event Color)
%     Only use this property when 'e_sst' property exists. 'event_color'
%     follows the same structure as 'e_sst'. If wave contains only one
%     waveform object, 'event_color' can be entered as a 1x3 array of RGB
%     values (between 0 and 1). For wave arguments longer than 1,
%     'event_color' should be entered as a 1xN cell array, each containing
%     a 1x3 array of RGB values. If multiple event sets exist for a single
%     waveform, 'event_color' must be specified accordingly
%     EXAMPLE: h = helicorder(w,'e_sst',A,'event_color',B)
%        A is the same as the previous example
%        B(1,1) is a 1x3 RGB array
%        A(1,2) is a 1x3 cell
%        A(1,2){1} is a 1x3 RGB array
%        A(1,2){2} is a 1x3 RGB array
%        A(1,2){3} is a 1x3 RGB array
%     DEFAULT = [] (No events)
%
%  'display'-->(Multi-Waveform Display Type)
%     'single' - used when there is only one waveform in wave. None of the
%        multi-waveform display types can be set unless multiple waveforms
%        are passed.
%     'stack' - plots multiple waveforms over top of each other.
%        (Multiple motivations could exist for wanting to display data in
%        this way such as overlaying VLP energy)
%     'alternate' - Alternate helicorder lines
%        i.e. Helicorder containing 3 waveforms would look like:
%        w1(t0 --> t1), w2(t0 --> t1), w3(t0 --> t1), w1(t1 -- t2),...
%     'group' - Plot all lines from individual waveforms in separate blocks
%        i.e. Helicorder containing 3 waveforms would look like:
%        w1(t0 --> t1), w1(t1 --> t2), w1(t2 --> t3),... w1(t_end-1 -->
%        t_end),w2(t0 --> t1), w2(t1 --> t2), w2(t2 --> t3),...
%     DEFAULT = 'single','alternate' (1 or multiple waveforms)
%
%INPUTS: wave     - a waveform object to be plotted on multiple helicorder
%                   trace lines
%        varargin - user-defined helicorder properties (argument pairs)
%
%OUTPUTS: h - helicorder object
%
%   See also HELICORDER/DISP, HELICORDER/BUILD, HELICORDER/GET,
%            HELICORDER/SET
%
% Author: Dane Ketner, Alaska Volcano Observatory
% $Date$
% $Revision$

%% EMPTY HELICORDER
if nargin == 0
   h.wave = waveform();
   h.e_sst = {};
   nw = 0;

%% ANOTHER HELICORDER OBJECT
elseif nargin==1 && isa(varargin{1},'helicorder')
    h = varargin{1};
    return

%% WAVEFORM OBJECT PROVIDED
elseif (nargin >= 1) && isa(varargin{1},'waveform')
   h.wave = varargin{1}; % Helicorder waveform data
   nw = numel(h.wave);
else
   error('helicorder: first argument must be a waveform')
end

%% DEFAULT PROPERTIES
for n = 1:nw,
   h.e_sst{n} = []; % Empty event start/stop time array
end

h.mpl = 10;   % 10 Minutes per line
h.ytick = 30; % Time Ticks every 30 Minutes (3 Lines)

if nw <= 1                     % Default trace color scheme
   h.trace_color{1} = [0 0 0]; % black
   h.event_color{1} = [1 0 0]; % red
elseif nw == 2
   h.trace_color{1} = [0 0 0]; % black
   h.trace_color{2} = [0 0 1]; % blue
else
   rr = zeros(1,nw);              % red values
   gg = round(linspace(-255,255,nw));
   gg = sign(sign(gg)+1).*gg/255; % green values
   bb = round(linspace(255,-255,nw));
   bb = sign(sign(bb)+1).*bb/255; % blue values
   for n = 1:nw
      h.trace_color{n} = [rr(n) gg(n) bb(n)];
   end
end

h.event_color = [];

if nw <= 1
   h.display = 'single';
else
   h.display = 'alternate';
end

h.scale = ones(1,nw); % Waveform amplitude scale is automatically set based
                      % on maximum values in wave. scale can be set after
                      % the fact. 1 is default, going higher or lower will
                      % scale wave proportionally.
h = class(h,'helicorder');

%% USER-DEFINED PROPERTIES
if (nargin > 1)
   v = varargin(2:end);
   nv = nargin-1;
   if ~rem(nv,2) == 0
      error(['helicorder: Arguments after wave must appear in ',...
             'property name/val pairs'])
   end
   for n = 1:2:nv-1
      name = lower(v{n});
      val = v{n+1};
      switch name
         case 'e_sst' % Event Start/Stop Times
            h = set(h,'e_sst',val);
            for n = 1:nw       % Default event color scheme (shades of red)
               if ~iscell(h.e_sst{n})
                  h.event_color{n} = [1 0 0];
               else
                  for m = 1:numel(h.e_sst{n})
                     rr = m/numel(h.e_sst{n});
                     h.event_color{n}{m} = [rr 0 0];
                  end
               end
            end
         case 'ytick' % Y-Tick Spacing
            h = set(h,'ytick',val);
         case 'mpl' % Minutes Per Line
            h = set(h,'mpl',val);
         case 'trace_color' % Minutes Per Line
            h = set(h,'trace_color',val);
         case 'event_color' % Minutes Per Line
            h = set(h,'event_color',val);
         case 'display' % Multiple Waveform Display Type
            h = set(h,'display',val);
         case 'scale' % Waveform Data Scaling Factor
            h = set(h,'scale',val);
         otherwise
            error('helicorder: Property name not recognized')
      end
   end
end

author: geoscience-community-codes

repo: GISMO


text:
function varargout=specsubm(varargin)
%
% For calling details please see v_specsubm.m
%
% This dummy routine is included for backward compatibility only
% and will be removed in a future release of voicebox. Please use
% v_specsubm.m in future and/or update with v_voicebox_update.m
%
%      Copyright (C) Mike Brookes 2018
%      Version: $Id: specsubm.m 10863 2018-09-21 15:39:23Z dmb $
%
if nargout
    varargout=cell(1,nargout);
    [varargout{:}]=v_specsubm(varargin{:});
else
    v_specsubm(varargin{:});
end

author: ImperialCollegeLondon

repo: sap-voicebox