fid=fopen('mah_smooth_mth-2d.bin','r','b');
mah=fread(fid,'float32');
mah=reshape(mah,360,160,12);
fclose(fid);
fid=fopen('songmiao_mineral_mth-2d.bin','r','b');
songmin=fread(fid,'float32');
songmin=reshape(songmin,360,160,12);
fclose(fid);
fid=fopen('songmiao_solubile_mth-2d.bin','r','b');
songsol=fread(fid,'float32');
songsol=reshape(songsol,360,160,12);
fclose(fid);

X=[.5:1:359.5];
Y=[-79.5:1:79.5];

mah_ann=squeeze(mean(mah,3));
songmin_ann=squeeze(mean(songmin,3));
songsol_ann=squeeze(mean(songsol,3));

filand=find(mah_ann==0);
conv=1e3*365*24*60*60

figure(1), hold off, subplot(1,1,1)
for it=1:3,
  if (it==1), 
    tmp=log10(.04*mah_ann*conv);
    titstr=['Mahowald, annual: 4% solubility (mmol Fe/m^2/y)'];
  end % if it 1
  if (it==2),
    tmp=log10(.04*songmin_ann*conv);
    titstr=['Fan, annual: 4% solubility (mmol Fe/m^2/y)'];
  end % if it 2
  if (it==3),
    tmp=log10(songsol_ann*conv);
    titstr=['Fan, annual: variable solubility (mmol Fe/m^2/y)'];
  end % if it 3
%%
  subplot(3,1,it)
  tmp(filand)=NaN;
  pcolor(X,Y,(tmp)'); 
  caxis([-3.5 .2]);
  shading flat, colorbar
  hold on, hp = gca ; hc = colorbar; axes(hc);
  set(gca,'Ytick',[-3 -2 -1 0 ]);
  set(gca,'YTickLabel',['1e-3';'1e-2';'1e-1';'1e0 ']); axes(hp),
%%
  title(titstr);
end % for it
