Skip to content

Commit

Permalink
replace for loop with vecrtorized approach
Browse files Browse the repository at this point in the history
  • Loading branch information
kkappler committed Jan 11, 2024
1 parent 25a60d4 commit 8c27a19
Showing 1 changed file with 8 additions and 18 deletions.
26 changes: 8 additions & 18 deletions aurora/transfer_function/weights/coherence_weights.py
Original file line number Diff line number Diff line change
Expand Up @@ -88,33 +88,23 @@ def jackknife_coherence_weights(
f"removing 'best' {n_obs - n_clip_high} of {n_obs} via jackknife coherence"
)

# Initialize a weight vector the length = num_observations
partial_coh = np.zeros(n_obs)
W = np.zeros(n_obs) # for example

c_xy = np.abs(x.conj() * y)
c_xx = np.real(x.conj() * x)
c_yy = np.real(y.conj() * y)

Check warning on line 93 in aurora/transfer_function/weights/coherence_weights.py

View check run for this annotation

Codecov / codecov/patch

aurora/transfer_function/weights/coherence_weights.py#L91-L93

Added lines #L91 - L93 were not covered by tests

# "Jackknife" matrix; all ones but zero on the diagonal
J = np.ones((n_obs, n_obs)) - np.eye(n_obs)
ccc = np.vstack([c_xy, c_xx, c_yy])

for i in range(n_obs):
partial_series = drop_column(ccc, i)
partial_coh[i] = coherence_from_fc_series(
partial_series[0, :],
partial_series[1, :],
partial_series[2, :],
)

# "Jackknife" matrix; all ones but zero on the dia may need some memory checks and a lazy approach with hdf5
# J = np.ones((n_obs, n_obs)) - np.eye(n_obs)
J = np.ones((n_obs, n_obs), dtype=np.int8) - np.eye(n_obs, dtype=np.int8)
p_xx = J @ c_xx.data
p_xy = J @ c_xy.data
p_yy = J @ c_yy.data
pc = p_xy**2 / (p_xx * p_yy)
assert np.isclose(pc, partial_coh).all()
p_xy**2 / (p_xx * p_yy)
partial_coh = p_xy**2 / (p_xx * p_yy)

Check warning on line 102 in aurora/transfer_function/weights/coherence_weights.py

View check run for this annotation

Codecov / codecov/patch

aurora/transfer_function/weights/coherence_weights.py#L97-L102

Added lines #L97 - L102 were not covered by tests

worst_to_best = np.argsort(partial_coh)
keepers = worst_to_best[n_clip_low:n_clip_high]

Check warning on line 105 in aurora/transfer_function/weights/coherence_weights.py

View check run for this annotation

Codecov / codecov/patch

aurora/transfer_function/weights/coherence_weights.py#L105

Added line #L105 was not covered by tests
# Initialize a weight vector the length = num_observations
W = np.zeros(n_obs) # for example

Check warning on line 107 in aurora/transfer_function/weights/coherence_weights.py

View check run for this annotation

Codecov / codecov/patch

aurora/transfer_function/weights/coherence_weights.py#L107

Added line #L107 was not covered by tests
W[keepers] = 1
return W

Expand Down

0 comments on commit 8c27a19

Please sign in to comment.