Skip to content

Commit caca96f

Browse files
committed
v496 - integrated EnKF fixes (M. Han from OPG)
-proper handling of windowsize to grab proper data in hydrograph, fix by M. Han (EnKF.cpp) -fixed indexing bug in CloseTimeStepOps (found by R. Chlumsky)
1 parent aeec3d6 commit caca96f

File tree

1 file changed

+10
-9
lines changed

1 file changed

+10
-9
lines changed

src/EnKF.cpp

Lines changed: 10 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -328,7 +328,7 @@ void CEnKFEnsemble::Initialize(const CModel* pModel,const optStruct &Options)
328328
{
329329
_aObsIndices[_nObs]=i;
330330
_nObs++;
331-
for(int nn=_nTimeSteps-_window_size;nn<_nTimeSteps;nn++) {
331+
for(int nn=_nTimeSteps-_window_size+1;nn<=_nTimeSteps;nn++) {
332332
obsval=pTSObs->GetSampledValue(nn);
333333
if(obsval!=RAV_BLANK_DATA) { _nObsDatapoints++; }
334334
}
@@ -380,7 +380,7 @@ void CEnKFEnsemble::Initialize(const CModel* pModel,const optStruct &Options)
380380
if (_pObsPerturbations[p]->state==sv){pPerturb=_pObsPerturbations[p];}
381381
}
382382

383-
for(int nn=_nTimeSteps-_window_size;nn<_nTimeSteps;nn++)
383+
for(int nn=_nTimeSteps-_window_size+1;nn<=_nTimeSteps;nn++)
384384
{
385385
obsval=pTSObs->GetSampledValue(nn);
386386
if(obsval!=RAV_BLANK_DATA) {
@@ -443,28 +443,29 @@ void CEnKFEnsemble::CloseTimeStepOps(CModel* pModel,optStruct& Options,const tim
443443
//----------------------------------------------------------
444444
int curr_nn=(int)((tt.model_time+TIME_CORRECTION)/Options.timestep);//current timestep index
445445

446-
if (curr_nn<_nTimeSteps-_window_size){return;} // haven't reached the data availability period yet
446+
if (curr_nn<_nTimeSteps-_window_size+1){return;} // haven't reached the data availability period yet
447447

448448
int j=0;
449449
double obsval;
450450
for(int ii=0;ii<_nObs;ii++) { //messy that we have to go through whole loop like this
451451
const CTimeSeriesABC* pTSObs=pModel->GetObservedTS (_aObsIndices[ii]);
452452
const CTimeSeriesABC* pTSMod=pModel->GetSimulatedTS(_aObsIndices[ii]);
453-
for(int nn=_nTimeSteps-_window_size;nn<_nTimeSteps;nn++) {
453+
for(int nn=_nTimeSteps-_window_size+1;nn<=_nTimeSteps;nn++) {
454454
obsval=pTSObs->GetSampledValue(nn);
455455
if(obsval!=RAV_BLANK_DATA) {
456456
if(nn==curr_nn) {
457457
_output_matrix[e][j]=pTSMod->GetSampledValue(nn);
458458
//cout<<"GRABBING OUTPUT ens: "<<e<<" nn="<<curr_nn<<" obs: "<<_output_matrix[e][j]<<" "<<obsval<<" duration "<<Options.duration<<endl;
459-
break;
459+
j++;
460+
break; //escape loop
461+
}
462+
else {
463+
j++;
460464
}
461-
j++;
465+
462466
}
463467
}
464468
}
465-
if (j!=_nObsDatapoints){
466-
ExitGracefully("Indexing issue in EnKF::CloseTimeStepOps: This shouldnt happen.",RUNTIME_ERR);
467-
}
468469
}
469470
//////////////////////////////////////////////////////////////////
470471
/// \brief the EnKF assimilation matrix calculations for updating states

0 commit comments

Comments
 (0)