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.
          """
          # extracting dictionary of coefficients specific to required
          # intensity measure type.
          C = self.COEFFS[imt]
          if isinstance(imt, PGA):
              imt_per = 0.0
          else:
              imt_per = imt.period
          # Fix site parameters for consistent dS2S application.
          sites.vs30 = np.array([250])
          sites.z1pt0 = np.array([330])
          # intensity on a reference soil is used for both mean
          # and stddev calculations.
          ln_y_ref = self._get_ln_y_ref(rup, dists, C)
          # exp1 and exp2 are parts of eq. 7
          exp1 = np.exp(C['phi3'] * (sites.vs30.clip(-np.inf, 1130) - 360))
          exp2 = np.exp(C['phi3'] * (1130 - 360))
          # v1 is the period dependent site term. The Vs30 above which, the
          # amplification is constant
          v1 = self._get_v1(imt)
          # Get log-mean from regular unadjusted model
          b13a_mean = self._get_mean(sites, C, ln_y_ref, exp1, exp2, v1)
          # Adjust mean and standard deviation
          mean = b13a_mean + self._get_dL2L(imt_per) + self._get_dS2S(imt_per)
          mean += convert_to_LHC(imt)
          stddevs = self._get_adjusted_stddevs(sites, rup, C, stddev_types,
                                               ln_y_ref, exp1, exp2, imt_per)
  
          return mean, stddevs