function adcp_ell_plots(theYears, dataDirectory, plotDirectory, siteArchive, theSite)
%
%   adcp_ell_plots(theYears, dataDirectory, plotDirectory, siteArchive, theSite)
%
%   Function to create ellipse plots of east and north components by month for
%   hourly averaged and low-passed ADCP measurements for theYears.
%   Output files are named adcp_ell_YYYY.

%Soupy Alexander, 11/6/03

theta=2*pi*(1:64)/64;
theta=[0 theta];

colList = {'r'; 'g'; 'b'; 'c'; 'm'; 'k'; 'r'; 'g'; 'b'; 'c'; 'k'; 'm'};
monthList = {'Jan'; 'Feb'; 'Mar'; 'Apr'; 'May'; 'Jun'; 'Jul'; 'Aug'; 'Sep'; ...
        'Oct'; 'Nov'; 'Dec'};
theSpacing = 10;
theLimit = 15;

%ADCP
timeVar = cell(1,3);
timeVar_lp = cell(1,3);
currVar = cell(1,3);
currVar_lp = cell(1,3);

if strcmp(lower(theSite), 'a')
    depthList = {'5 mbs'; '15 mbs'; 'Deepest Bin'};
else
    depthList = {'5 mbs'; '10 mbs'; 'Deepest Bin'};
end

%Load in the site A data
for indexFile = 1:length(siteArchive)
    if rem(indexFile,10) == 0
        fprintf(['On file number ' num2str(indexFile, '%.0f') ' of ' ...
                num2str(length(siteArchive), '%.0f') '.\n']);
    end
    fileName = siteArchive(indexFile).fileList;
    if isempty(strfind(fileName, 'adc'))
        continue
    end
    currFileA = netcdf(fullfile(dataDirectory, siteArchive(indexFile).fileList), 'nowrite');
    fileDepth = currFileA{'depth'}(:);
    fileTime = singlejd(currFileA{'time'}(:), currFileA{'time2'}(:));
    fileU = currFileA{'u_1205'}(:);
    fileV = currFileA{'v_1206'}(:);
    bads = find(fileU > 1e3 | fileV > 1e3);
    fileU(bads) = NaN;
    fileV(bads) = NaN;
    
    theHit = value2Index(fileDepth, 5);
    theU = fileU(:, theHit);
    try
        [theU_lp, fileTime_lp] = plfilt(theU, fileTime, 6); 
    catch
        continue
    end
    theV = fileV(:, theHit);
    [theV_lp, fileTime_lp] = plfilt(theV, fileTime, 6); 
        
    fileCurr = theU + i.*theV;
    fileCurr_lp = theU_lp + i.*theV_lp;
    
    oldTime = timeVar{1};
    timeVar{1} = [oldTime(:); fileTime(:)]; 
    oldTime_lp = timeVar_lp{1};
    timeVar_lp{1} = [oldTime_lp(:); fileTime_lp(:)];
    
    oldCurr = currVar{1};
    currVar{1} = [oldCurr(:); fileCurr(:)]; 
    oldCurr_lp = currVar_lp{1};
    currVar_lp{1} = [oldCurr_lp(:); fileCurr_lp(:)];
    
    
    if strcmp(lower(theSite), 'a')
        theHit = value2Index(fileDepth, 15);
    elseif strcmp(lower(theSite), 'b')
        theHit = value2Index(fileDepth, 10);
    end
    theU = fileU(:, theHit);
    try
        [theU_lp, fileTime_lp] = plfilt(theU, fileTime, 6); 
    catch
        continue
    end
    theV = fileV(:, theHit);
    [theV_lp, fileTime_lp] = plfilt(theV, fileTime, 6); 
        
    fileCurr = theU + i.*theV;
    fileCurr_lp = theU_lp + i.*theV_lp;
    
    oldTime = timeVar{2};
    timeVar{2} = [oldTime(:); fileTime(:)]; 
    oldTime_lp = timeVar_lp{2};
    timeVar_lp{2} = [oldTime_lp(:); fileTime_lp(:)];
    
    oldCurr = currVar{2};
    currVar{2} = [oldCurr(:); fileCurr(:)]; 
    oldCurr_lp = currVar_lp{2};
    currVar_lp{2} = [oldCurr_lp(:); fileCurr_lp(:)];
    
    
    theHit = find(fileDepth == max(fileDepth));
        theU = fileU(:, theHit);
    try
        [theU_lp, fileTime_lp] = plfilt(theU, fileTime, 6); 
    catch
        continue
    end
    theV = fileV(:, theHit);
    [theV_lp, fileTime_lp] = plfilt(theV, fileTime, 6); 
        
    fileCurr = theU + i.*theV;
    fileCurr_lp = theU_lp + i.*theV_lp;
    
    oldTime = timeVar{3};
    timeVar{3} = [oldTime(:); fileTime(:)]; 
    oldTime_lp = timeVar_lp{3};
    timeVar_lp{3} = [oldTime_lp(:); fileTime_lp(:)];
    
    oldCurr = currVar{3};
    currVar{3} = [oldCurr(:); fileCurr(:)]; 
    oldCurr_lp = currVar_lp{3};
    currVar_lp{3} = [oldCurr_lp(:); fileCurr_lp(:)];
