Home > atmlab > datasets > +satreaders > modis.m

modis

PURPOSE ^

SATREADERS.MODIS reads modis data

SYNOPSIS ^

function S = modis(file, varargin)

DESCRIPTION ^

 SATREADERS.MODIS reads modis data

 Read MODIS data and output the data in the format common to all
 satreaders.<dataset>.m readers in atmlab. Geodata and time data are
 always retrieved from the data file.

 For info on the common format, see <a href="matlab:help SatDataset/reader">SatDataset/reader</a>.

 IN

   file    string      Path to hdf file
   extra   cell array (optional) extra fields.

 OUT

   data    struct  With fields:
                   time    time in seconds since 00:00 UT
                   lat     latitude in degrees, one column per viewing angle
                   lon     longitude in [-180, 180] degrees, colums as for lat

 FORMAT

   S = satreaders.modis(file, varargin)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

DOWNLOAD ^

modis.m

SOURCE CODE ^

0001 function S = modis(file, varargin)
0002 
0003 % SATREADERS.MODIS reads modis data
0004 %
0005 % Read MODIS data and output the data in the format common to all
0006 % satreaders.<dataset>.m readers in atmlab. Geodata and time data are
0007 % always retrieved from the data file.
0008 %
0009 % For info on the common format, see <a href="matlab:help SatDataset/reader">SatDataset/reader</a>.
0010 %
0011 % IN
0012 %
0013 %   file    string      Path to hdf file
0014 %   extra   cell array (optional) extra fields.
0015 %
0016 % OUT
0017 %
0018 %   data    struct  With fields:
0019 %                   time    time in seconds since 00:00 UT
0020 %                   lat     latitude in degrees, one column per viewing angle
0021 %                   lon     longitude in [-180, 180] degrees, colums as for lat
0022 %
0023 % FORMAT
0024 %
0025 %   S = satreaders.modis(file, varargin)
0026  
0027 % $Id: modis.m 8720 2013-10-21 20:41:39Z gerrit $
0028 % Created by Salomon Eliasson
0029 
0030 core_fields   = {'Scan_Start_Time','Latitude','Longitude'};
0031 extra_fields  = optargs(varargin, {{}});
0032 all_fields    = [core_fields(:); extra_fields(:)];
0033 
0034 for F = all_fields'
0035     S.(F{1}) = hdfread(file,F{1});
0036 end
0037 
0038 % Geodata and time need interpolating to make them the same size as the rest of
0039 % the data
0040 S = interpolate_geodataTime(S);
0041 
0042 % CONSTRUCT time axis with same dimensions
0043 
0044 S.lat = double(S.Latitude);
0045 S.lon = double(S.Longitude);
0046 
0047 % get the verion directly from the filename
0048 D = datasets;
0049 info = D.modis_aqua_L2.find_info_from_granule(file);
0050 date = dayofyear_inverse(str2double(info.year), str2double(info.doy));
0051 S.epoch = round(date2unixsecs(date.year, date.month, date.day));
0052 
0053 %I need to do this after flagging, or the dimensions will be wrong.
0054 % data starts at 1993, I need an offset. Only need one time per scan line
0055 %(first 4 are NaN due to interpolation)
0056 S.time = S.Scan_Start_Time(:,5) + date2unixsecs(1993) -S.epoch;
0057 
0058 S = rmfield(S,core_fields);
0059 
0060 S.version = info.version;
0061 S.path = file;
0062 
0063 S = MaskInvalidGeoTimedataWithNaN(S);
0064 
0065 end
0066 %%%%%%%%%%%%%
0067 % SUBFUNCTION
0068 
0069 function S = interpolate_geodataTime(S)
0070 
0071 % NOTE: Latitude and Longitudes must be interpolated to the right
0072 % dimensions. This is tested to be accurate compared to the MYD03 full
0073 % geodata to within 1e-4 degrees. The alternative is to download the MYD03
0074 % dataset, but that is about 3TB per year.
0075 %
0076 % NOTE: This will destroy all data outside the bounds of the sampled data!
0077 %
0078 % Warning: Slight dodgyness
0079 %
0080 % - I am using 'linear' interpolation on a geoid.
0081 % - When LONGITUDES contain lon < -175 & lon > 175 the satellite clearly
0082 %   passes the  dateline. For such swaths I temporarily switch to 0:360
0083 %   representation, interpolate, and then switch back to -180:180
0084 %   representation.
0085 % - The incorrect interpolation of LATITUDES if the modis passes the poles is deemed
0086 %   small and ignored here since the sampled footprints are about 5km
0087 %   apart. e.g. in the worst case (for this footprint resolution):
0088 %   [89.95 89.95] is interpolated to [89.95 89.95 89.95 89.95 89.95]
0089 %   instead of [89.95 89.98 90 89.98 89.95]
0090 %
0091 
0092 % Use the sampling as described in variable documentation
0093 %Longitude:Cell_Along_Swath_Sampling = 3, 2028, 5 ;
0094 %Longitude:Cell_Across_Swath_Sampling = 3, 1348, 5 ;
0095 
0096 
0097 [X,Y]   = meshgrid(3:5:size(S.Scan_Start_Time,2)*5-2,3:5:size(S.Scan_Start_Time,1)*5-2);
0098 [Xi,Yi] = meshgrid(1:1354,1:size(S.Scan_Start_Time,1)*5);
0099 
0100 fields2DealWith = {'Latitude','Longitude','Scan_Start_Time'};
0101 
0102 lonflag = false;
0103 for F = fields2DealWith
0104     Z = S.(F{1});
0105     
0106     % dateline patch
0107     if strcmp(F{1},'Longitude')
0108         Z(Z<-180) = NaN; % get rid of fill values
0109         if any(Z(:)<-175 ) && any(Z(:)>175) % Fill = -999
0110             % CONVERT LONGITUDES to [0 360] regime if the dateline is crossed
0111             Z=Z+(Z < 0)*360;
0112             lonflag=true;
0113         end
0114     end
0115     
0116     S.(F{1}) = interp2(X,Y,Z,Xi,Yi);
0117     
0118     if strcmp(F{1},'Longitude') && lonflag
0119         S.(F{1}) = S.(F{1})-(S.(F{1}) > 180)*360;
0120     end
0121 end
0122 
0123 end

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005