% Convert a GrADS-formatted QuikSCAT climatology kept by Jan Hafner at IPRC
% to netcdf. Makes use of parse_ctl.m read_subvol.m, and 
%
% Recalculated time with new parse_ctl.m on 10/23/06,
% which should not interfere with u,v. SPdeS

%fname='/home/helios/xie/sapporo/CLIM/orog/qscatC.ctl';
fname='/home/honuea/hafner/xie/CLIM/qscatC.ctl'

[nlon,lon,nlat,lat,nlev,lev,ntime,time,dset,dtitle,...
 machineformat,undef,nvars,name,vlev,units,sequential]...
    =parse_ctl(fname);
if sequential==1
  disp('Warning: file reads in read_subvol must be modified for sequential binary.')
end

% open the file for reading in the correct architecture
fid = fopen(dset,'r',machineformat);
nbyte=4; % GrADS is single-precision by default

fivar=1:length(name); % read all

% initialize arrays
for m=fivar
  eval([name{m} '=NaN+zeros(nlon,nlat,vlev(m),ntime);'])
end
frewind(fid);
% read variables in order:
% (fastest varying) lon,lat,lev,variable,time (slowest varying)
for l=1:ntime % loop time
  for m=1:length(name)
    vol=NaN+zeros(nlon,nlat,vlev(m));
    if ismember(m,fivar)
      vol(:)=read_subvol(fid,1,nlon,1,nlat,1,vlev(m),nlon,nlat,vlev(m),'single');
      vol(vol==undef)=NaN;
      eval([name{m} '(:,:,:,l)=vol;']);
    else
      status=fseek(fid,nbyte*nlat*nlon*vlev(m),0); % skip over the 3D volume
    end
  end
end

% write QuikSCAT climatology  netcdf
nc=netcdf('qscatC.nc','clobber');

nc('lon') = nlon;
nc{'lon'} = {'lon'};
nc{'lon'}.long_name = 'longitude';
nc{'lon'}.units = 'degrees_east';
nc{'lon'}(:)=lon;
nc('lat') = nlat;
nc{'lat'} = {'lat'};
nc{'lat'}.long_name = 'latitude';
nc{'lat'}.units = 'degrees_north';
nc{'lat'}(:)=lat;
nc('lev') = nlev;
nc{'lev'} = {'lev'};
nc{'lev'}.long_name = 'vertcal level';
nc{'lev'}.units = 'surface';
nc{'lev'}(:)=0;
nc('time') = ntime;
nc{'time'} = {'time'};
nc{'time'}.units='calendar month';
nc{'time'}.first_month='January';
nc{'time'}(:) = 1:12;

nc{'u'}=ncfloat('lon','lat','lev','time');
nc{'u'}.long_name = 'u wind component';
nc{'u'}.units = 'm s-1';
nc{'u'}(:)=u;
nc{'v'}=ncfloat('lon','lat','lev','time');
nc{'v'}.long_name = 'v wind component';
nc{'v'}.units = 'm s-1';
nc{'v'}(:)=v;
nc{'w'}=ncfloat('lon','lat','lev','time');
nc{'w'}.long_name = 'wind speed';
nc{'w'}.units = 'm s-1';
nc{'w'}(:)=sp;
nc{'taux'}=ncfloat('lon','lat','lev','time');
nc{'taux'}.long_name = 'u pseudostress component';
nc{'taux'}.units = 'm2 s-2';
nc{'taux'}(:)=tx;
nc{'tauy'}=ncfloat('lon','lat','lev','time');
nc{'tauy'}.long_name = 'v pseudostress component';
nc{'tauy'}.units = 'm2 s-2';
nc{'tauy'}(:)=ty;

nc.title= 'QuikSCAT (August1999-August2006) climatological surface wind and stress for the global ocean from RSS, F. Wentz)';
nc.resolution='0.25 x 0.25 degree';
nc.history= 'Created by make_qscat_nc.m, which converts a QuikSCAT climatology GrADS file by Jan Hafner on 9/21/2006. SPdeS';

nc=close(nc);