end

for indexYear = 1:length(theYears)
    
    targetYear = theYears(indexYear);
    
    f = figure; orient landscape
    
    for indexVar = 1:3
        
        theTime = timeVar{indexVar};
        theU = real(currVar{indexVar});
        theV = imag(currVar{indexVar});
        
        bads = find(isnan(theU + theV));
        theTime(bads) = [];
        theU(bads) = [];
        theV(bads) = [];
        
        gregTime = gregorian(theTime);
        
        if isempty(find(gregTime(:,1) == targetYear))
            close
            clear f
            continue
        end
        
        subplot(3,1,indexVar), hold on
        
        for indexMonth = 1:12
            inYear = find(gregTime(:,1) == targetYear & ...
                gregTime(:,2) == indexMonth);
            
            theInTime = theTime(inYear);
            theInU = theU(inYear);
            theInV = theV(inYear);
            
            bads = find(isnan(theInU + theInV));
            theInTime(bads) = [];
            theInU(bads) = [];
            theInV(bads) = [];
            
            
            if length(theInU) > 10*24
                [majAx, majAz, minAx, minAz, theEllip] = pcaben(theInU, theInV);
                x2=(majAx*cos(theta))';
                y2=(minAx*sin(theta))';
                theAngle=-1*majAz*pi/180+pi/2;
                xx2 = x2(:).*cos(theAngle) - y2.*sin(theAngle);
                yy2 = x2(:).*sin(theAngle) + y2.*cos(theAngle);
                
                meanU = mean(theInU);
                meanV = mean(theInV);
                
                q = quiver(theSpacing*indexMonth, 0, meanU, meanV, colList{indexMonth});
                hold on
                plot(theSpacing*indexMonth+xx2 + meanU, yy2 + meanV, colList{indexMonth});
                
                l = line([-10 theSpacing*15], [0 0]);
                set(l, 'color', 'k')
            end  
        end
        
        
        axis equal
        axis([0 theSpacing*13 -15 15])
        set(gca, 'xtick', [theSpacing:theSpacing:12*theSpacing]);
        set(gca, 'xticklabel', monthList)
        set(gca, 'ytick', [-15:5:15]);
        ylabel('cm/s')
        box('on')
        
        title([depthList{indexVar}])
    end
    
    outTitle = {['Site ' upper(theSite) '- ' num2str(targetYear, '%.0f') ': Monthly Mean Current and Variability']; 
        'ADCP Hour-Averaged Data'};
    suptitle(outTitle)
    
    printTitle = [ lower(theSite) '_ell_h_' num2str(targetYear, '%.0f')];
    print('-dpdf', fullfile(plotDirectory, printTitle))
    
    close
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    
    figure, orient landscape
    
    for indexVar = 1:3
        
        theTime = timeVar_lp{indexVar};
        theU = real(currVar_lp{indexVar});
        theV = imag(currVar_lp{indexVar});
        
        gregTime = gregorian(theTime);
       
        
        subplot(3,1,indexVar), hold on
        
        for indexMonth = 1:12
            inYear = find(gregTime(:,1) == targetYear & ...
                gregTime(:,2) == indexMonth);
            
            theInTime = theTime(inYear);
            theInU = theU(inYear);
            theInV = theV(inYear);
            
            bads = find(isnan(theInU + theInV));
            theInTime(bads) = [];
            theInU(bads) = [];
            theInV(bads) = [];
            
            
            if length(theInU) > (10*4)
                [majAx, majAz, minAx, minAz, theEllip] = pcaben(theInU, theInV);
                x2=(majAx*cos(theta))';
                y2=(minAx*sin(theta))';
                theAngle=-1*majAz*pi/180+pi/2;
                xx2 = x2(:).*cos(theAngle) - y2.*sin(theAngle);
                yy2 = x2(:).*sin(theAngle) + y2.*cos(theAngle);
                
                meanU = mean(theInU);
                meanV = mean(theInV);
                
                q = quiver(theSpacing*indexMonth, 0, meanU, meanV, colList{indexMonth});
                hold on
                plot(theSpacing*indexMonth+xx2 + meanU, yy2 + meanV, colList{indexMonth});
                
                l = line([-10 theSpacing*15], [0 0]);
                set(l, 'color', 'k')
            end
        end
        
        
        axis equal
        axis([0 theSpacing*13 -15 15])
        set(gca, 'xtick', [theSpacing:theSpacing:12*theSpacing]);
        set(gca, 'xticklabel', monthList)
        set(gca, 'ytick', [-15:5:15]);
        ylabel('cm/s')
        box('on')
        
        title([depthList{indexVar}])
    end
    
    outTitle = {['Site ' upper(theSite) '- ' num2str(targetYear, '%.0f') ': Monthly Mean Current and Variability']; 
        'ADCP Low-Passed Data'};
    suptitle(outTitle)
    
    printTitle = [lower(theSite) '_ell_lp_' num2str(targetYear, '%.0f')];
    print('-dpdf', fullfile(plotDirectory, printTitle))
    
    close
end