Documentation of interp_oi2CCM

Global Index (short | long) | Local contents | Local Index (short | long)

Help text

  interp2CCM.m

  This file should load the sstoi data in, interpolate it to 
  the CCM grid, redefine the time axis, and git movin'

Cross-Reference Information

Listing of script interp_oi2CCM

clear
cd /home/disk/hayes/dvimont/nco

%  Open files, and read variables

filoi  = ['/home/disk/hayes/dvimont/nco/sstoi8298.nc'];
filccm = ['/home/disk/hayes/dvimont/nco/sstccmnew.nc'];
filccm_orig = ['/home/disk/hayes/dvimont/nco/sst_jan1978-sep1993.nc'];

ccmid = ncmex('OPEN', filccm, 'nowrite');
  [latccm, status] = ncmex('VARGET', ccmid, 'lat', 0, -1, 0);
  [lonccm, status] = ncmex('VARGET', ccmid, 'lon', 0, -1, 0);
  [sstccm, status] = ncmex('VARGET', ccmid, 'SST', [0 0 0], [-1 -1 -1], 0);
ccmid = ncmex('CLOSE', ccmid);

%  Fix time (NOTE:  The year_days_leap makes the time axis consistent
%     for FERRET, but I think one should use a 365 day year for CCM)

time_new = [];
year_days = [1 32 60 91 121 152 182 213 244 274 305 335];

%year_days_leap = [year_days, 365+year_days 730+[year_days(1:2) 1+year_days(3:12)] 1096+year_days];
%for i = 1:4;
%  time_new = [time_new, [(365*3+366)*(i-1) + year_days_leap]];

for i = 1:16;
   time_new = [time_new, 365*(i-1)+year_days];
end;
time_new = [time_new, [365*16 + year_days(1:7)]];

%  Redo date axis on netcdf file:

date_new = [];
day =  [16 14 16 15 16 15 16 16 15 16 15 16];
for i = 1982:1997
  for j = 1:12
    date_new = [date_new; (10000*i+100*j+day(j))];
  end
end

for j = 1:7
  date_new = [date_new; (19980000+100*j+day(j))];
end

date_new = [198112116; date_new];

nbdate_new = 811216;

%  Write new data to netcdf file

filccm = ['sst_nov1981-jul1998.nc']
ccmid = ncmex('OPEN', filccm, 'NC_WRITE');

[timeid, rcode] = ncmex('VARID', ccmid, 'time');
[time, status] = ncmex('VARGET', ccmid, timeid, 0, -1, 0);
status = ncmex('VARPUT', ccmid, timeid, 0, -1, time_new, 0);

[dateid, rcode] = ncmex('VARID', ccmid, 'date');
[date, status] = ncmex('VARGET', ccmid, dateid, 0, -1, 0);
status = ncmex('VARPUT', ccmid, dateid, 0, -1, date_new', 0);

%  Change beginning day for dataset

[datatype, len, status] = ncmex('ATTINQ', ccmid, 'time', 'units');
[value, status] = ncmex('ATTGET', ccmid, 'time', 1)
value_new = ['days since 0081-11-14 00:00:00']
status = ncmex('ATTPUT', ccmid, 'time', 'units', datatype, len+1, value_new)

ccmid = ncmex('CLOSE', ccmid);

nc = netcdf(filccm, 'write');
nc{'nbdate'}(:) = nbdate_new;
nc = close(nc);

%  Load big SST file from sstoi8298.nc, and interpolate

oiid = ncmex('OPEN', filoi, 'nowrite');
  [sstid, rcode] = ncmex('VARID', oiid, 'data');
  [sstoi, status] = ncmex('VARGET', oiid, sstid, [0 0 0], [-1 -1 -1], 0);
  [add_offset, status] = ncmex('ATTGET', oiid, 'data', 'add_offset');
  [scale_factor, status] = ncmex('ATTGET', oiid, 'data', 'scale_factor');
  [latoi, status] = ncmex('VARGET', oiid, 'lat', 0, -1, 0);
  [lonoi, status] = ncmex('VARGET', oiid, 'lon', 0, -1, 0);
  status
oiid = ncmex('CLOSE', oiid);

sstoi = sstoi(:,:,1:199) * scale_factor + add_offset;

[nx, ny, ntm] = size(sstoi);

%  Flip the y=axis upside down so the grid matches sstccm

if (latoi(1) > latoi(length(latoi)))
  latoi = fliplr(latoi);
  for i = 1:ntm;
    sstoi(:,:,i) = reshape(fliplr(squeeze(sstoi(:,:,i))),nx, ny, 1);
  end
end 

%  Add a few degrees on either side for better interpolation

sstoi = reshape(sstoi, nx, ny*ntm);
sstoi = [sstoi((nx-4):nx,:); sstoi; sstoi(1:5,:)];
sstoi = reshape(sstoi, (nx+10), ny, ntm);
lonoi = [(lonoi((nx-4):nx)-360) lonoi (lonoi(1:5)+360)];

sstoi = shiftdim(sstoi, 2);
for i = 1:ntm;
  hello(['iteration ' num2str(i)])
  [temp, n] = shiftdim(sstoi(i,:,:));
  temp = interp2(latoi, lonoi, temp, latccm', lonccm);
  sstccm_new(i,:,:) = shiftdim(temp, -n);
end
sstoi = shiftdim(sstccm_new, 1);

[nx, ny, ntm] = size(sstoi);

%  Change ice flag to -1.8 on sstoi:

sstoi = shiftdim(sstoi, 2);
cur_flag = min(min(min(sstoi)));
new_flag = -1.8;
for i = 1:ntm
  [temp, n] = shiftdim(sstoi(i,:,:));
  temp = reshape(temp, nx*ny, 1);
  flag = find(temp == cur_flag);
  temp(flag) = new_flag .* ones(length(flag), 1);
  temp = reshape(temp, nx, ny);
  temp = shiftdim(temp, -n);
  sstoi(i,:,:) = temp;
end
sstoi = shiftdim(sstoi, 1);

%  Write sstoi to the new dataset

ccmid = ncmex('OPEN', filccm, 'NC_WRITE');
  [sstid, rcode] = ncmex('VARID', ccmid, 'SST');
  [status] = ncmex('VARPUT', ccmid, sstid, [0 0 0], [-1 -1 -1], sstoi, 0);
  status
ccmid = ncmex('CLOSE', ccmid);