Skip to content

Commit

Permalink
Require phase_sd in degrees rather than radians
Browse files Browse the repository at this point in the history
  • Loading branch information
gavinmacaulay committed Dec 9, 2024
1 parent 05fe2cf commit ea396cf
Show file tree
Hide file tree
Showing 2 changed files with 4 additions and 4 deletions.
6 changes: 3 additions & 3 deletions src/echosms/dwbamodel.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

from .scattermodelbase import ScatterModelBase
from .utils import wavenumber, as_dict
from math import log10, cos, acos, pi, isclose
from math import log10, cos, acos, pi, isclose, radians
from cmath import exp
from scipy.spatial.transform import Rotation as R
import numpy as np
Expand Down Expand Up @@ -94,7 +94,7 @@ def calculate_ts_single(self, medium_c, medium_rho, theta, phi, f, target_c, tar
phase_sd : float
If non-zero, this model becomes the SDWBA (stochastic DWBA). A random phase is
applied to each term in the DWBA integral, obtained from a Gaussian distribution
centered on 0 with standard deviation of `phase_sd` [rad].
centered on 0 with standard deviation of `phase_sd` [°].
num_runs : int
The number of times to run the SDWBA model. The mean TS (calculated in the
linear domain) is returned. Intended to be used in conjunction with `phase_sd`.
Expand Down Expand Up @@ -168,7 +168,7 @@ def calculate_ts_single(self, medium_c, medium_rho, theta, phi, f, target_c, tar

for run in range(int(num_runs)): # is only ever 1 run for the DWBA
if do_sdwba:
phase_factors = np.exp(1j*self.rng.normal(scale=phase_sd, size=len(a)))
phase_factors = np.exp(1j*self.rng.normal(scale=radians(phase_sd), size=len(a)))

integral = 0.0
for a_, r_pos, d_r_pos, r_tan, p in zip(a, rv_pos, d_rv_pos, rv_tan, phase_factors):
Expand Down
2 changes: 1 addition & 1 deletion src/example_code.py
Original file line number Diff line number Diff line change
Expand Up @@ -264,7 +264,7 @@ def plot_compare_angle(theta1, ts1, label1, theta2, ts2, label2, title):
dwba_ts = mod.calculate_ts(m)

# and then a SDWBA version of the same
m |= {'phase_sd': 20*np.pi/180, 'num_runs': 100}
m |= {'phase_sd': 20, 'num_runs': 100}
sdwba_ts = mod.calculate_ts(m)

plt.plot(m['theta'], sdwba_ts, label='sdwba')
Expand Down

0 comments on commit ea396cf

Please sign in to comment.