Skip to content

Commit

Permalink
minor bug fixes in oi_ensi_lr. var name "_debug_level" changed becaus…
Browse files Browse the repository at this point in the history 
…e variables in R can't begin with "_" and this was creating issues with the R bindings
  • Loading branch information
Cristian Lussana committed Sep 3, 2024
1 parent 29e1890 commit 03372e5
Show file tree
Hide file tree
Showing 3 changed files with 51 additions and 20 deletions.
35 changes: 33 additions & 2 deletions include/gridpp.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -314,7 +314,38 @@ namespace gridpp {
const vec3& background_l,
const vec3& background_L,
const Points& obs_points,
const vec& obs,
const vec2& obs,
const vec2& pbackground_r,
const vec2& pbackground_R,
const StructureFunction& structure,
float var_ratios_or,
float std_ratios_lr,
float weigth,
int max_points,
bool allow_extrapolation=true);

/** Optimal interpolation for an ensemble gridded field (alternative version)
* Work in progress
* @param bgrid Grid of background field
* @param background 3D vector of (left) background values to update (Y, X, E)
* @param background 3D vector of (LEFT) background values (Y, X, E) used to compute correlations
* @param obs_points observation points
* @param obs 2D vector of perturbed observations (S, E)
* @param background 3D vector of (right) background values used to compute innovations (Y, X, E)
* @param background 3D vector of (RIGHT) background values (Y, X, E) used to compute correlations
* @param structure Structure function
* @param variance_ratio (ratio of observation to right background error variance)
* @param standard deviation ratio (ratio of left to right background error standard deviation)
* @param weight given to the analysis increment
* @param max_points Maximum number of observations to use inside localization zone; Use 0 to disable
* @param allow_extrapolation Allow OI to extrapolate increments outside increments at observations
* @returns 3D vector of analised values (Y, X, E)
*/
vec2 optimal_interpolation_ensi_lr(const Points& bpoints,
const vec2& background_l,
const vec2& background_L,
const Points& obs_points,
const vec2& obs,
const vec2& pbackground_r,
const vec2& pbackground_R,
const StructureFunction& structure,
Expand Down Expand Up @@ -1256,7 +1287,7 @@ namespace gridpp {
*/
void set_debug_level(int level);

static int _debug_level = 0;
static int debug_level = 0;

/** Get the currently set level of debug messagess
* @returns Currently set debug level
Expand Down
4 changes: 2 additions & 2 deletions src/api/gridpp.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -68,11 +68,11 @@ int gridpp::get_omp_threads() {
}

void gridpp::set_debug_level(int level) {
gridpp::_debug_level = level;
gridpp::debug_level = level;
}

int gridpp::get_debug_level() {
return gridpp::_debug_level;
return gridpp::debug_level;
}
gridpp::not_implemented_exception::not_implemented_exception() :
std::logic_error("Function not yet implemented") {
Expand Down
32 changes: 16 additions & 16 deletions src/api/oi_ensi_lr.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,32 +66,32 @@ vec3 gridpp::optimal_interpolation_ensi_lr(const gridpp::Grid& bgrid,
throw std::invalid_argument("Both background grid and observations points must be of same coordinate type (lat/lon or x/y)");
}
// Check ensembles have consistent sizes
int nE = background_l[0][0].size();
if(background_l.size() != nY || background_l[0].size() != nX) {
std::stringstream ss;
ss << "Input left field (" << background_l.size() << "," << background_l[0].size() << "," background_l[0][0].size() << ") is not the same size as the grid (" << nY << "," << nX << "," << nE << ")";
ss << "Input left field (" << background_l.size() << "," << background_l[0].size() << "," << background_l[0][0].size() << ") is not the same size as the grid (" << nY << "," << nX << "," << nE << ")";
throw std::invalid_argument(ss.str());
}
int nE = background_l[0][0].size();
if(background_L.size() != nY || background_L[0].size() != nX || background_L[0][0].size() != nE) {
std::stringstream ss;
ss << "Input LEFT field (" << background_L.size() << "," << background_L[0].size() << "," background_L[0][0].size() << ") is not the same size as the grid (" << nY << "," << nX << "," << nE << ")";
ss << "Input LEFT field (" << background_L.size() << "," << background_L[0].size() << "," << background_L[0][0].size() << ") is not the same size as the grid (" << nY << "," << nX << "," << nE << ")";
throw std::invalid_argument(ss.str());
}
if(pbackground_r.size() != nS || pbackground_r[0].size() != nE) {
std::stringstream ss;
ss << "Input right field at observation location (" << pbackground_r.size() << "," << pbackground_r[0].size() << ") and points (" << nS "," << nE << ") size mismatch";
ss << "Input right field at observation location (" << pbackground_r.size() << "," << pbackground_r[0].size() << ") and points (" << nS << "," << nE << ") size mismatch";
throw std::invalid_argument(ss.str());
}

if(pbackground_R.size() != nS || pbackground_R[0].size() != nE) {
std::stringstream ss;
ss << "Input RIGHT field at observation location (" << pbackground_R.size() << "," << pbackground_R[0].size() << ") and points (" << nS "," << nE << ") size mismatch";
ss << "Input RIGHT field at observation location (" << pbackground_R.size() << "," << pbackground_R[0].size() << ") and points (" << nS << "," << nE << ") size mismatch";
throw std::invalid_argument(ss.str());
}
// Check observations have consistent size
if(pobs.size() != nS || pobs[0].size() != nE) {
std::stringstream ss;
ss << "Observations (" << pobs.size() << "," << pobs[0].size() << ") and points (" << nS "," << nE << ") size mismatch";
ss << "Observations (" << pobs.size() << "," << pobs[0].size() << ") and points (" << nS << "," << nE << ") size mismatch";
throw std::invalid_argument(ss.str());
}

Expand Down Expand Up @@ -195,8 +195,8 @@ vec2 gridpp::optimal_interpolation_ensi_lr(const gridpp::Points& bpoints,
numInvalid++;
}
for(int i = 0; i < nS; i++) {
float value_r = pbackground_r[y][e];
float value_R = pbackground_R[y][e];
float value_r = pbackground_r[i][e];
float value_R = pbackground_R[i][e];
if(!gridpp::is_valid(value_r) || !gridpp::is_valid(value_R))
numInvalid++;
}
Expand Down Expand Up @@ -267,9 +267,9 @@ vec2 gridpp::optimal_interpolation_ensi_lr(const gridpp::Points& bpoints,
vec rhos = structure.corr_background(p1, p2);
for(int i = 0; i < lLocIndices0.size(); i++) {
int index = lLocIndices0[i];
if(gridpp::is_valid(pobs[index])) {
if(rhos_static_lr[i] > 0) {
lRhos0.push_back(std::pair<float,int>(rhos_static_lr[i], i));
if(gridpp::is_valid(pobs[index][0])) {
if(rhos[i] > 0) {
lRhos0.push_back(std::pair<float,int>(rhos[i], i));
}
}
}
Expand Down Expand Up @@ -322,7 +322,7 @@ vec2 gridpp::optimal_interpolation_ensi_lr(const gridpp::Points& bpoints,
float std = gridpp::calc_statistic(backgroundValid_L, gridpp::Std);
if(gridpp::is_valid(mean) && gridpp::is_valid(std) && std != 0) {
for(int e = 0; e < nValidEns; e++)
lX_L[0][e] = 1 / sqrt(nValidEns-1) * (backgroundValid_L[e] - mean) / std;
lX_L(0,e) = 1 / sqrt(nValidEns-1) * (backgroundValid_L[e] - mean) / std;
}
// lZ_R: used to compute ensemble-based background correlations i) between yth gridpoint and observations ii) among observations
mattype lZ_R(lS, nValidEns);
Expand Down Expand Up @@ -372,13 +372,13 @@ vec2 gridpp::optimal_interpolation_ensi_lr(const gridpp::Points& bpoints,
// Anti-extrapolation filter //
///////////////////////////////
if(!allow_extrapolation) {
vectype MaxIncEns = arma::max(lInnov);
vectype MinIncEns = arma::min(lInnov);
vectype maxIncEns = arma::max(lInnov);
vectype minIncEns = arma::min(lInnov);

for(int e = 0; e < nValidEns; e++) {
float increment = dx[e];
float MaxInc = MaxIncEns[e];
float MinInc = MinIncEns[e];
float maxInc = maxIncEns[e];
float minInc = minIncEns[e];
if(maxInc > 0 && increment > maxInc) {
increment = maxInc;
}
Expand Down

0 comments on commit 03372e5

Please sign in to comment.