function [figNo] = plotReconsEmission(hDir, info, periods, USBndryFile, ...
                                      topcorner, botcorner, emissionsTrue2D, ...
                                      emissionsRecons2D, figNo, ...
                                      visualize, saveStuff, oDir, fs, ...
                                      marker, ms, mfc, lw)
%plotReconsEmission() Function to plot reconstructed emissions along with 
%sensor locations
%   [figNo] = plotReconsEmission(hDir, periods, USBndryFile, ...
%                                topcorner, botcorner, emissionsTrue2D, ...
%                                emissionsRecons2D, visualize, saveStuff, oDir)
%
%  Function to plot emissions in the Lower48, along with US bndry and
%  sensor locations.
%
%  INPUTS :
%   hDir        ; directory where info about tower locations are kept
%   info        : struct with info
%   periods     : [start, end] of periods of interest
%   USBndryFile : name of file with US bndry polygon
%   topcorner, botcorner : [lat, lon] of topleft, bot-right corner of
%                          domain
%   emissionsTrue2D : 3D matrix, of emissions on a 64^2 grid
%   emissionsRecons2D : ditto, but for estimated emissions
%   figNo        : figiNo to start numbering
%   visualize    : plto stuff?
%   saveStuff    : save to file?
%   oDir         : output directory
%   fs           : font size
%   marker, ms, mfc : symbol for sensor, marker size & marker fill color
%   lw           : linewidth
%
%  OUTPUTS :
%   figNo   : figure number to use downstream
%
%  Jaideep Ray, 6/27/2012
% -------------------------------------------------------------------------

% ---- Load the tower locations
[header, LatLonOfTowers] = hdrload([hDir 'MeasuredCO2/LatLongOfTowers_Lower48.dat']);
LatLonOfTowers = LatLonOfTowers(:, 2:3) ;

% ---- Scale the LatLon so that we can plot it on a scale that goes 1:29, 1:65
LatLonOfTowers(:, 1) = topcorner(1) - LatLonOfTowers(:, 1) + 1 ; % lats
LatLonOfTowers(:, 2) = LatLonOfTowers(:, 2) - topcorner(2) + 1 ; % longs
periodList = periods(1) : periods(2) ;

% ---- Read in US bndry
usbndry = load(USBndryFile) ;
usbndry(:, 1) = usbndry(:, 1) - topcorner(2) + 1 ;
usbndry(:, 2) = topcorner(1) - usbndry(:, 2) + 1 ;

% Plot the reconstructed emissions
clim = [-5.0e-2  0.5] ;  

if (visualize == true) 
    
    % ---- Loop over all periods
    for i = info.buffer + 1 : (info.nPeriods - info.buffer)
        
        % ---- True emissions
        figNo = figNo + 1 ; figure(figNo) ;
        imagesc(flipud(emissionsTrue2D(:, :, i)), clim) ; hold on ; colorbar ;
        periodNo = periodList(i - info.buffer) ;
        str = sprintf('True emissions in 8-day period %d [microMoles/m^2/sec]', ...
            periodNo) ;
        title(str, 'FontSize', fs) ;
        plot(LatLonOfTowers(:,2), LatLonOfTowers(:, 1), marker, ...
            'MarkerSize',ms, 'MarkerFaceColor', mfc) ; hold on ;
        
        % ---- Plot the US bndry
        plot(usbndry(:, 1), usbndry(:, 2), '-w', 'LineWidth', lw) ;
        
        % Save file to disk
        if (saveStuff == true)
            str = sprintf('TrueEmissions_P%02d.jpg', periodNo);
            print('-djpeg', [oDir '/' str]);
        end
        
        % Increase fonts
        set(gca, 'FontSize', fs);
        
        % ----- Now do reconstructed emissions
        figNo = figNo + 1 ; figure(figNo) ;
        imagesc(flipud(emissionsRecons2D(:, :, i)), clim) ;  hold on ; colorbar ;
        str = sprintf('Reconstructed emissions in 8-day period %d [microMoles/m^2/sec] ', ...
            periodNo) ;
        title(str, 'FontSize', fs) ;
        plot(LatLonOfTowers(:,2), LatLonOfTowers(:, 1), marker, ...
            'MarkerSize',ms, 'MarkerFaceColor', mfc) ; hold on ;
        
        % ---- Plot the US bndry
        plot(usbndry(:, 1), usbndry(:, 2), '-w', 'LineWidth', lw) ;
        
        % Save file to disk
        if (saveStuff == true)
            str = sprintf('ReconsEmissions_P%02d.jpg', periodNo);
            print('-djpeg', [oDir '/' str]);
        end
        
        % Increase fonts
        set(gca, 'FontSize', fs);
    end

end

end

