Skip to content

Commit

Permalink
ensure noisy gap filling is correctly applied
Browse files Browse the repository at this point in the history
- propagated fillval from the downloaded CDF to fill_gaps
- made sure fill_gaps' output data wasn't overwriting the input
  • Loading branch information
drsteve authored and jhaiduce committed Apr 24, 2020
1 parent 100ef78 commit d4603c5
Show file tree
Hide file tree
Showing 2 changed files with 11 additions and 4 deletions.
10 changes: 8 additions & 2 deletions advect_imf.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,10 @@ def load_acedata(tstart, tend, noise=True, proxy=None):
for dataset, local_name, cdaweb_name in [(mag_data, 'b', 'BGSM'),
(swepam_data, 'u', 'V_GSM'),
(swepam_data, '', 'SC_pos_GSM')]:
t, values = dataset['Epoch'], fill_gaps(dataset[cdaweb_name][:, i], noise=noise)
t = dataset['Epoch']
values = fill_gaps(dataset[cdaweb_name][:, i],
fillval=dataset[cdaweb_name].attrs['FILLVAL'],
noise=noise)

# Grab the appropriate component from
# VALIDMIN and VALIDMAX attributes
Expand Down Expand Up @@ -99,7 +102,10 @@ def load_dscovr(tstart, tend, noise=True, proxy=None):
(mag_data, 'b', 'B1GSE', 'Epoch1'),
(plasma_data, 'u', 'V_GSE', 'Epoch'),
(orbit_data, '', 'GSE_POS', 'Epoch')]:
t, values = dataset[cdaweb_time_var], fill_gaps(dataset[cdaweb_name][:, i], noise=noise)
t = dataset[cdaweb_time_var]
values = fill_gaps(dataset[cdaweb_name][:, i],
fillval=dataset[cdaweb_name].attrs['FILLVAL'],
noise=noise)

for attr in ('VALIDMIN', 'VALIDMAX'):

Expand Down
5 changes: 3 additions & 2 deletions missing.py
Original file line number Diff line number Diff line change
Expand Up @@ -60,11 +60,12 @@ def fill_gaps(data, fillval=9999999, sigma=5, winsor=0.05, noise=False, constrai
# draw fluctuations from CDF and apply to linearly filled gaps
for gap in gaps:
for i in range(gap[1]-gap[0]+1):
data[gap[0]+i] += dx[p.searchsorted(random.random())]
series[gap[0]+i] += dx[p.searchsorted(random.random())]

# cap variable if it should be strictly positive (e.g. number density)
# use lowest measured value as floor
if constrain and series.min() > 0.0:
data[data < series.min()] = series.min()
series[series < series.min()] = series.min()
return series

return data

0 comments on commit d4603c5

Please sign in to comment.