Newer
Older
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% template_get_adcp_data.m
% -------------------------------
% Author : Pierre ROUSSELOT - IRD (pierre.rousselot@ird.fr)
% Jeremie HABASQUE - IRD (jeremie.habasque@ird.fr)
% -------------------------------
% INPUTS:
% - binary raw file with .000 extension
% OUTPUTS:
% - U and V fields interpolated on a regulard grid, filtered and subsampled
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
close all; clear
addpath(genpath('../ADCP_mooring_data_processing'));
addpath('/home/proussel/Documents/OUTILS/TOOLS/nansuite'); % NaNSuitePath
rawfile = 'FR31/FR30_000.000'; % binary file with .000 extension
fpath_output = './FR31/0-0/'; % Output directory
cruise.name = 'PIRATA'; % cruise name
mooring.name = '0W0N'; % '0N10W'
mooring.lat = '0000.176'; % latitude [°']
mooring.lon = '-0003.920'; % longitude [°']
adcp.sn = 22546; % ADCP serial number
adcp.type = '150 khz Quartermaster'; % Type : Quartermaster, longranger
adcp.direction = 'up'; % upward-looking 'up', downward-looking 'dn'
adcp.instr_depth = 300; % nominal instrument depth
instr = 1; % this is just for name convention and sorting of all mooring instruments
offset_pres_sens = 1; % offset in m between pressure sensore and ADCP
% NCFILE info
d_fillvalue = -9999;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
latDegInd = strfind(mooring.lat,'');
lonDegInd = strfind(mooring.lon,'');
if strfind(mooring.lon,'-')
mooring.lat = str2double(mooring.lat(1:latDegInd-1))-str2double(mooring.lat(latDegInd+1:end-1))/60;
else
mooring.lat = str2double(mooring.lat(1:latDegInd-1))+str2double(mooring.lat(latDegInd+1:end-1))/60;
end
if strfind(mooring.lon,'-')
mooring.lon = str2double(mooring.lon(1:lonDegInd-1))-str2double(mooring.lon(lonDegInd+1:end-1))/60;
else
mooring.lon = str2double(mooring.lon(1:lonDegInd-1))+str2double(mooring.lon(lonDegInd+1:end-1))/60;
end
clock_drift = clock_drift/3600; % convert into hrs
%% Read rawfile
disp('****')
raw_file = [fpath_output, mooring.name '_' num2str(adcp.sn) '_instr_' sprintf('%02d',instr) '_raw.mat'];
if exist(raw_file)
fprintf('Read %s\n', raw_file);
load(raw_file)
else
fprintf('Read %s\n', rawfile);
raw = read_os3(rawfile,'all');
save(raw_file,'raw','-v7.3');
end
%% Correct clock drift
time0 = julian(raw.juliandate);
clockd = linspace(0, clock_drift, length(time0));
raw.juliandate = raw.juliandate - clockd / 24;
disp('****')
disp('Correct clock drift')
%% Plot pressure and temperature sensor
figure;
subplot(2,1,1)
plot(raw.pressure);
title('Pressure sensor');
ylabel('Depth [m]');
xlabel('Time index');
subplot(2,1,2)
plot(raw.temperature);
title('Temperature sensor');
pierre.rousselot_ird.fr
committed
saveas(gcf,[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Pressure_Temp_sensor'],'fig')
%% Determine first and last indiced when instrument was at depth (you can do this by plotting 'raw.pressure' for example
disp('****')
first = input('Determine first indice when instrument was at depth (with pres/temp plot): ');
disp('****')
last = input('Determine last indice when instrument was at depth (with pres/temp plot): ');
%% Correct velocity data with external T/S sensor
%% Extract data
freq = raw.config.sysconfig.frequency;
u2 = squeeze(raw.vel(:,1,first:last));
v2 = squeeze(raw.vel(:,2,first:last));
w = squeeze(raw.vel(:,3,first:last));
vel_err = squeeze(raw.vel(:,4,first:last));
pierre.rousselot_ird.fr
committed
corr = squeeze(mean(raw.cor(:,4,first:last),2));
ea = squeeze(mean(raw.amp(:,:,first:last),2));
pierre.rousselot_ird.fr
committed
pg = squeeze(raw.pg(:,4,first:last));
ang = [raw.pitch(first:last) raw.roll(first:last) raw.heading(first:last)];
soundspeed = raw.soundspeed(first:last);
temp = raw.temperature(first:last);
press = raw.pressure(first:last);
if press < 0
press = -press;
end
nbin = raw.config.ncells; % number of bins
bin = 1:nbin; % bin number
blen = raw.config.cell; % bin length
tlen = raw.config.pulse; % pulse length
lag = raw.config.lag; % transmit lag distance
blnk = raw.config.blank; % blank distance after transmit
fbind = raw.config.bin1distance;% middle of first bin distance
dt = (time(2)-time(1))*24; % Sampling interval in hours
EA0 = round(mean(ea(nbin,:)));
%% Calculate Magnetic deviation values
[a,~] = gregorian(raw.juliandate(1));
magnetic_deviation_ini = magdev(mooring.lat,mooring.lon,0,a+(raw.juliandate(1)-julian(a,1,1,0,0,0))/365.25);
[a,~] = gregorian(raw.juliandate(end));
magnetic_deviation_end = magdev(mooring.lat,mooring.lon,0,a+(raw.juliandate(end)-julian(a,1,1,0,0,0))/365.25);
rot = (magnetic_deviation_ini+magnetic_deviation_end)/2;
mag_dev = linspace(magnetic_deviation_ini, magnetic_deviation_end, length(time0));
pierre.rousselot_ird.fr
committed
mag_dev = mag_dev(first:last);
disp('****')
disp('Correct magnetic deviation')
M(1,1) = cos(mag_dev(ii));
M(1,2) = -sin(mag_dev(ii));
M(2,1) = sin(mag_dev(ii));
M(2,2) = cos(mag_dev(ii));
vvel(1,:) = u2(:,ii);
vvel(2,:) = v2(:,ii);
vvvel = M * vvel;
u(:,ii) = vvvel(1,:);
v(:,ii) = vvvel(2,:);
%% Correct percent good: Exclude data with percent good below prct_good
figure;
colormap jet;
pcolor(pg);
shading flat;
pierre.rousselot_ird.fr
committed
title('Percent good of the bins'); colorbar;
ylabel('Bins');
pierre.rousselot_ird.fr
committed
saveas(gcf,[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Perceent_good'],'fig')
disp('****')
prct_good = input('Determine percent good threshold (generally 20): ');
igap = find(pg<prct_good); % Exclude data with percent good below prct_good
u(igap) = NaN;
v(igap) = NaN;
w(igap) = NaN;
vel_err(igap) = NaN;
pierre.rousselot_ird.fr
committed
ea(igap) = NaN;
corr(igap) = NaN;
%% Correction data with high attitude (pitch/roll)
figure;
subplot(2,1,1)
plot(abs(ang(:,1)));
hold on
plot(10*ones(1,length(ang(:,1))),'--r');
hold off
title('Attitude sensor');
pierre.rousselot_ird.fr
committed
xlabel('Time index');
axis([0 length(ang(:,1)) 0 20])
pierre.rousselot_ird.fr
committed
subplot(2,1,2)
plot(abs(ang(:,2)));
hold on
plot(10*ones(1,length(ang(:,1))),'--r');
hold off
pierre.rousselot_ird.fr
committed
xlabel('Time index');
axis([0 length(ang(:,1)) 0 20])
pierre.rousselot_ird.fr
committed
saveas(gcf,[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Attitude'],'fig')
disp('****')
disp('Delete high attitude ADCP data (>=10)');
pierre.rousselot_ird.fr
committed
high_pitch = find(abs(ang(:,1))>=10);
high_roll = find(abs(ang(:,2))>=10);
if ~isempty(high_pitch) || ~isempty(high_roll)
high_att = input('Do you want to delete them? 1/0 [1]');
pierre.rousselot_ird.fr
committed
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
if isempty(high_att)
high_att = 1;
end
if high_att == 1
high_att = union(high_pitch,high_roll);
offset(high_att) = NaN;
ang(high_att,1) = NaN;
ang(high_att,2) = NaN;
ang(high_att,3) = NaN;
mag_dev(high_att) = NaN;
u(high_att) = NaN;
v(high_att) = NaN;
w(high_att) = NaN;
vel_err(high_att) = NaN;
ea(high_att) = NaN;
corr(high_att) = NaN;
pg(high_att) = NaN;
time(high_att) = NaN;
z_bins(high_att) = NaN;
depth(high_att) = NaN;
temp(high_att) = NaN;
soundspeed(high_att) = NaN;
end
end
%% amplitude of the bins / Correction ADCP's depth
figure;
colormap jet;
pcolor(ea);
shading flat;
pierre.rousselot_ird.fr
committed
title('Amplitude of the bins'); colorbar;
ylabel('Bins');
pierre.rousselot_ird.fr
committed
saveas(gcf,[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Amplitude_bins'],'fig')
%% Add external pressure sensor choice
pierre.rousselot_ird.fr
committed
disp('****')
pressure_data = input('Do you have external pressure data? 1/0 [0]');
if isempty(pressure_data)
pressure_data = 0;
if pressure_data == 1
[filename, path, ~] = uigetfile('*.cnv', 'Select external pressure file');
press_file = fullfile(path, filename);
press_file = read_sbe(press_file);
real_depth = press_file.depSM(first:last);
dpt = sw_dpth(press,mooring.lat)'; % convert pressure to depth, press needs to have dimension (n x 1)
plot(real_depth, 'b')
grid on
plot(dpt, 'r')
legend('External sensor depth','ADCP sensor depth');
xlabel('Time Index')
exsens_data = input('Do you want to use external sensor data? 1/0 [0]');
if isempty(exsens_data)
exsens_data = 1;
if exsens_data
dpt = real_depth-offset_pres_sens;
dpt1 = repmat(dpt,nbin,1);
binmat = repmat((1:nbin)',1,length(dpt1));
z(1,:) = dpt1(1,:)-fbind;
for ii = 2:length(binmat(:,1))
z(ii,:) = z(1,:)-(binmat(ii,:)-1.5)*blen;
end
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
end
end
%% If upward looking: determine range of surface bins used for instrument depth correction below!
surface_bins = 0;
if pressure_data == 0
disp('****')
surface_bins = input('Do you want to apply an offset using surface reflection? 1/0 [1]');
if isempty(surface_bins)
surface_bins = 1;
end
if surface_bins == 1
sbins = input('Determine range of surface bins used for instrument depth correction (with aplitude plot, ie. 30:35): ');
%% Calculate depth of each bin:
dpt = sw_dpth(press,mooring.lat)'; % convert pressure to depth, press needs to have dimension (n x 1)
dpt1 = repmat(dpt,nbin,1);
binmat = repmat((1:nbin)',1,length(dpt1));
% If ADCP is upward-looking a depth correction can be inferred from the surface reflection, which is done in adcp_surface_fit
disp('****')
if strcmp(adcp.direction,'up')
[z,dpt1,offset,xnull] = adcp_surface_fit(dpt,ea,sbins,blen,tlen,lag,blnk,nbin,fbind);
elseif strcmp(adcp.direction,'dn')
z(1,:) = dpt1(1,:)+(tlen+blen+lag)*0.5+blnk;
for ii = 2:length(binmat(:,1))
z(ii,:) = z(1,:)+(binmat(ii,:)-1.5)*blen;
end
else
error('Bin depth calculation: unknown direction!');
saveas(figure(5),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Histdiff_depth'],'fig')
saveas(figure(6),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Offset_depth'],'fig')
saveas(figure(7),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Amplitude_bins_2'],'fig')
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
dpt = sw_dpth(press,mooring.lat)'; % convert pressure to depth, press needs to have dimension (n x 1)
dpt1 = repmat(dpt,nbin,1);
z = dpt;
figure
plot(dpt)
hold on
plot(adcp.instr_depth*ones(length(dpt),1),'--r')
hold off
grid on
axis([0 length(dpt) min(min(dpt,adcp.instr_depth*ones(1,length(dpt))))-50 max(max(dpt,adcp.instr_depth*ones(1,length(dpt))))+50])
xlabel('Ensembles')
ylabel('Depth [m]')
legend('Depth','Nominal Depth')
offset = input('Do you want to apply a manual offset? 1/0 [0]');
if isempty(offset)
offset = 0;
elseif offset == 1
offset = input('->Enter new offset:');
dpt = dpt + offset;
dpt1 = dpt1 + offset;
z = z + offset;
hold on
plot(dpt,'g')
hold off
legend('Depth', 'Nominal Depth', 'Corrected depth')
elseif offset == 0
offset = 0;
end
binmat = repmat((1:nbin)',1,length(dpt1));
if strcmp(adcp.direction,'up')
z(1,:) = dpt1(1,:)-(tlen+blen+lag)*0.5-blnk;
for ii = 2:length(binmat(:,1))
z(ii,:) = z(1,:)-(binmat(ii,:)-1.5)*blen;
end
elseif strcmp(adcp.direction,'dn')
z(1,:) = dpt1(1,:)+(tlen+blen+lag)*0.5+blnk;
for ii = 2:length(binmat(:,1))
z(ii,:) = z(1,:)+(binmat(ii,:)-1.5)*blen;
end
else
error('Bin depth calculation: unknown direction!');
end
u1=u; v1=v; w1=w; vel_err1=vel_err; ea1=ea; corr1=corr; z1=z;
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
if surface_bins == 1
disp('****')
disp('Remove bad data near surface due to sidelobe');
if strcmp(adcp.direction,'up')
for i = 1:length(time)
sz_dpt(i) = adcp_shadowzone(dpt1(i),raw.config.sysconfig.angle); % depending on the instrument depth and the beam angle the shadow zone, i.e. the depth below the surface which is contaminated by the surface reflection is determined
iz(i) = find(z(:,i)>sz_dpt(i),1,'last');
sbin(i) = bin(iz(i));
%sbin(i)=30; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% here a manual criterion should be hard-coded if
% adcp_check_surface (below) shows bad velocities close to the
% surface
u1(sbin(i)+1:end,i) = NaN;
v1(sbin(i)+1:end,i) = NaN;
w1(sbin(i)+1:end,i) = NaN;
vel_err1(sbin(i)+1:end,i) = NaN;
ea1(sbin(i)+1:end,i) = NaN;
corr1(sbin(i)+1:end,i) = NaN;
z1(sbin(i)+1:end,i) = NaN;
end
if 1
bins = nmedian(sbin)-4:nmedian(sbin)+4;
adcp_check_surface(bins,u,u1,v,v1,corr,corr1,time,bin,z,z1,sz_dpt);
% here the closest bins below the surface are plotted that are supposed to have good velocities, if there are still bad velocities a manual criterion needs to be found
reply_sidelobe = input('Do you want to apply manual criterion? 1/0 [0]:');
if isempty(reply_sidelobe)
reply_sidelobe = 0;
end
if reply_sidelobe
bin_cutoff = input('->Enter new bin cutoff value:');
u1=u; v1=v; w1=w; vel_err1=vel_err; ea1=ea; corr1=corr; z1=z;
u1(bin_cutoff+1:end,:) = NaN;
v1(bin_cutoff+1:end,:) = NaN;
w1(bin_cutoff+1:end,:) = NaN;
vel_err1(bin_cutoff+1:end,:) = NaN;
ea1(bin_cutoff+1:end,:) = NaN;
corr1(bin_cutoff+1:end,:) = NaN;
z1(bin_cutoff+1:end,:) = NaN;
adcp_check_surface(bins,u,u1,v,v1,corr,corr1,time,bin,z,z1,sz_dpt);
end
end
saveas(figure(7),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','Meridional_zonal_velocity'],'fig')
if pressure_data == 0
adcp.z_offset = offset;
end
data.u = u1;
data.v = v1;
data.w = w1;
data.e = vel_err1;
data.ea = ea1;
data.pg = pg;
data.time = time;
data.z_bins = z1;
pierre.rousselot_ird.fr
committed
data.depth = dpt1;
data.temp = temp;
data.sspd = soundspeed;
data.lat = mooring.lat;
data.lon = mooring.lon;
disp('****')
reply_ts = input('Do you want to calculate Target Strength? (CTD profile needed) 1/0 [0]:');
if isempty(reply_ts)
reply_ts = 0;
end
if reply_ts == 1
% EA0: noisefloor from ADCP electronic noise; a scalar
disp('->Calculate Target Strength')
disp(['->EA0=' num2str(EA0) '(noisefloor from ADCP electronic noise) ; pH=8.1; (seawater pH)'])
file_ctd = input('->Path to CTD .nc file [''*.nc'']:');
if exist(file_ctd, 'file')
ctd.depth = ncread(file_ctd, 'PRES');
ctd.temp = ncread(file_ctd, 'TEMP');
ctd.sal = ncread(file_ctd, 'PSAL');
ctd.tempR = ones(size(z));
ctd.salR = ones(size(z));
% for each ADCP bin, find salinity and temperature values
for i_bin = 1:nbin
[~,ind_depth] = min(abs(z(i_bin,:)-ctd.depth(:)));
ctd.salR(i_bin,:) = ctd.sal(ind_depth);
ctd.tempR(i_bin,:) = ctd.temp(ind_depth);
end
% Compute TS
[sspd,~] = meshgrid(soundspeed,1:nbin);
data.ts = target_strength(data.ea,EA0,ctd.salR,ctd.tempR,sspd,temp',soundspeed',z,adcp.config.cell,adcp.config.blank,adcp.config.sysconfig.angle,freq);
figure;
colormap jet;
pcolor(data.ts.ts);
shading flat;
title('Target Strength'); colorbar;
ylabel('Bins');
xlabel('Time index');
saveas(gcf,[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','target_strength'],'fig')
else
disp('CTD file doesn''t exist')
reply_ts = 0;
end
end
%% Save raw data
disp('****')
disp('Saving raw data')
save([fpath_output, mooring.name '_' num2str(adcp.sn) '_instr_' sprintf('%02d',instr) '.mat'],'adcp','mooring','data','-v7.3');
if strcmp(adcp.direction,'up')
Z = fliplr(blen/2:blen:max(z(:))+blen);
Zmax = max(Z);
u_interp = NaN(length(time),length(Z));
v_interp = NaN(length(time),length(Z));
if reply_ts == 1
ts_interp = NaN(length(time),length(Z));
end
for i=1:length(time)
% indice correspondant sur la grille finale Z
ind = round((Zmax-z(1,i))/blen)+1;
% filling the grid
npts = min([length(Z)+ind+1 length(bin)]);
u_interp(i,ind:ind+npts-1) = u1(1:npts,i);
v_interp(i,ind:ind+npts-1) = v1(1:npts,i);
if reply_ts == 1
ts_interp(i,ind:ind+npts-1) = data.ts.ts(1:npts,i);
end
end
else
Z = (round(min(data.depth))+blen/2:blen:max(z(:))+blen);
Zmax = min(Z);
u_interp = NaN(length(time),length(Z));
v_interp = NaN(length(time),length(Z));
if reply_ts == 1
ts_interp = NaN(length(time),length(Z));
end
for i=1:length(time)
% indice correspondant sur la grille finale Z
ind = -(round((Zmax-z(1,i))/blen)-1);
% filling the grid
npts = min([length(Z)-ind-1 length(bin)]);
u_interp(i,ind:ind+npts-1) = u1(1:npts,i);
v_interp(i,ind:ind+npts-1) = v1(1:npts,i);
if reply_ts == 1
ts_interp(i,ind:ind+npts-1) = data.ts.ts(1:npts,i);
end
end
end
if reply_ts == 0
ts_interp = NaN;
disp('****')
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
horz_interp = input('Do you want to do an horizontal interpolation? 1/0 [1]:');
if isempty(horz_interp)
horz_interp = 1;
end
if horz_interp
filt = input('Do you want to filter tide? 1/0 [1]:');
if isempty(filt)
filt = 1;
end
sub = input('Do you want to subsampling data daily? 1/0 [1]:');
if isempty(filt)
sub = 1;
end
if filt
if sub
disp('Horizontal interpolation, Filtering & Subsampling')
else
disp('Horizontal interpolation & Filtering')
end
disp('Filtering using 40h Hamming filter (other disponible filter: FFT filter, tide model)')
else
if sub
disp('Horizontal interpolation & Subsampling')
else
disp('Horizontal interpolation')
end
end
else
filt = 0;
sub = 0;
end
[uintfilt,vintfilt,tsintfilt,inttim,utid_baro,vtid_baro] = adcp_filt_sub(data,u_interp',v_interp',ts_interp',1:length(Z),reply_ts, filt, sub);
saveas(figure(5),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','data_raw_filt_subsampled_1'],'fig')
saveas(figure(6),[fpath_output,mooring.name,'_',num2str(adcp.sn),'_instr_',num2str(instr),'_','data_raw_filt_subsampled_2'],'fig')
if horz_interp == 0
uintfilt = u_interp';
vintfilt = v_interp';
inttim = data.time;
end
disp('****')
bin_start = input('Determine first bin indice with good interpolated data: '); % bin indice where good interpolated data for the whole dataset start
disp('****')
bin_end = input('Determine last bin indice with good interpolated data: '); % bin indice where good interpolated data for the whole dataset start
data.uintfilt = uintfilt(bin_start:bin_end,:);
data.vintfilt = vintfilt(bin_start:bin_end,:);
if reply_ts == 1
data.tsintfilt = tsintfilt(bin_start:bin_end,:);
end
data.Z = Z(bin_start:bin_end);
data.inttim = inttim;
save([fpath_output, mooring.name '_' num2str(adcp.sn) '_instr_' sprintf('%02d',instr) '_int_filt_sub.mat'],'adcp','mooring','data');
close all
[C,h] = contourf(inttim,Z(bin_start:bin_end),uintfilt(bin_start:bin_end,:),niv_u);
set(h,'LineColor','none');
caxis(niv_u([1 end]));
h=colorbar;
ylabel(h,'U [m s^-^1]');
set(gca,'ydir', 'reverse');
ylabel('Depth (m)');
if filt
title({[mooring.name, ' - ZONAL VELOCITY - RDI ',num2str(freq),' kHz (filtered from tide)']});
else
title({[mooring.name, ' - ZONAL VELOCITY - RDI ',num2str(freq),' kHz']});
end
[C,h] = contourf(inttim,Z(bin_start:bin_end),vintfilt(bin_start:bin_end,:),niv_v);
set(h,'LineColor','none');
ylabel(h,'V [m s^-^1]');
set(gca,'ydir', 'reverse');
ylabel('Depth (m)');
if filt
title({[mooring.name, ' - MERIDIONAL VELOCITY - RDI ',num2str(freq),' kHz (filtered from tide)']});
else
title({[mooring.name, ' - MERIDIONAL VELOCITY - RDI ',num2str(freq),' kHz']});
end
graph_name = [fpath_output, mooring.name '_', num2str(instr), '_U_V_int_filt_sub'];
set(hf,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)]);
if reply_ts == 1
hf=figure('position', [0, 0, 1400, 500]);
colormap jet
pcolor(inttim,Z(bin_start:bin_end),tsintfilt(bin_start:bin_end,:)); shading interp;
h = colorbar;
ylabel(h,'Target Strength [dB1]');
set(gca,'ydir', 'reverse');
ylabel('Depth (m)');
ylim([0,round(max(max(z,adcp.instr_depth)))]);
%change figure label in HH:MM
gregtick;
if filt
title({[mooring.name, ' - TARGET STRENGTH - RDI ',num2str(freq),' kHz (filtered from tide)']});
else
title({[mooring.name, ' - TARGET STRENGTH - RDI ',num2str(freq),' kHz']});
end
graph_name = [fpath_output, mooring.name, '_', num2str(instr), '_TS_int_filt_sub'];
set(hf,'Units','Inches');
pos = get(hf,'Position');
set(hf,'PaperPositionMode','Auto','PaperUnits','Inches','PaperSize',[pos(3), pos(4)]);
saveas(hf,graph_name,'jpeg');
end
% %% Plot tide
% hf=figure('position', [0, 0, 1400, 1000]);
% niv_tide = (-5:0.2:5);
% utid_baro = utid_baro * 100;
% vtid_baro = vtid_baro * 100;
% %u
% subplot(2,1,1);
% colormap jet
% [C,h] = contourf(inttim,Z(bin_start:bin_end),utid_baro(bin_start:bin_end,:),niv_tide);
% set(h,'LineColor','none');
% caxis(niv_tide([1 end]));
% h=colorbar;
% ylabel(h,'U [cm s^-^1]');
% set(gca,'ydir', 'reverse');
% ylabel('Depth (m)');
% ylim([0,adcp.instr_depth]);
% %change figure label in HH:MM
% gregtick;
% title({[mooring.name, ' - ZONAL TIDE VELOCITY - RDI ',num2str(freq),' kHz']});
% [C,h] = contourf(inttim,Z(bin_start:bin_end),vtid_baro(bin_start:bin_end,:),niv_tide);
% set(h,'LineColor','none');
% caxis(niv_tide([1 end]));
% h = colorbar;
% ylabel(h,'V [cm s^-^1]');
% set(gca,'ydir', 'reverse');
% ylabel('Depth (m)');
% ylim([0,adcp.instr_depth]);
% %change figure label in HH:MM
% gregtick;
% title({[mooring.name, ' - MERIDIONAL TIDE VELOCITY - RDI ',num2str(freq),' kHz']});
disp('****')
disp('Creating .nc file')
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
% Input parameters for NETCDF Global Attributes
tc_globAttFilename = fullfile('tools/input_GlobalAttrParameters.xls'); % JLL 2020/12 Il serait judicieux de remonter cette valeur en début du script template_get_adcp_data.m
%% Prepare informations and variables required to create NETCDF file %%
[yr_start , ~, ~] = gregorian(data.inttim(1));
[yr_end, ~, ~] = gregorian(data.inttim(length(data.inttim)));
% Read inputs metadata required for NETCDF Global Attributes
[~,~,cell_ncAttributes] = xlsread(tc_globAttFilename);
% Complete output path and filename
tc_ncFilenam_out = fullfile(fpath_output,['ADCP_',mooring.name,'_',num2str(yr_start),'_',num2str(yr_end),'_1d.nc']);
% Assign a "4D-size" (TIME,DEPTH,LATITUDE,LONGITUDE) to the ADCP current variables : UINTFILT, VINTFILT
td_uADCP = ones(numel(data.inttim),numel(data.Z),numel(data.lat),numel(data.lon)) * d_fillvalue;
td_uADCP(:,:,1,1) = data.uintfilt';
td_vADCP = ones(numel(data.inttim),numel(data.Z),numel(data.lat),numel(data.lon)) * d_fillvalue;
td_vADCP(:,:,1,1) = data.vintfilt';
% Flip for convention
data.Z = fliplr(data.Z);
td_uADCP = fliplr(td_uADCP);
td_vADCP = fliplr(td_vADCP);
% Group general ADCP mooring informations and ADCP data to be written in NETCDF file format
struct_dataADCP = struct('mooringName', mooring.name, 'mooringLat', mooring.lat,...
'mooringLon', mooring.lon, 'time', data.inttim, 'depth', data.Z,...
'u', td_uADCP, 'v', td_vADCP);
%% Write netcdf file %%
f_w_ADCP_ncOS(tc_ncFilenam_out,cell_ncAttributes,struct_dataADCP,d_fillvalue);
disp('****')