docs: def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
          """
          See :meth:`superclass method
          <.base.GroundShakingIntensityModel.get_mean_and_stddevs>`
          for spec of input and result values.
          """
  
          # assign constant
          log10e = np.log10(np.e)
  
          # Distance term
          R = np.sqrt(dists.rjb ** 2 + 11.29 ** 2)
          # Magnitude term
          M = rup.mag - 6
  
          # Site term only distinguishes between lava and ash;
          # since ash sites have Vs30 in the range 60-200m/s,
          # we use this upper value as class separator
          S = np.zeros(R.shape)
          S[sites.vs30 <= 200] = 1
  
          # Mean ground motion (natural log)
          # call super
          mean, stddevs = super().get_mean_and_stddevs(sites, rup, dists,
                                                       imt, stddev_types)
  
          if rup.mag > 7. and rup.mag <= 7.7:
              mean = (0.171 * (1 - M)) / log10e + mean
  
          elif rup.mag > 7.7:
              mean = (0.1512 + 0.387 * (1 - M)) / log10e + mean
  
          # define natural log of SA 0.3 sec and 0.2 sec
          if isinstance(imt, SA):
              if imt.period == 0.3:
                  mean = np.log(2.2) + mean
  
              if imt.period == 0.2:
                  mean = np.log(2.5) + mean
  
          return mean, stddevs