Skip to content

Commit 3f1e68c

Browse files
cabanesccabanesc
authored
Merge pull request #16 from delphinedobler/theta_levels-record-for-the-plot-function
Record and Use Theta levels for the plot function and a few enhancements
2 parents 5a48c93 + c28b279 commit 3f1e68c

File tree

9 files changed

+884
-435
lines changed

9 files changed

+884
-435
lines changed

‎matlab_codes/calculate_piecewisefit.m‎

Lines changed: 46 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -8,6 +8,8 @@ function calculate_piecewisefit( pn_float_dir, pn_float_name, po_system_configur
88

99
%
1010
% Cecile Cabanes, June 2013 : calculate off-diagonal terms for error estimate: add horizontal covariance to track changes: see "change config 129"
11+
% Delphine Dobler (DD), September 2024 :
12+
% 4 - save computed theta levels in a file
1113

1214
%pn_float_dir='testfloats/';
1315
%pn_float_name='robbins4900178';
@@ -156,10 +158,15 @@ function calculate_piecewisefit( pn_float_dir, pn_float_name, po_system_configur
156158
return
157159
end
158160

161+
% DD (2024/09-4): initialisation of arrays to be saved in a file
162+
Theta_nseq=NaN*zeros(10,n_seq);
163+
Index=NaN*zeros(10,n);
164+
Plevels=NaN*zeros(10,n_seq);
165+
159166
% loop through sequences of calseries, ie if time series is broken into segments --
160-
for i=1:n_seq
167+
for iseq=1:n_seq
161168

162-
calindex = find(calseries==unique_cal(i));
169+
calindex = find(calseries==unique_cal(iseq));
163170
k = length(calindex);
164171

165172
% choose 10 float theta levels to use in the piecewise linear fit --------
@@ -180,6 +187,31 @@ function calculate_piecewisefit( pn_float_dir, pn_float_name, po_system_configur
180187

181188
[Theta, P, index, var_s_th, th] =...
182189
find_10thetas( unique_SAL, unique_PTMP, unique_PRES, unique_la_ptmp, use_theta_gt, use_theta_lt, use_pres_gt, use_pres_lt, use_percent_gt);
190+
191+
% DD (2024/09-4): Save Theta-related information for plotting functions
192+
Index(:,calindex)=index;
193+
Plevels(:,iseq)=P;
194+
[mvth,~]=size(var_s_th);
195+
Theta_nseq(:,iseq)=Theta;
196+
if iseq==1
197+
Var_s_Thetas=NaN.*ones(mvth,n_seq);
198+
Thetas=NaN.*ones(mvth,n_seq);
199+
Var_s_Thetas(:,iseq)=var_s_th;
200+
Thetas(:,iseq)=th;
201+
else
202+
[mVth,~]=size(Var_s_Thetas);
203+
tmp1=Var_s_Thetas;
204+
tmp2=Thetas;
205+
mm=max(mVth,mvth);
206+
Var_s_Thetas=NaN.*ones(mm,n_seq);
207+
Thetas=NaN.*ones(mm,n_seq);
208+
Var_s_Thetas(1:mVth,:)=tmp1;
209+
Var_s_Thetas(1:mvth,iseq)=var_s_th;
210+
Thetas(1:mVth,:)=tmp2;
211+
Thetas(1:mvth,iseq)=th;
212+
end
213+
%
214+
183215

184216
pp=find(isnan(index)==0);
185217
if(isempty(pp)==0) % only proceed when there are valid levels ----
@@ -229,18 +261,18 @@ function calculate_piecewisefit( pn_float_dir, pn_float_name, po_system_configur
229261
if isempty(breaks)
230262
[xfit(calindex), pcond_factor(calindex), pcond_factor_err(calindex), time_deriv(calindex), ...
231263
time_deriv_err(calindex), sta_mean(calindex), sta_rms(calindex), NDF, fitcoef, fitbreaks] = ...
232-
fit_cond(x, y, err, covariance, 'max_no_breaks', max_breaks(i));
264+
fit_cond(x, y, err, covariance, 'max_no_breaks', max_breaks(iseq));
233265
else
234-
breaks_in = breaks(i,:);
266+
breaks_in = breaks(iseq,:);
235267
breaks_in = breaks_in(find(isfinite(breaks_in)));
236-
if isempty(max_breaks(i))
268+
if isempty(max_breaks(iseq))
237269
[xfit(calindex), pcond_factor(calindex), pcond_factor_err(calindex), time_deriv(calindex), ...
238270
time_deriv_err(calindex), sta_mean(calindex), sta_rms(calindex), NDF, fitcoef, fitbreaks] = ...
239271
fit_cond(x, y, err, covariance, 'breaks', breaks_in);
240272
else
241273
[xfit(calindex), pcond_factor(calindex), pcond_factor_err(calindex), time_deriv(calindex), ...
242274
time_deriv_err(calindex), sta_mean(calindex), sta_rms(calindex), NDF, fitcoef, fitbreaks] = ...
243-
fit_cond(x, y, err, covariance, 'breaks', breaks_in, 'max_no_breaks', max_breaks(i));
275+
fit_cond(x, y, err, covariance, 'breaks', breaks_in, 'max_no_breaks', max_breaks(iseq));
244276
end
245277
end
246278

@@ -266,16 +298,22 @@ function calculate_piecewisefit( pn_float_dir, pn_float_name, po_system_configur
266298
sta_SAL1(:,calindex) = sw_salt( (sta_COND(:,calindex)+sta_COND_err(:,calindex))/sw_c3515, unique_PTMP, 0 );
267299
sta_SAL_err(:,calindex) = abs(sta_SAL(:,calindex)-sta_SAL1(:,calindex));
268300

269-
fcoef(i,1:length(fitcoef)) = fitcoef;
301+
fcoef(iseq,1:length(fitcoef)) = fitcoef;
270302
if ~isempty(fitbreaks)
271-
fbreaks(i,1:length(fitbreaks)) = fitbreaks;
303+
fbreaks(iseq,1:length(fitbreaks)) = fitbreaks;
272304
end
273305

274306
end
275307

276308
end %if there are valid levels
277309
end %for each unique_cal
278310

311+
% DD (2024/09-4) : save theta levels information inside a file.
312+
Theta=Theta_nseq;
313+
ls_theta = strcat( po_system_configuration.FLOAT_CALIB_DIRECTORY, ...
314+
pn_float_dir, 'selected_theta_', pn_float_name,'.mat');
315+
save(ls_theta, 'Theta','Index','Plevels','Var_s_Thetas','Thetas')
316+
% end of DD (2024/09-4)
279317

280318
% save calibration data --------------------------------
281319

‎matlab_codes/create_la_wmo_boxes_file.m‎

Lines changed: 111 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,111 @@
1+
function create_la_wmo_boxes_file( po_system_configuration )
2+
3+
% function create_la_wmo_boxes_file( po_system_configuration )
4+
%
5+
% Delphine Dobler, August 2024
6+
7+
clim_ctd_=strcat(po_system_configuration.HISTORICAL_DIRECTORY,po_system_configuration.HISTORICAL_CTD_PREFIX);
8+
clim_argo_=strcat(po_system_configuration.HISTORICAL_DIRECTORY,po_system_configuration.HISTORICAL_ARGO_PREFIX);
9+
clim_bottle_=strcat(po_system_configuration.HISTORICAL_DIRECTORY,po_system_configuration.HISTORICAL_BOTTLE_PREFIX);
10+
la_wmo_boxes_file=strcat( po_system_configuration.CONFIG_DIRECTORY, po_system_configuration.CONFIG_WMO_BOXES);
11+
12+
la_wmo_boxes=zeros(648,4);
13+
la_wmo_boxes(:,1)=[1800;1700;1600;1500;1400;1300;1200;1100;1000;3000;3100;...
14+
3200;3300;3400;3500;3600;3700;3800;1801;1701;1601;1501;...
15+
1401;1301;1201;1101;1001;3001;3101;3201;3301;3401;3501;...
16+
3601;3701;3801;1802;1702;1602;1502;1402;1302;1202;1102;...
17+
1002;3002;3102;3202;3302;3402;3502;3602;3702;3802;1803;...
18+
1703;1603;1503;1403;1303;1203;1103;1003;3003;3103;3203;...
19+
3303;3403;3503;3603;3703;3803;1804;1704;1604;1504;1404;...
20+
1304;1204;1104;1004;3004;3104;3204;3304;3404;3504;3604;...
21+
3704;3804;1805;1705;1605;1505;1405;1305;1205;1105;1005;...
22+
3005;3105;3205;3305;3405;3505;3605;3705;3805;1806;1706;...
23+
1606;1506;1406;1306;1206;1106;1006;3006;3106;3206;3306;...
24+
3406;3506;3606;3706;3806;1807;1707;1607;1507;1407;1307;...
25+
1207;1107;1007;3007;3107;3207;3307;3407;3507;3607;3707;...
26+
3807;1808;1708;1608;1508;1408;1308;1208;1108;1008;3008;...
27+
3108;3208;3308;3408;3508;3608;3708;3808;1809;1709;1609;...
28+
1509;1409;1309;1209;1109;1009;3009;3109;3209;3309;3409;...
29+
3509;3609;3709;3809;1810;1710;1610;1510;1410;1310;1210;...
30+
1110;1010;3010;3110;3210;3310;3410;3510;3610;3710;3810;...
31+
1811;1711;1611;1511;1411;1311;1211;1111;1011;3011;3111;...
32+
3211;3311;3411;3511;3611;3711;3811;1812;1712;1612;1512;...
33+
1412;1312;1212;1112;1012;3012;3112;3212;3312;3412;3512;...
34+
3612;3712;3812;1813;1713;1613;1513;1413;1313;1213;1113;...
35+
1013;3013;3113;3213;3313;3413;3513;3613;3713;3813;1814;...
36+
1714;1614;1514;1414;1314;1214;1114;1014;3014;3114;3214;...
37+
3314;3414;3514;3614;3714;3814;1815;1715;1615;1515;1415;...
38+
1315;1215;1115;1015;3015;3115;3215;3315;3415;3515;3615;...
39+
3715;3815;1816;1716;1616;1516;1416;1316;1216;1116;1016;...
40+
3016;3116;3216;3316;3416;3516;3616;3716;3816;1817;1717;...
41+
1617;1517;1417;1317;1217;1117;1017;3017;3117;3217;3317;...
42+
3417;3517;3617;3717;3817;7817;7717;7617;7517;7417;7317;...
43+
7217;7117;7017;5017;5117;5217;5317;5417;5517;5617;5717;...
44+
5817;7816;7716;7616;7516;7416;7316;7216;7116;7016;5016;...
45+
5116;5216;5316;5416;5516;5616;5716;5816;7815;7715;7615;...
46+
7515;7415;7315;7215;7115;7015;5015;5115;5215;5315;5415;...
47+
5515;5615;5715;5815;7814;7714;7614;7514;7414;7314;7214;...
48+
7114;7014;5014;5114;5214;5314;5414;5514;5614;5714;5814;...
49+
7813;7713;7613;7513;7413;7313;7213;7113;7013;5013;5113;...
50+
5213;5313;5413;5513;5613;5713;5813;7812;7712;7612;7512;...
51+
7412;7312;7212;7112;7012;5012;5112;5212;5312;5412;5512;...
52+
5612;5712;5812;7811;7711;7611;7511;7411;7311;7211;7111;...
53+
7011;5011;5111;5211;5311;5411;5511;5611;5711;5811;7810;...
54+
7710;7610;7510;7410;7310;7210;7110;7010;5010;5110;5210;...
55+
5310;5410;5510;5610;5710;5810;7809;7709;7609;7509;7409;...
56+
7309;7209;7109;7009;5009;5109;5209;5309;5409;5509;5609;...
57+
5709;5809;7808;7708;7608;7508;7408;7308;7208;7108;7008;...
58+
5008;5108;5208;5308;5408;5508;5608;5708;5808;7807;7707;...
59+
7607;7507;7407;7307;7207;7107;7007;5007;5107;5207;5307;...
60+
5407;5507;5607;5707;5807;7806;7706;7606;7506;7406;7306;...
61+
7206;7106;7006;5006;5106;5206;5306;5406;5506;5606;5706;...
62+
5806;7805;7705;7605;7505;7405;7305;7205;7105;7005;5005;...
63+
5105;5205;5305;5405;5505;5605;5705;5805;7804;7704;7604;...
64+
7504;7404;7304;7204;7104;7004;5004;5104;5204;5304;5404;...
65+
5504;5604;5704;5804;7803;7703;7603;7503;7403;7303;7203;...
66+
7103;7003;5003;5103;5203;5303;5403;5503;5603;5703;5803;...
67+
7802;7702;7602;7502;7402;7302;7202;7102;7002;5002;5102;...
68+
5202;5302;5402;5502;5602;5702;5802;7801;7701;7601;7501;...
69+
7401;7301;7201;7101;7001;5001;5101;5201;5301;5401;5501;...
70+
5601;5701;5801;7800;7700;7600;7500;7400;7300;7200;7100;...
71+
7000;5000;5100;5200;5300;5400;5500;5600;5700;5800];
72+
73+
for i_hist_type=2:4
74+
if i_hist_type==2
75+
clim_dir=clim_ctd_;
76+
clim_type="ctd_";
77+
prefix=po_system_configuration.HISTORICAL_CTD_PREFIX;
78+
end
79+
if i_hist_type==4
80+
clim_dir=clim_argo_;
81+
clim_type="argo_";
82+
prefix=po_system_configuration.HISTORICAL_ARGO_PREFIX;
83+
end
84+
if i_hist_type==3
85+
clim_dir=clim_bottle_;
86+
clim_type="bot_";
87+
prefix=po_system_configuration.HISTORICAL_BOTTLE_PREFIX;
88+
end
89+
90+
prefix_str=split(prefix,"/");
91+
clim_dir=strcat(po_system_configuration.HISTORICAL_DIRECTORY,strjoin(prefix_str(1:end-1),"/"));
92+
93+
listing=dir(clim_dir);
94+
n_files=size(listing,1);
95+
for ifile=1:n_files
96+
l_name=listing(ifile).name;
97+
if strlength(l_name) < strlength(clim_type)
98+
continue
99+
end
100+
if l_name(1:strlength(clim_type))== clim_type
101+
box_no=strrep(l_name,clim_type,'');
102+
box_no=strrep(box_no,'.mat','');
103+
box_no=str2double(box_no);
104+
i_box=(la_wmo_boxes(:,1)==box_no);
105+
la_wmo_boxes(i_box,i_hist_type)=1;
106+
end
107+
108+
end
109+
end
110+
111+
save(la_wmo_boxes_file,'la_wmo_boxes')

‎matlab_codes/get_region_ow.m‎

Lines changed: 18 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11

2-
function [ pa_grid_lat, pa_grid_long, pa_grid_dates ] = get_region_ow(pa_wmo_numbers, pa_float_name, po_config_data) ;
2+
function [ pa_grid_lat, pa_grid_long, pa_grid_dates, pa_grid_box_id ] = get_region_ow(pa_wmo_numbers, pa_float_name, po_config_data)
33

44
% function [ pa_grid_lat, pa_grid_long, pa_grid_dates ] = get_region_ow( pa_wmo_numbers, pa_float_name, po_config_data ) ;
55
%
@@ -18,12 +18,16 @@
1818
% Annie Wong, December 2007
1919
% Breck Owens, December 2006
2020
% C.Cabanes (2015) to save time: load only long, lat, dates and source data
21+
% D.Dobler (DD), August 2024: to save more time later in retr_region_ow: also output
22+
% corresponding boxes ids. Additional minor modifications: comment date_hist as it is unused, and
23+
% correct indentation twice.
2124

22-
pa_grid_lat = [ ] ;
23-
pa_grid_long = [ ] ;
24-
pa_grid_dates = [ ] ;
25+
pa_grid_lat = [ ] ;
26+
pa_grid_long = [ ] ;
27+
pa_grid_dates = [ ] ;
28+
pa_grid_box_id = [ ] ;
2529

26-
[m,n]=size(pa_wmo_numbers);
30+
[m,~]=size(pa_wmo_numbers);
2731

2832

2933

@@ -40,7 +44,7 @@
4044
%toc
4145

4246
not_use=[];
43-
date_hist = changedates(lo_box_data.dates);
47+
%date_hist = changedates(lo_box_data.dates);
4448
lo_box_data.lat(not_use)=[];
4549
lo_box_data.long(not_use)=[];
4650
lo_box_data.dates(not_use)=[];
@@ -55,22 +59,24 @@
5559
not_use=[];
5660
for i=1:length(lo_box_data.lat)
5761
profile=lo_box_data.source{i};
58-
jj=findstr(profile,'_');
62+
jj=findstr(profile,'_');
5963
ref_float=profile(1:jj-1);
6064
kk=findstr(pa_float_name, ref_float);
6165
if(isempty(kk)==0)
62-
not_use=[not_use,i];
66+
not_use=[not_use,i];
6367
end
6468
end
6569
lo_box_data.lat(not_use)=[];
6670
lo_box_data.long(not_use)=[];
6771
lo_box_data.dates(not_use)=[];
6872
%-----------------------------------------------------------------------
6973
end
70-
71-
pa_grid_lat = [ pa_grid_lat, lo_box_data.lat ] ;
72-
pa_grid_long = [ pa_grid_long, lo_box_data.long ] ;
73-
pa_grid_dates = [ pa_grid_dates, lo_box_data.dates ] ;
74+
75+
n_observations = size(lo_box_data.lat,2);
76+
pa_grid_box_id = [ pa_grid_box_id, pa_wmo_numbers(ln_index,1)*ones(1,n_observations)] ;
77+
pa_grid_lat = [ pa_grid_lat, lo_box_data.lat ] ;
78+
pa_grid_long = [ pa_grid_long, lo_box_data.long ] ;
79+
pa_grid_dates = [ pa_grid_dates, lo_box_data.dates ] ;
7480
end
7581
end
7682
fclose('all');

‎matlab_codes/map_data_grid.m‎

Lines changed: 9 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -39,6 +39,10 @@
3939
% Bretherton, etal, 1975) is removed. The error estimate includes
4040
% the contributions from both the mapping and the mean
4141
% estimatation.
42+
%
43+
% Delphine Dobler (DD), August 2024:
44+
% 3.3 - Performance: Inside map_data_grid, comment calculation
45+
% of vdataerror as this is never used afterwards.
4246
% **********************************************************************
4347

4448
[m,n]=size(v);
@@ -72,8 +76,11 @@
7276
Cmd = s*covarxyt_pv(posdata, posdata, gs, ts, q, phi, pv);
7377
vdata = Cmd*wght + vm;
7478
% include error in mean (eqn 24 of Bretherton, etal)
75-
vdataerror = repmat(sqrt(s-diag(Cmd*Cdd*Cmd') ...
76-
+(1- sum(Cmd*Cdd,2)).^2/sum_Cdd),1,1);
79+
%vdataerror = repmat(sqrt(s-diag(Cmd*Cdd*Cmd') ...
80+
% +(1- sum(Cmd*Cdd,2)).^2/sum_Cdd),1,1);
81+
% DD (2024/08-3.3): never used afterwards, save some time (Uncomment to
82+
% compute)
83+
vdataerror = [];
7784

7885
% map to posgrid -----
7986

0 commit comments

Comments
 (0)