From ef7b57f80ef2aec6289a617b9e32665266a57858 Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Thu, 3 Feb 2022 13:32:59 -0700 Subject: [PATCH 1/8] #1824 Added MLCAPE --- met/data/config/PB2NCConfig_default | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/met/data/config/PB2NCConfig_default b/met/data/config/PB2NCConfig_default index 82750e9912..615bea33af 100644 --- a/met/data/config/PB2NCConfig_default +++ b/met/data/config/PB2NCConfig_default @@ -121,7 +121,8 @@ obs_prefbufr_map = [ { key = "D_MIXR"; val = "MIXR"; }, { key = "D_PRMSL"; val = "PRMSL"; }, { key = "D_PBL"; val = "PBL"; }, - { key = "D_CAPE"; val = "CAPE"; } + { key = "D_CAPE"; val = "CAPE"; }, + { key = "D_MLCAPE"; val = "MLCAPE"; } ]; //////////////////////////////////////////////////////////////////////////////// From a6ba1abf807d0ae193aea1fc2e6e3771979c034a Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Thu, 3 Feb 2022 13:33:24 -0700 Subject: [PATCH 2/8] #1824 Cleanup --- met/src/tools/other/pb2nc/calcape.f | 28 ---------------------------- 1 file changed, 28 deletions(-) diff --git a/met/src/tools/other/pb2nc/calcape.f b/met/src/tools/other/pb2nc/calcape.f index e8d3a60b9d..d833246671 100644 --- a/met/src/tools/other/pb2nc/calcape.f +++ b/met/src/tools/other/pb2nc/calcape.f @@ -130,8 +130,6 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, CINS(I,J) = 0.0 THESP(I,J)= 0.0 IEQL(I,J) = LM+1 -C print*,' DEBUG HS calcape IEQL(I,J): by LM+1: ',IEQL(I,J) -C print*,' DEBUG HS calcape LCL(I,J) to 0 ' LCL(I,J)=0 ENDDO @@ -233,11 +231,8 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, APESPK=(H10E5/TPSPK)**CAPA TTHESK=TTHBTK*EXP(ELOCP*QBTK*APESPK/TTHBTK) C--------------CHECK FOR MAXIMUM THETA E-------------------------------- -C print*,' DEBUG HS calcape PTBL(IQ ,IT )',PTBL(IQ ,IT ), C 2 PTBL(IQ+1 ,IT ),PTBL(IQ ,IT+1 ),PTBL(IQ+1 ,IT+1 ) -C print*,' DEBUG HS calcape TTHESK.GT.THESP(I,J:',TTHESK,THESP(I,J) IF(TTHESK.GT.THESP(I,J)) THEN -C print*,' DEBUG HS calcape Update PSP:',TPSPK,' LMM:',LMM PSP (I,J)=TPSPK THESP(I,J)=TTHESK ENDIF @@ -253,9 +248,7 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, DO J=1,JM DO I=1,IM IF (L.LT.LMH(I,J)) THEN -C print*,' DEBUG HS calcape L.LT.LMH(I,J):',L,LMH(I,J) PKL = P(I,J,L) -C print*,' DEBUG HS calcape L,PKL.LT.PSP(I,J)',L,PKL,PSP(I,J) IF (PKL.LT.PSP(I,J)) THEN LCL(I,J)=L+1 PLCL(I,J)=P(I,J,L+1) @@ -276,7 +269,6 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, DO J=1,JM DO I=1,IM PKL = P(I,J,L) -C print*,' DEBUG HS calcape L.LE.LCL(I,J):',L,LCL(I,J) IF(L.LE.LCL(I,J)) THEN IF(PKL.LT.PLQ)THEN KNUML=KNUML+1 @@ -311,23 +303,19 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, ENDIF C------------SEARCH FOR EQ LEVEL---------------------------------------- -C print*,' DEBUG HS calcape SEARCH FOR EQ LEVEL, KNUMH:',KNUMH DO N=1,KNUMH I=IHRES(N) J=JHRES(N) IF(TPAR(I,J,L).GT.T(I,J,L)) THEN IEQL(I,J)=L -C print*,' DEBUG HS calcape 111 IEQL(I,J):',IEQL(I,J) PEQL(I,J)=P(I,J,L) ENDIF ENDDO -C print*,' DEBUG HS calcape SEARCH FOR EQ LEVEL, KNUML:',KNUML DO N=1,KNUML I=ILRES(N) J=JLRES(N) IF(TPAR(I,J,L).GT.T(I,J,L)) THEN IEQL(I,J)=L -C print*,' DEBUG HS calcape 222 IEQL(I,J):',IEQL(I,J) PEQL(I,J)=P(I,J,L) ENDIF ENDDO @@ -339,7 +327,6 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, DO I=1,IM LCLK=LCL(I,J) IEQK=IEQL(I,J) -C print*,' DEBUG HS calcape IEQK,LCLK:',IEQK,LCLK DO L=IEQK,LCLK c print*,'l=',l c print*,'p(i,j,l)=',p(i,j,l) @@ -369,12 +356,10 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, C print*,'virtual thetap=',thetap C print*,'virtual thetaa=',thetaa endif -C print*,' DEBUG HS calcape before CAPE(I,J)=',CAPE(I,J) IF (THETAP.LT.THETAA) & CINS(I,J)=CINS(I,J)+G*(ALOG(THETAP)-ALOG(THETAA))*DZKL IF (THETAP.GT.THETAA) & CAPE(I,J)=CAPE(I,J)+G*(ALOG(THETAP)-ALOG(THETAA))*DZKL -C print*,' DEBUG HS calcape after CAPE(I,J)=',CAPE(I,J) ENDDO ENDDO ENDDO @@ -389,7 +374,6 @@ SUBROUTINE CALCAPE(ivirt,itype,T,Q,P,p1d,t1d,q1d,PINT,LMH, DO I = 1,IM CAPE(I,J) = AMAX1(0.0,CAPE(I,J)) CINS(I,J) = AMIN1(CINS(I,J),0.0) -C print*,' DEBUG HS calcape final CAPE(I,J)=',CAPE(I,J) ENDDO ENDDO c if(itype.eq.2) print*,'CAPE,CINS=',CAPE,CINS @@ -791,19 +775,7 @@ SUBROUTINE SPLINE(JTB,NOLD,XOLD,YOLD,Y2,NNEW,XNEW,YNEW,P,Q) 600 K1=K1+1 IF(K1.LE.NNEW) GO TO 300 C----------------------------------------------------------------------- -C print*,' DEBUG HS calcape NOLD,NNEW: ',NOLD,NNEW -CC print*,' DEBUG HS calcape XOLD: ',XOLD -CC print*,' DEBUG HS calcape YOLD: ',YOLD -C print*,' DEBUG HS calcape XOLD(1)->: XNEW(1)',XOLD(1),XNEW(1) -C print*,' DEBUG HS calcape XOLD(2)->: XNEW(2)',XOLD(2),XNEW(2) -CC print*,' DEBUG HS calcape XNEW: ',XNEW -CC print*,' DEBUG HS calcape YNEW: ',YNEW -C print*,' DEBUG HS calcape YOLD(1)->: YNEW(1)',YOLD(1),YNEW(1) -C print*,' DEBUG HS calcape YOLD(2)->: YNEW(2)',YOLD(2),YNEW(2) -C print*,' DEBUG HS calcape YOLD(3)->: YNEW(3)',YOLD(3),YNEW(3) -C print*,' DEBUG HS calcape YOLD(-2)->: YNEW(-2)', C 2 YOLD(NOLD-1),YNEW(NNEW-1) -C print*,' DEBUG HS calcape YOLD(-1)->: YNEW(-1)', C 2 YOLD(NOLD),YNEW(NNEW) RETURN END From 5dffde227fb6bb9249c2b1f1edeeb9eb38bdedb5 Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Thu, 3 Feb 2022 13:51:42 -0700 Subject: [PATCH 3/8] #1824 Support mixed layer cape --- met/src/tools/other/pb2nc/pb2nc.cc | 207 ++++++++++++++++++++++------- 1 file changed, 160 insertions(+), 47 deletions(-) diff --git a/met/src/tools/other/pb2nc/pb2nc.cc b/met/src/tools/other/pb2nc/pb2nc.cc index 9190b4106f..b063b613ce 100644 --- a/met/src/tools/other/pb2nc/pb2nc.cc +++ b/met/src/tools/other/pb2nc/pb2nc.cc @@ -260,8 +260,11 @@ static const char *airnow_aux_vars = "TPHR QCIND"; // Pick the latter one if exists multiuple variables static const char *bufr_avail_sid_names = "SID SAID RPID"; static const char *bufr_avail_latlon_names = "XOB CLON CLONH YOB CLAT CLATH"; +static const char *derived_mlcape = "D_MLCAPE"; static const char *derived_cape = "D_CAPE"; static const char *derived_pbl = "D_PBL"; +static const float MLCAPE_INTERVAL = 3000.; +static bool use_CAPE_BIN = false; static double bufr_obs[mxr8lv][mxr8pm]; static double bufr_obs_extra[mxr8lv][mxr8pm]; @@ -511,6 +514,7 @@ void initialize() { prepbufr_derive_vars.add("D_MIXR"); prepbufr_derive_vars.add("D_PRMSL"); prepbufr_derive_vars.add("D_CAPE"); + prepbufr_derive_vars.add("D_MLCAPE"); prepbufr_derive_vars.add("D_PBL"); for (int idx=0; idx<(sizeof(hdr) / sizeof(hdr[0])); idx++) { @@ -982,10 +986,13 @@ void process_pbfile(int i_pb) { int itype = 1; // itype 1: where a parcel is lifted from the ground // itype 2: Where the "best cape" in a number of parcels int cape_code = -1; + int mlcape_code = -1; float p1d,t1d,q1d; int IMM, JMM; - int cape_level=0, cape_count=0, cape_cnt_too_big=0, cape_cnt_surface_msgs=0; - int cape_cnt_no_levels=0, cape_cnt_missing_values=0, cape_cnt_zero_values=0; + int cape_level, cape_count, cape_cnt_too_big, cape_cnt_surface_msgs; + int cape_cnt_no_levels, cape_cnt_missing_values, cape_cnt_zero_values; + int mlcape_count, mlcape_cnt_too_big; + int mlcape_cnt_missing_values, mlcape_cnt_zero_values; float cape_p, cape_h; float cape_qm = bad_data_float; @@ -997,8 +1004,9 @@ void process_pbfile(int i_pb) { bool has_pbl_data; bool do_pbl = false; - bool cal_cape = bufr_obs_name_arr.has(derived_cape, cape_code, false); bool cal_pbl = bufr_obs_name_arr.has(derived_pbl, pbl_code, false); + bool cal_cape = bufr_obs_name_arr.has(derived_cape, cape_code, false); + bool cal_mlcape = bufr_obs_name_arr.has(derived_mlcape, mlcape_code, false); bool is_same_header; unixtime prev_hdr_vld_ut = (unixtime) 0; @@ -1010,6 +1018,10 @@ void process_pbfile(int i_pb) { // Initialize prev_hdr_lat = prev_hdr_lon = prev_hdr_elv = bad_data_double; + cape_level = cape_count = cape_cnt_too_big = 0; + cape_cnt_no_levels = cape_cnt_missing_values = cape_cnt_zero_values = 0; + mlcape_count = mlcape_cnt_too_big = 0; + mlcape_cnt_missing_values = mlcape_cnt_zero_values = 0; if (cal_pbl) { is_same_header = false; @@ -1327,7 +1339,7 @@ void process_pbfile(int i_pb) { } } - if (cal_cape) { + if (cal_cape || cal_mlcape) { cape_level = 0; } @@ -1471,7 +1483,7 @@ void process_pbfile(int i_pb) { obs_arr[4] += c_to_k; } - if (cal_cape) { + if (cal_cape || cal_mlcape) { if(grib_code == pres_grib_code) { is_cape_input = true; if (cape_level < MAX_CAPE_LEVEL) cape_data_pres[cape_level] = obs_arr[4]; @@ -1567,7 +1579,7 @@ void process_pbfile(int i_pb) { } // end for i } - if (cal_cape) { + if (cal_cape || cal_mlcape) { if (cape_member_cnt >= 3) cape_level++; } @@ -1603,7 +1615,7 @@ void process_pbfile(int i_pb) { } } // end for lv - if (cal_cape) { + if (cal_cape || cal_mlcape) { if (1 < cape_level) { bool reverse_levels; float cape_val, cin_val, PLCL,PEQL; @@ -1646,14 +1658,6 @@ void process_pbfile(int i_pb) { } } - //p1d = cape_p; - //t1d = cape_data_temp[cape_level-1]; - //q1d = cape_data_spfh[cape_level-1]; - calcape_(&ivirt,&itype, cape_data_temp, cape_data_spfh, cape_data_pres, - &p1d,&t1d,&q1d, static_dummy_201, - &cape_level, &IMM,&JMM, &cape_level, - &cape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); - if(mlog.verbosity_level() >= 7) { mlog << Debug(7) << method_name << " index,P,T,Q to compute CAPE from " << i_read << "-th message\n" ; @@ -1661,37 +1665,146 @@ void process_pbfile(int i_pb) { mlog << Debug(7) << method_name << " " << idx << ", " << cape_data_pres[idx] << ", " << cape_data_temp[idx] << ", " << cape_data_spfh[idx] << "\n"; } - mlog << Debug(7) << method_name - << " calcape_(" << ivirt << "," << itype << ") cape_val: " - << cape_val << " cape_level: " << cape_level - << ", cin_val: " << cin_val - << ", PLCL: " << PLCL << ", PEQL: " << PEQL - << " lat: " << hdr_lat << ", lon: " << hdr_lon - << " valid_time: " << unix_to_yyyymmdd_hhmmss(hdr_vld_ut) - << " " << hdr_typ << " " << hdr_sid - << "\n\n" ; } - if (cape_val > MAX_CAPE_VALUE) { - cape_cnt_too_big++; - mlog << Debug(5) << method_name - << " Ignored cape_value: " << cape_val << " cape_level: " << cape_level - << ", cin_val: " << cin_val - << ", PLCL: " << PLCL << ", PEQL: " << PEQL << "\n"; + if (cal_cape) { + calcape_(&ivirt,&itype, cape_data_temp, cape_data_spfh, cape_data_pres, + &p1d,&t1d,&q1d, static_dummy_201, + &cape_level, &IMM,&JMM, &cape_level, + &cape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); + + if(mlog.verbosity_level() >= 7) { + mlog << Debug(7) << method_name + << " calcape_(" << ivirt << "," << itype << ") cape_val: " + << cape_val << " cape_level: " << cape_level + << ", cin_val: " << cin_val + << ", PLCL: " << PLCL << ", PEQL: " << PEQL + << " lat: " << hdr_lat << ", lon: " << hdr_lon + << " valid_time: " << unix_to_yyyymmdd_hhmmss(hdr_vld_ut) + << " " << hdr_typ << " " << hdr_sid + << "\n\n" ; + } + + if (cape_val > MAX_CAPE_VALUE) { + cape_cnt_too_big++; + mlog << Debug(5) << method_name + << " Ignored cape_value: " << cape_val << " cape_level: " << cape_level + << ", cin_val: " << cin_val + << ", PLCL: " << PLCL << ", PEQL: " << PEQL << "\n"; + } + else if (cape_val >= 0) { + obs_arr[1] = cape_code; + obs_arr[2] = cape_p; + obs_arr[3] = cape_h; + obs_arr[4] = cape_val; // observation value + addObservation(obs_arr, (string)hdr_typ, (string)hdr_sid, hdr_vld_ut, + hdr_lat, hdr_lon, hdr_elv, cape_qm, + OBS_BUFFER_SIZE); + cape_count++; + n_derived_obs++; + if (is_eq(cape_val, 0.)) cape_cnt_zero_values++; + } + else cape_cnt_missing_values++; } - else if (cape_val >= 0) { - obs_arr[1] = cape_code; - obs_arr[2] = cape_p; - obs_arr[3] = cape_h; - obs_arr[4] = cape_val; // observation value - addObservation(obs_arr, (string)hdr_typ, (string)hdr_sid, hdr_vld_ut, - hdr_lat, hdr_lon, hdr_elv, cape_qm, - OBS_BUFFER_SIZE); - cape_count++; - n_derived_obs++; - if (is_eq(cape_val, 0.)) cape_cnt_zero_values++; + if (cal_mlcape) { + int input_levels; + float mlcape_val = bad_data_float; + float *p_values, *t_values, *q_values; + float *P_bins = 0; + float *T_bins = 0; + float *Q_bins = 0; + + input_levels = cape_level; + if (use_CAPE_BIN) { + int p_levels = (cape_data_pres[cape_level-1] - cape_data_pres[0]) / MLCAPE_INTERVAL + 1; + if (p_levels > 1) { + float p_upper; + int idx_mlcape, counts[p_levels]; + P_bins = new float [p_levels]; + T_bins = new float [p_levels]; + Q_bins = new float [p_levels]; + + idx_mlcape = 0; + p_upper = MLCAPE_INTERVAL; + for (int idx=0; idx p_upper) { + for (int idx2=1; idx2 0) { + for (int idx=0; idx 0) { + P_bins[idx] /= counts[idx]; + T_bins[idx] /= counts[idx]; + Q_bins[idx] /= counts[idx]; + } + } + } + p_values = P_bins; + t_values = T_bins; + q_values = Q_bins; + } + } + else { + p_values = (float *)cape_data_pres; + t_values = (float *)cape_data_temp; + q_values = (float *)cape_data_spfh; + } + + if(p_values != 0) { + // The second last seems to be better than the average of last two or three + p1d = p_values[input_levels-2]; + t1d = t_values[input_levels-2]; + q1d = q_values[input_levels-2]; + ivirt = 0; + itype = 2; + calcape_(&ivirt,&itype, t_values, q_values, p_values, + &p1d,&t1d,&q1d, static_dummy_201, + &input_levels, &IMM,&JMM, &input_levels, + &cape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); + ivirt = 1; + calcape_(&ivirt,&itype, t_values, q_values, p_values, + &p1d,&t1d,&q1d, static_dummy_201, + &input_levels, &IMM,&JMM, &input_levels, + &mlcape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); + if (mlcape_val>0 && cape_val>0) mlcape_val = (mlcape_val + cape_val) / 2; + } + + if (mlcape_val > MAX_CAPE_VALUE) { + mlcape_cnt_too_big++; + mlog << Debug(5) << method_name + << " Ignored cape_value: " << cape_val << " cape_level: " << cape_level + << ", cin_val: " << cin_val + << ", PLCL: " << PLCL << ", PEQL: " << PEQL << "\n"; + } + else if (mlcape_val >= 0) { + obs_arr[1] = mlcape_code; + obs_arr[2] = cape_p; + obs_arr[3] = cape_h; + obs_arr[4] = mlcape_val; // observation value + addObservation(obs_arr, (string)hdr_typ, (string)hdr_sid, hdr_vld_ut, + hdr_lat, hdr_lon, hdr_elv, cape_qm, + OBS_BUFFER_SIZE); + mlcape_count++; + n_derived_obs++; + if (is_eq(mlcape_val, 0.)) mlcape_cnt_zero_values++; + } + else mlcape_cnt_missing_values++; } - else cape_cnt_missing_values++; } else if (1 < buf_nlev) cape_cnt_no_levels++; else cape_cnt_surface_msgs++; @@ -1948,13 +2061,13 @@ void process_pbfile(int i_pb) { << "Total observations retained or derived\t= " << (n_file_obs + n_derived_obs) << "\n"; - if (cal_cape) { - mlog << Debug(3) << "\nDerived CAPE = " << cape_count - << "\tZero = " << cape_cnt_zero_values + if (cal_cape || cal_mlcape) { + mlog << Debug(3) << "\nDerived CAPE = " << (cape_count + mlcape_count) + << "\tZero = " << (cape_cnt_zero_values + mlcape_cnt_zero_values) << "\n\tnot derived: No cape inputs = " << (cape_cnt_no_levels) << "\tNo vertical levels = " << (cape_cnt_surface_msgs) - << "\n\tfiltered: " << cape_cnt_missing_values << ", " - << cape_cnt_too_big + << "\n\tfiltered: " << (cape_cnt_missing_values + mlcape_cnt_missing_values) << ", " + << (cape_cnt_too_big + mlcape_cnt_too_big) << "\n"; } From 254ebf7e2310881a483d315b814f503bc88a832d Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Thu, 3 Feb 2022 14:00:43 -0700 Subject: [PATCH 4/8] #1824 Support mixed layer cape --- met/src/tools/other/pb2nc/pb2nc.cc | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/met/src/tools/other/pb2nc/pb2nc.cc b/met/src/tools/other/pb2nc/pb2nc.cc index b063b613ce..5e57694bf7 100644 --- a/met/src/tools/other/pb2nc/pb2nc.cc +++ b/met/src/tools/other/pb2nc/pb2nc.cc @@ -1804,6 +1804,10 @@ void process_pbfile(int i_pb) { if (is_eq(mlcape_val, 0.)) mlcape_cnt_zero_values++; } else mlcape_cnt_missing_values++; + if (P_bins != 0) delete [] P_bins; + if (T_bins != 0) delete [] T_bins; + if (Q_bins != 0) delete [] Q_bins; + } } else if (1 < buf_nlev) cape_cnt_no_levels++; From 432f7e0956d8edb115ebadaaf1d1841020d9ce11 Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Thu, 17 Feb 2022 17:24:24 -0700 Subject: [PATCH 5/8] #1824 ci-run-test Apply ivirt 1 and itype 2 only for MLCAPE and cleanup --- met/src/tools/other/pb2nc/pb2nc.cc | 101 +++++------------------------ 1 file changed, 16 insertions(+), 85 deletions(-) diff --git a/met/src/tools/other/pb2nc/pb2nc.cc b/met/src/tools/other/pb2nc/pb2nc.cc index b7b841df05..90f15b40b0 100644 --- a/met/src/tools/other/pb2nc/pb2nc.cc +++ b/met/src/tools/other/pb2nc/pb2nc.cc @@ -264,7 +264,6 @@ static const char *derived_mlcape = "D_MLCAPE"; static const char *derived_cape = "D_CAPE"; static const char *derived_pbl = "D_PBL"; static const float MLCAPE_INTERVAL = 3000.; -static bool use_CAPE_BIN = false; static double bufr_obs[mxr8lv][mxr8pm]; static double bufr_obs_extra[mxr8lv][mxr8pm]; @@ -1706,90 +1705,25 @@ void process_pbfile(int i_pb) { } else cape_cnt_missing_values++; } + if (cal_mlcape) { - int input_levels; float mlcape_val = bad_data_float; - float *p_values, *t_values, *q_values; - float *P_bins = 0; - float *T_bins = 0; - float *Q_bins = 0; - - input_levels = cape_level; - if (use_CAPE_BIN) { - int p_levels = (cape_data_pres[cape_level-1] - cape_data_pres[0]) / MLCAPE_INTERVAL + 1; - if (p_levels > 1) { - float p_upper; - int idx_mlcape, counts[p_levels]; - P_bins = new float [p_levels]; - T_bins = new float [p_levels]; - Q_bins = new float [p_levels]; - - idx_mlcape = 0; - p_upper = MLCAPE_INTERVAL; - for (int idx=0; idx p_upper) { - for (int idx2=1; idx2 0) { - for (int idx=0; idx 0) { - P_bins[idx] /= counts[idx]; - T_bins[idx] /= counts[idx]; - Q_bins[idx] /= counts[idx]; - } - } - } - p_values = P_bins; - t_values = T_bins; - q_values = Q_bins; - } - } - else { - p_values = (float *)cape_data_pres; - t_values = (float *)cape_data_temp; - q_values = (float *)cape_data_spfh; - } - - if(p_values != 0) { - // The second last seems to be better than the average of last two or three - p1d = p_values[input_levels-2]; - t1d = t_values[input_levels-2]; - q1d = q_values[input_levels-2]; - ivirt = 0; - itype = 2; - calcape_(&ivirt,&itype, t_values, q_values, p_values, - &p1d,&t1d,&q1d, static_dummy_201, - &input_levels, &IMM,&JMM, &input_levels, - &cape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); - ivirt = 1; - calcape_(&ivirt,&itype, t_values, q_values, p_values, - &p1d,&t1d,&q1d, static_dummy_201, - &input_levels, &IMM,&JMM, &input_levels, - &mlcape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); - if (mlcape_val>0 && cape_val>0) mlcape_val = (mlcape_val + cape_val) / 2; - } + ivirt = 1; + itype = 2; + // The second last seems to be better than the average of last two or three + p1d = cape_data_pres[cape_level-2]; + t1d = cape_data_temp[cape_level-2]; + q1d = cape_data_spfh[cape_level-2]; + calcape_(&ivirt,&itype, cape_data_temp, cape_data_spfh, cape_data_pres, + &p1d,&t1d,&q1d, static_dummy_201, + &cape_level, &IMM,&JMM, &cape_level, + &mlcape_val, &cin_val, &PLCL, &PEQL, static_dummy_200); if (mlcape_val > MAX_CAPE_VALUE) { mlcape_cnt_too_big++; mlog << Debug(5) << method_name - << " Ignored cape_value: " << cape_val << " cape_level: " << cape_level - << ", cin_val: " << cin_val - << ", PLCL: " << PLCL << ", PEQL: " << PEQL << "\n"; + << " Ignored ML_cape: " << mlcape_val << " cape_level: " + << cape_level << "\n"; } else if (mlcape_val >= 0) { obs_arr[1] = mlcape_code; @@ -1804,10 +1738,6 @@ void process_pbfile(int i_pb) { if (is_eq(mlcape_val, 0.)) mlcape_cnt_zero_values++; } else mlcape_cnt_missing_values++; - if (P_bins != 0) delete [] P_bins; - if (T_bins != 0) delete [] T_bins; - if (Q_bins != 0) delete [] Q_bins; - } } else if (1 < buf_nlev) cape_cnt_no_levels++; @@ -3118,9 +3048,10 @@ float compute_pbl(map pqtzuv_map_tq, << " from TQ (" << tq_count << ") and UV (" << uv_count << ")\n"; - if (pbl_level <= 0) { + if (pbl_level <= 1) { mlog << Debug(4) << method_name - << "Skip CALPBL because an empty list after combining TQZ and UV\n"; + << "Skip CALPBL because of " << (pbl_level == 1) ? "only one record " : "an empty list" + << " after combining TQZ and UV\n"; } else { // Order all observations by pressure from bottom to top From 656a103c41f6a6ff1cb856a42e2a1af8e274451a Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Fri, 18 Feb 2022 08:48:12 -0700 Subject: [PATCH 6/8] #1824 Added MLCAPE to obs_bufr_map --- met/docs/Users_Guide/reformat_point.rst | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/met/docs/Users_Guide/reformat_point.rst b/met/docs/Users_Guide/reformat_point.rst index 7a3fce768f..0953afdd7c 100644 --- a/met/docs/Users_Guide/reformat_point.rst +++ b/met/docs/Users_Guide/reformat_point.rst @@ -227,7 +227,7 @@ _____________________ obs_bufr_var = [ 'QOB', 'TOB', 'ZOB', 'UOB', 'VOB' ]; -Each PrepBUFR message will likely contain multiple observation variables. The **obs_bufr_var** variable is used to specify which observation variables should be retained or derived. The variable name comes from BUFR file which includes BUFR table. The following BUFR names may be retained: QOB, TOB, ZOB, UOB, and VOB for specific humidity, temperature, height, and the u and v components of winds. The following BUFR names may be derived: D_DPT, D_WIND, D_RH, D_MIXR, D_PRMSL, D_PBL, and D_CAPE for dew point, wind speed, relative humidity, mixing ratio, pressure reduced to MSL, planetary boundary layer height, and convective available potential energy. This configuration replaces **obs_grib_code**. If the list is empty, all BUFR variables are retained. +Each PrepBUFR message will likely contain multiple observation variables. The **obs_bufr_var** variable is used to specify which observation variables should be retained or derived. The variable name comes from BUFR file which includes BUFR table. The following BUFR names may be retained: QOB, TOB, ZOB, UOB, and VOB for specific humidity, temperature, height, and the u and v components of winds. The following BUFR names may be derived: D_DPT, D_WIND, D_RH, D_MIXR, D_PRMSL, D_PBL, D_CAPE, and D_MLCAPE for dew point, wind speed, relative humidity, mixing ratio, pressure reduced to MSL, planetary boundary layer height, convective available potential energy, and mixed layer convective available potential energy. This configuration replaces **obs_grib_code**. If the list is empty, all BUFR variables are retained. _____________________ @@ -248,6 +248,7 @@ _____________________ { key = 'D_PRMSL'; val = 'PRMSL'; }, { key = 'D_PBL'; val = 'PBL'; }, { key = 'D_CAPE'; val = 'CAPE'; } + { key = 'D_MLCAPE'; val = 'MLCAPE'; } ]; From bdf70b75bb5e12fea4425859629347c67cf0168e Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Fri, 18 Feb 2022 08:50:29 -0700 Subject: [PATCH 7/8] #1824 ci-run-test FIx for the Docker build --- met/src/tools/other/pb2nc/pb2nc.cc | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/met/src/tools/other/pb2nc/pb2nc.cc b/met/src/tools/other/pb2nc/pb2nc.cc index 90f15b40b0..0b914b3198 100644 --- a/met/src/tools/other/pb2nc/pb2nc.cc +++ b/met/src/tools/other/pb2nc/pb2nc.cc @@ -3048,10 +3048,13 @@ float compute_pbl(map pqtzuv_map_tq, << " from TQ (" << tq_count << ") and UV (" << uv_count << ")\n"; - if (pbl_level <= 1) { + if (pbl_level <= 0) { mlog << Debug(4) << method_name - << "Skip CALPBL because of " << (pbl_level == 1) ? "only one record " : "an empty list" - << " after combining TQZ and UV\n"; + << "Skip CALPBL because of an empty list after combining TQZ and UV\n"; + } + else if (pbl_level == 1) { + mlog << Debug(4) << method_name + << "Skip CALPBL because of only one available record after combining TQZ and UV\n"; } else { // Order all observations by pressure from bottom to top From f8d6888cb7e24dad594f0e7f6c455b6cfe02aec0 Mon Sep 17 00:00:00 2001 From: Howard Soh Date: Fri, 18 Feb 2022 11:12:43 -0700 Subject: [PATCH 8/8] #1824 Added MLCAPE output into the existing unit test pb2nc_compute_pbl_capepb2nc_compute_pbl_capepb2nc_compute_pbl_cape --- test/config/PB2NCConfig_pbl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/test/config/PB2NCConfig_pbl b/test/config/PB2NCConfig_pbl index 924c06a54c..ac48b5850f 100644 --- a/test/config/PB2NCConfig_pbl +++ b/test/config/PB2NCConfig_pbl @@ -92,7 +92,7 @@ level_category = [0, 1, 4, 5, 6]; // Use obs_bufr_map to rename variables in the output. // If empty, process all available variables. // -obs_bufr_var = ["D_CAPE", "D_PBL"]; +obs_bufr_var = ["D_CAPE", "D_PBL", "D_MLCAPE"]; //////////////////////////////////////////////////////////////////////////////// // @@ -120,6 +120,7 @@ obs_prepbufr_map = [ { key = "D_PBL"; val = "HPBL"; }, { key = "D_PRMSL"; val = "PRMSL"; }, { key = "D_CAPE"; val = "CAPE"; }, + { key = "D_MLCAPE"; val = "MLCAPE"; }, { key = "TDO"; val = "DPT"; }, { key = "PMO"; val = "PRMSL"; }, { key = "TOCC"; val = "TCDC"; },