diff --git a/src/Science/QSPCount.cpp b/src/Science/QSPCount.cpp index 4ccfc085a8c5eba6596638c338cfd9886cece1af..358eeaaa04dc7e8b0324737447acfcd2d67a03c0 100644 --- a/src/Science/QSPCount.cpp +++ b/src/Science/QSPCount.cpp @@ -265,18 +265,22 @@ void QSPCount(QSPOptions qsp_options, QSPData qsp_data) { return; } - TArrayn::DTArray *temp_ptr; - if (S1_ptr) { // If S1 is not ke - temp_ptr = S1_ptr; - } else { // If S1 is ke we know u must exist - temp_ptr = qsp_data.u; + int i_low, j_low, k_low, i_high, j_high, k_high; + + { + TArrayn::DTArray *temp_ptr; + if (S1_ptr) { // If S1 is not ke + temp_ptr = S1_ptr; + } else { // If S1 is ke we know u must exist + temp_ptr = qsp_data.u; + } + i_low = temp_ptr->lbound(blitz::firstDim); + j_low = temp_ptr->lbound(blitz::secondDim); + k_low = temp_ptr->lbound(blitz::thirdDim); + i_high = temp_ptr->ubound(blitz::firstDim); + j_high = temp_ptr->ubound(blitz::secondDim); + k_high = temp_ptr->ubound(blitz::thirdDim); } - int i_low = S1_ptr->lbound(blitz::firstDim); - int j_low = S1_ptr->lbound(blitz::secondDim); - int k_low = S1_ptr->lbound(blitz::thirdDim); - int i_high = S1_ptr->ubound(blitz::firstDim); - int j_high = S1_ptr->ubound(blitz::secondDim); - int k_high = S1_ptr->ubound(blitz::thirdDim); double double_max = std::numeric_limits<double>::max(); if (qsp_options.T1_max == double_max || qsp_options.S1_max == double_max || @@ -285,8 +289,10 @@ void QSPCount(QSPOptions qsp_options, QSPData qsp_data) { j_low, k_low, i_high, j_high, k_high); } - double hS = (qsp_options.S1_max - qsp_options.S1_min) / (double)qsp_options.NS; - double hT = (qsp_options.T1_max - qsp_options.T1_min) / (double)qsp_options.NT; + double hS = + (qsp_options.S1_max - qsp_options.S1_min) / (double)qsp_options.NS; + double hT = + (qsp_options.T1_max - qsp_options.T1_min) / (double)qsp_options.NT; double hS_inv = 1 / hS; double hT_inv = 1 / hT; @@ -375,10 +381,10 @@ void QSPCount(QSPOptions qsp_options, QSPData qsp_data) { break; } - int idxS = floor((Sval - S1_min) * hS_inv); - int idxT = floor((Tval - T1_min) * hT_inv); - idxS = std::max(std::min(idxS, qsp_options.NS), 0); - idxT = std::max(std::min(idxT, qsp_options.NT), 0); + int idxS = floor((Sval - qsp_options.S1_min) * hS_inv); + int idxT = floor((Tval - qsp_options.T1_min) * hT_inv); + idxS = std::max(std::min(idxS, qsp_options.NS - 1), 0); + idxT = std::max(std::min(idxT, qsp_options.NT - 1), 0); double volume_weight; if (qsp_data.mapped) { diff --git a/src/cases/derivatives/derivatives.cpp b/src/cases/derivatives/derivatives.cpp index b0c5a99d3f0ab62676d4bb6e220527a7f0ae6176..9fd3b0db1004198ee3f364d79674f5f6e57c2254 100644 --- a/src/cases/derivatives/derivatives.cpp +++ b/src/cases/derivatives/derivatives.cpp @@ -445,7 +445,7 @@ class userControl : public BaseCase { QSPOptions qsp_opts; qsp_opts.S1_name = S1_name; qsp_opts.T1_name = T1_name; - qsp_opts.filename = filename; + qsp_opts.filename = QSP_filename; qsp_opts.NS = NS; qsp_opts.NT = NT; qsp_opts.S1_max = S1_max; @@ -475,11 +475,11 @@ class userControl : public BaseCase { qsp_data.salinity = temp1; } if (T1_name.compare("temp") == 0 || S1_name.compare("temp") == 0) { - init_tracer_restart("T", *temp2); + init_tracer_restart("t", *temp2); qsp_data.temp = temp2; } if (T1_name.compare("rho") == 0 || S1_name.compare("rho") == 0) { - init_tracer_restart("t", *temp3); + init_tracer_restart("rho", *temp3); qsp_data.rho = temp3; } @@ -563,8 +563,8 @@ int main(int argc, char ** argv) { add_option("do_Q_and_R",&do_Q_and_R,false,"Calculate Q and R?"); add_option("do_lambda2",&do_lambda2,false,"Calculate Lambda2?"); add_option("do_hist",&do_hist,false,"Create QSP Data?"); - add_option("T1",&T1_name,"u", "Name of tracer 1 for QSP. Valid values are t (for rho),u,v,w,T (for temp) or k for K.E."); - add_option("S1",&S1_name,"w", "Name of tracer 2 for QSP. Valid values are t (for rho),u,v,w,T (for temp) or k for K.E."); + add_option("T1",&T1_name,"u", "Name of tracer 1 for QSP. Valid values are rho,u,v,w,temp or ke"); + add_option("S1",&S1_name,"w", "Name of tracer 2 for QSP. Valid values are rho,u,v,w,temp or ke"); add_option("T1_max",&T1_max,std::numeric_limits<double>::max(), "Maximum explicit bin for T1 in QSP."); add_option("T1_min",&T1_min,std::numeric_limits<double>::min(), "Minimum explicit bin for T1 in QSP."); add_option("S1_max",&S1_max,std::numeric_limits<double>::max(), "Maximum explicit bin for S1 in QSP.");