source: altifloat/matlab_toolbox/script_prepare_exp_modele.m @ 199

%% STEP ONE : produce the meshfile

%meshgrid de Dan

meshfile='../exp_forw_dan/meshgrid_dan.dat';
make_meshgrid(meshfile,33,31,36.92572,36.22729,489,438);

%% STEP TWO : transform and interpolate aviso file to YAO file
dir_uvfiles='../../data/MODEL_CYCO/';
aviso_netcdf={[dir_uvfiles 'uveta_cyco1_2013090112.nc']...
              [dir_uvfiles 'uveta_cyco1_2013090212.nc'];};

it=1:24:1+24*(length(aviso_netcdf)-1);

%meshmask='../exp_forw_dan/meshgrid_dan.dat';
outfile='../exp_forw_dan/uv_bck_dan.dat';
outmask='../exp_forw_dan/mask_dan.dat';

dt=1;

interp_bck(aviso_netcdf,it,dt,meshfile,outfile,outmask);

%% Lissage
outfile_s='../exp_forw_dan/uv_bck_dan_s.dat';


uvb=load(outfile);
meshg=load(meshfile);
lon=unique(meshg(:,1));
lat=unique(meshg(:,2));
nlon=length(lon);
nlat=length(lat);
ntime=length(unique(uvb(:,1)));
Ub=reshape(uvb(:,end-1),nlat,nlon,ntime);
Vb=reshape(uvb(:,end),nlat,nlon,ntime);

Ub_s=nan*ones(size(Ub));
Vb_s=nan*ones(size(Vb));

K=ones(95,95);
K=K./sum(K(:));

for it=1:ntime
Ub_s(:,:,it)=conv2(Ub(:,:,it),K,'same');
Vb_s(:,:,it)=conv2(Vb(:,:,it),K,'same');
end

uvb_s=uvb;
uvb_s(:,end-1) = Ub_s(:);
uvb_s(:,end)   = Vb_s(:);
dlmwrite(outfile_s,uvb_s,'\t');

%% STEP THREE : Prepare init
drifterdir='../altifloat/';
initfile='init.dat';
fname={'a300234011043340.mat','a300234011120770.mat','a300234011249440.mat'};
lon_drift=nan*ones(1,length(fname));
lat_drift=lon_drift;
MM=[];
stime_num=datenum('20130901T120000','yyyymmddTHHMMSS');
dt=1/24; %one hour


for j=1:length(fname)
       load([drifterdir fname{j}]);
       [x,y,t,M]=obs2yao(lon,lat,stime_num,stime_num,dt,meshfile,j);
       MM=[MM;M];
    
end % for j
dlmwrite(['../exp_forw_dan/' initfile],MM,'\t');

%% STEP OBS : modify obs file
init=load(['../exp_forw_dan/' initfile]);
obs=load(['../exp_forw_dan/rfloat_total.dat']);

ilon=find(init(:,1)==1);
ilat=find(init(:,1)==2);

init=[init(ilon,2)-1 init(ilon,3)-1 init(ilat,4) init(ilon,4)];

istep=5;
I=find(obs(:,end)>1e-8);
obs=obs(I,:);

obs=[init;obs];
obs=sortrows(obs,[2 1]);

t=unique(obs(:,1));
tech=t(1:istep:end);

I=find(ismember(obs(:,1),tech));
obs=obs(I,:);
dlmwrite(['../exp_forw_dan/obs_jum.dat'],obs,'\t');

%% Create bck
aviso_netcdf={[dir_uvfiles 'uveta_cyco1_2013092812.nc']...
              [dir_uvfiles 'uveta_cyco1_2013092912.nc'];};

it=1:24:1+24*(length(aviso_netcdf)-1);

%meshmask='../exp_forw_dan/meshgrid_dan.dat';
outfile='../exp_forw_dan/uv_bck2_dan.dat';
outmask='../exp_forw_dan/mask_dan2.dat';

dt=1;

interp_bck(aviso_netcdf,it,dt,meshfile,outfile,outmask);