%*******************  DIsjoint STA/LTA  ******************
%	Returns the disjoint STA/LTA r for a time series v,  *
%      and the index of to the first element of r, j1.   *
%   ws and wl are numbers or (row) arrays specifying     *
%      the short and long weigths, respectively,         * 
% Examples:                                              *
%    [r,j1]=dishlo(v,3,8) uses uniform weight arrays 3   *
%                                 and 8 elements long.   *
%    [r,j1]=dishlo(v,[3 1],[8,1]) uses weight arrays 3   *
%       and 8 elements long diminishing as [3 2 1] and   *
%       [8 7 ... 2 1] and normalized                     *
%    [r,j1]=dishlo(v,[ws(1) ... ws(ns)],...              *
%                                 [wl(1) ... wl(nl)]])   *
%        employs user-specified averaging weights.       *
%                                                        *
%                                       02/Jul/2024 Nava *
%*********************************************************
function [r,j1]=dishlo(v,ws,wl)
ns=length(ws);
if ns<3  % Explicit short weights array should have ns>=3 
    nss=ws(1);  % First member of array is length
    % Second member (if any) >0 indicates diminishing
    %                weigths, default=uniform weigths
    if ns==2 && ws(2)>0 
        ws=nss:-1:1; 
    else
        ws=ones(1,nss); 
    end
    ns=nss;
end
ws=ws/sum(ws); % specified or not, weigths are normalized
    
nl=length(wl);
if nl<3  % Long weights array not specified explicitly
    nll=wl(1);   % First member of array is  length
    % Second member (if any) >0 indicates diminishing
    %                weigths, default=uniform weigths
    if nl==2 &&  wl(2)>0
        wl=nll:-1:1; 
    else
        wl=ones(1,nll); 
    end
    nl=nll;
end
wl=wl/sum(wl); % specified or not, weigths are normalized  

[nr,n]=size(v);   % n is length of input time series
if nr>n   % if series is a column vector
    n=nr; ws=ws'; wl=wl';
end

if min(v)<0, v=abs(v); end  % v should be non-negative

nr=n-nl-ns+1;  r=zeros(1,nr);  % nr is length of sta/lta
j1=nl+ns;  
il1=1; il2=nl; is1=il2+1; is2=il2+ns;
for j=1:nr  
    rsh=sum(v(is1:is2).*ws);
    rlo=sum(v(il1:il2).*wl);  % this should never be =0
    r(j)=rsh/rlo;
    il1=il1+1; il2=il2+1; is1=is1+1; is2=is2+1;
end
