  function [bfdata]=read_rttov_emis(bfemis_dir,year,month,ich)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  function to read the UW/CIMSS global IR land s
%  surface BF emissivity data and location files
%  
%
%    CREATED by:  Eva Borbas UW/CIMSS  - july 2009
%
%    Note: This algorithm was released to be useful, 
%          but without any warranty.
%
%
%    INPUT: bfemis_dir  - directory of the emissivity database
%           year       -  number of  year of the data
%           month      -  number of the month
%           ich        -  number of the hinge point (channel)
%
%    OUTPUT: data       - structure includes emissivity, quality flag, wavenumbers, 
%			   latitude, longitude
%
% 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%addpath /home/evab/mfiles/
%bfemis_dir='/data4/emis/global_data/cdf/';
%year=2009
%month=7
%ich=1

% scale factors from the netcdf file
emis_sf=0.001;
flag_sf=1;
wavenum_sf=0.1;


nemis=10;

craft='aqua';

day_str={'001', '032', '060', '091', '121', '152', '182', '213', '244', '274', '305', '335'};
day_str_leap={'001', '032', '061', '092', '122', '153', '183', '214', '245', '275', '306', '336'};

if mod(year,4) == 0 
nday = day_str_leap;
else
nday= day_str;
end

%----------------
% read emis netcdf file
%----------------

bf_filename = ['global_emis_inf10_monthFilled_MYD11C3.A'  num2str(year) ...
                  char(nday(month)) '.041.nc']

bfem = read_netcdf([char(bfemis_dir) char(bf_filename)]);
disp(['Read BF File: ' bf_filename]);

%----------------
% emis data
%----------------

% can be looped here for the all 10 channels
%for ich=1:10

eval(sprintf('bfdata.emis%g=double(bfem.emis%g).*emis_sf;',ich,ich));
eval(sprintf('bfdata.emis%g(find(bfdata.emis%g < 0))=NaN;',ich,ich));
eval(sprintf('bfdata.emis%g;',ich,ich));

%figure;clf;
%eval(sprintf('imagesc(bfdata.emis%g);hold on;colorbar;',ich));

%end

%----------------
% quality flag
%----------------
bfdata.emis_flag=double(bfem.emis_flag).*flag_sf;
%figure);clf;imagesc(bfdata.emis_flag);colorbar;

%----------------
% wavenumbers
%----------------
bfdata.wavenum=double(bfem.wavenumber).*wavenum_sf;

clear bfem


%----------------
% read emis location  netcdf file
%----------------

bfloc_filename = ['global_emis_inf10_location.nc'];

bfdata_loc = read_netcdf([bfemis_dir bfloc_filename]);
disp(['Read BF File: ' bfloc_filename]);

%----------------
% lat/lon data
%----------------

bfdata.lat=double(bfdata_loc.lat);
bfdata.lon=double(bfdata_loc.lon);


