Q model

sourcespec/config_files/configspec.conf

@@ -546,6 +546,24 @@ t_star_min_max = float_list(default=list(0.001, 0.25))
 # If you want to search also negative t_star values:
 #   Qo_min_max = -10, 10
 Qo_min_max = float_list(default=None)
+# Prior knowledge of the attenuation model for your region can be specified
+# through Q_model as either a fixed value or a frequency-dependent function.
+# Use 'f' to represent frequency in the expression.
+# The frequency dependence is commonly expressed as:
+#   Q_model = Qo * (f/fo)**alpha
+# where:
+#   - Qo is the quality factor at the reference frequency fo
+#   - fo is the reference frequency (typically 1 Hz), to keep Q_model unitless
+#   - alpha is the frequency exponent
+# Any functional form can be used, provided Q_model remains unitless.
+# When Q_model is specified, t_star will not be inverted and all t_star-related
+# inversion parameters above will be ignored.
+# Examples:
+#   Fixed value:
+#       Q_model = 500
+#   Frequency-dependent (assuming fo=1 Hz):
+#       Q_model = 529 * f**0.42
+Q_model = string(default=None)
 # -------- INVERSION PARAMETERS
 
 # POST-INVERSION CONSTRAINTS --------

sourcespec/ssp_data_types.py

@@ -25,14 +25,18 @@ def __init__(self, Mw_0=None, fc_0=None, t_star_0=None):
 
     def __str__(self):
         """String representation."""
+        t_star_0 = (
+            'None' if self.t_star_0 is None else round(self.t_star_0, 4))
         return (
             f'Mw_0: {round(self.Mw_0, 4)}; '
             f'fc_0: {round(self.fc_0, 4)}; '
-            f't_star_0: {round(self.t_star_0, 4)}'
+            f't_star_0: {t_star_0}'
         )
 
     def get_params0(self):
         """Get initial values as a tuple."""
+        if self.t_star_0 is None:
+            return (self.Mw_0, self.fc_0)
         return (self.Mw_0, self.fc_0, self.t_star_0)
 
 
@@ -60,8 +64,11 @@ def __str__(self):
             *[round(x, 4) if x is not None else x for x in self.bounds[0]])
         s += 'fc: {}, {}; '.format(
             *[round(x, 4) if x is not None else x for x in self.bounds[1]])
-        s += 't_star: {}, {}'.format(
-            *[round(x, 4) if x is not None else x for x in self.bounds[2]])
+        try:
+            s += 't_star: {}, {}'.format(
+                *[round(x, 4) if x is not None else x for x in self.bounds[2]])
+        except IndexError:
+            s += 't_star: None, None'
         return s
 
     def _set_fc_min_max(self, config):
@@ -117,10 +124,8 @@ def _fix_initial_values_t_star(self):
     def __call__(self, **kwargs):
         """Interface for basin-hopping."""
         params = kwargs['x_new']
-        params_min = np.array(
-            (self.Mw_min, self.fc_min, self.t_star_min)).astype(float)
-        params_max = np.array(
-            (self.Mw_max, self.fc_max, self.t_star_max)).astype(float)
+        params_min = np.array([b[0] for b in self.bounds]).astype(float)
+        params_max = np.array([b[1] for b in self.bounds]).astype(float)
         params_min[np.isnan(params_min)] = 1e-99
         params_max[np.isnan(params_min)] = 1e+99
         tmin = bool(np.all(params >= params_min))
@@ -130,9 +135,13 @@ def __call__(self, **kwargs):
     @property
     def bounds(self):
         """Get bounds for minimize() as sequence of (min, max) pairs."""
-        self._bounds = ((self.Mw_min, self.Mw_max),
-                        (self.fc_min, self.fc_max),
-                        (self.t_star_min, self.t_star_max))
+        if self.t_star_min is None and self.t_star_max is None:
+            self._bounds = ((self.Mw_min, self.Mw_max),
+                            (self.fc_min, self.fc_max))
+        else:
+            self._bounds = ((self.Mw_min, self.Mw_max),
+                            (self.fc_min, self.fc_max),
+                            (self.t_star_min, self.t_star_max))
         return self._bounds
 
     @bounds.setter