path=('/datadir/gsat/');
co2_month=[];
lat_month=[];
lon_month=[];
files=dir('*h5');
monthsStr={'01';'02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12'};
monstring={'January','February','March','April','May','June','July','August','September','October','November','December'};
for ii=1:numel(files)
for year=2009:2010
for jj=1:numel(monthsStr)
nameFound=strfind(files(ii).name,['GOSATTFTS' num2str(year) monthsStr{jj}]);
if nameFound
latitude=hdf5read(strcat(path,files(ii).name),'/Data/geolocation/latitude');
longitude=double(hdf5read(strcat(path,files(ii).name),'/Data/geolocation/longitude'));
xco2=double(hdf5read(strcat(path,files(ii).name),'/Data/mixingRatio/XCO2'));
time =hdf5read(strcat(path,files(ii).name),'/scanAttribute/time');
times=length(latitude);
d_deg=1.5;
match=find(latitude < (45+d_deg/2) & latitude >=(15-d_deg/2) & longitude < (-75)+d_deg/2 & longitude >=(-105)-d_deg/2);
lat_G=double(latitude(match));
lon_G=double(longitude(match));
xco2_G=double(xco2(match));
if ~isempty(xco2_G)
co2_month=[co2_month;xco2_G];
lat_month=[lat_month;lat_G];
lon_month=[lon_month;lon_G];
end
d_deg=0.5;
latgrid=15:d_deg:45;
longrid=(-105):d_deg:(-75);
latind=round(lat_month/d_deg-15/d_deg+1);
lonind=round(lon_month/d_deg -(-105/d_deg)+1);
sum_co2=zeros(length(latgrid),length(longrid));
cnt_co2=zeros(length(latgrid),length(longrid));
for j=1:length(co2_month)
sum_co2(latind(j),lonind(j))=sum_co2(latind(j),lonind(j))+co2_month(j);
cnt_co2(latind(j),lonind(j))=cnt_co2(latind(j),lonind(j))+1;
end
grid_co2=sum_co2./cnt_co2;
end
end
end
end
figure(1)
load coast
h1=axesm('MapProjection','eqdcylin',...
'Grid','off',...
'MapLatLimit',[15-d_deg/2 45+d_deg/2],...
'MapLonLimit',[(-105-d_deg/2) (-75+d_deg/2)],...
'flatlimit',[15-d_deg/2 45+d_deg/2],...
'flonlimit',[(-105-d_deg/2) (-75+d_deg/2)],...
'MeridianLabel','off',...
'Frame','on');
plotm(lat,long,'k')
hold on
tightmap on
colormap(Jet);
grid on
h=surfacem(latgrid,longrid,grid_co2);