Commit d9bb5d76 authored by Peters, Wouter's avatar Peters, Wouter
Browse files

KG matrix redimensioned to use less memory

parent 584136a8
......@@ -79,7 +79,8 @@ class Optimizer(object):
self.speciesmask = {}
# Kalman Gain matrix
self.KG = np.zeros((self.nlag * self.nparams, self.nobs,), float)
#self.KG = np.zeros((self.nlag * self.nparams, self.nobs,), float)
self.KG = np.zeros((self.nlag * self.nparams,), float)
def state_to_matrix(self, statevector):
allsites = [] # collect all obs for n=1,..,nlag
......@@ -318,10 +319,12 @@ class Optimizer(object):
self.flags[n] = 2
continue
logging.debug('Proceeding to assimilate observation %s, %i' % (self.sitecode[n], self.obs_ids[n]))
PHt = 1. / (self.nmembers - 1) * np.dot(self.X_prime, self.HX_prime[n, :])
self.HPHR[n] = 1. / (self.nmembers - 1) * (self.HX_prime[n, :] * self.HX_prime[n, :]).sum() + self.R[n]
self.KG[:, n] = PHt / self.HPHR[n]
self.KG[:] = PHt / self.HPHR[n]
if self.may_localize[n]:
logging.debug('Trying to localize observation %s, %i' % (self.sitecode[n], self.obs_ids[n]))
......@@ -331,10 +334,10 @@ class Optimizer(object):
alpha = np.double(1.0) / (np.double(1.0) + np.sqrt((self.R[n]) / self.HPHR[n]))
self.x[:] = self.x + self.KG[:, n] * res
self.x[:] = self.x + self.KG[:] * res
for r in range(self.nmembers):
self.X_prime[:, r] = self.X_prime[:, r] - alpha * self.KG[:, n] * (self.HX_prime[n, r])
self.X_prime[:, r] = self.X_prime[:, r] - alpha * self.KG[:] * (self.HX_prime[n, r])
#WP !!!! Very important to first do all obervations from n=1 through the end, and only then update 1,...,n. The current observation
#WP should always be updated last because it features in the loop of the adjustments !!!!
......
......@@ -69,7 +69,7 @@ class CO2Optimizer(Optimizer):
corr = np.corrcoef(self.HX_prime[n, :], self.X_prime[r, :].squeeze())[0, 1]
prob = corr / np.sqrt((1.000000001 - corr ** 2) / (self.nmembers - 2))
if abs(prob) < self.tvalue:
self.KG[r, n] = 0.0
self.KG[r] = 0.0
count_localized = count_localized + 1
logging.debug('Localized observation %i, %i%% of values set to 0' % (self.obs_ids[n],count_localized*100/(self.nlag * self.nparams)))
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment