year=2001:2014;
lat=-89.5:1:89.5;
lon=-179.5:1:179.5;
month=1:168;

AOD_temp=ncread('/Users/yangyang/Aerosol_Variaiton/IRUN/AOD/ncfile/AOD_MODIS.nc','AOD_MONTH');

AOD_MONTH=zeros(360,180,168);

AOD_temp(AOD_temp<0)=nan;

% Deseasonal
for i=1:360
    for j=1:180
        for monthID=1:12
            temp=nanmean(AOD_temp(i,j,monthID,:),4);
            for yearID=1:14
                AOD_temp(i,j,monthID,yearID)=AOD_temp(i,j,monthID,yearID)-temp;
            end
        end
    end
end

% Change size
for i=1:360
    for j=1:180
        n=0;
        for yearID=1:length(year)
            for monthID=1:12
                n=n+1;
                AOD_MONTH(i,j,n)=AOD_temp(i,j,monthID,yearID);
            end
        end
    end
end




%for i=1:360
%    for j=1:180
%        for n=1:168
%            if AOD_MONTH(i,j,n)<0
%                AOD_MONTH(i,j,:)=-999;
%            end
%        end
%    end
%end
                



save('AOD_MODIS_monthly.mat','AOD_MONTH')

nccreate('AOD_MODIS_monthly.nc','Month','Dimensions',{'time',Inf},'Format','classic')
ncwrite('AOD_MODIS_monthly.nc','Month', month)
ncwriteatt('AOD_MODIS_monthly.nc','Month','long_name','Month')
nccreate('AOD_MODIS_monthly.nc','lat','Dimensions',{'lat',180},'Format','classic')
ncwrite('AOD_MODIS_monthly.nc','lat', lat)
ncwriteatt('AOD_MODIS_monthly.nc','lat','units','degrees_north')
ncwriteatt('AOD_MODIS_monthly.nc','lat','long_name','Latitude')
nccreate('AOD_MODIS_monthly.nc','lon','Dimensions',{'lon',360},'Format','classic')
ncwrite('AOD_MODIS_monthly.nc','lon', lon)
ncwriteatt('AOD_MODIS_monthly.nc','lon','units','degrees_east')
ncwriteatt('AOD_MODIS_monthly.nc','lon','long_name','Longitude')
nccreate('AOD_MODIS_monthly.nc','AOD_MONTH','Dimensions',{'lon',360,'lat',180,'month',length(month)},'Format','classic')
ncwrite('AOD_MODIS_monthly.nc','AOD_MONTH', AOD_MONTH)
ncdisp('AOD_MODIS_monthly.nc')