Returns the mean ground motion in terms of log10 m/s/s, implementing
        equation 2 (page 502)
544395
def get_mean(self, C, rup, sites, dists):
        
        # W2 needs to be a 1 by 5 matrix (not a vector
        w_2 = np.array([
            [C["W_21"], C["W_22"], C["W_23"], C["W_24"], C["W_25"]]
            ])
        # Gets the style of faulting dummy variable
        sof = self._get_sof_dummy_variable(rup.rake)
        # Get the normalised coefficients
        p_n = self.get_pn(rup, sites, dists, sof)
        mean = np.zeros_like(dists.rjb)
        # Need to loop over sites - maybe this can be improved in future?
        # ndenumerate is used to allow for application to 2-D arrays
        for idx, rval in np.ndenumerate(p_n[0]):
            # Place normalised coefficients into a single array
            p_n_i = np.array([rval, p_n[1], p_n[2][idx], p_n[3], p_n[4]])
            # Executes the main ANN model
            mean_i = np.dot(w_2, np.tanh(np.dot(self.W_1, p_n_i) + self.B_1))
            mean[idx] = (0.5 * (mean_i + C["B_2"] + 1.0) *
                         (C["tmax"] - C["tmin"])) + C["tmin"]
        return mean
9367
def RMS_energy(frames):
    """"""
    f = frames.flatten()
    return N.sqrt(N.mean(f * f))
132122
def _get_stddev_deep_soil(self, mag, imt):
        
        # footnote from table 4 says that stderr for magnitudes over 7
        # is equal to one of magnitude 7.
        if mag > 7:
            mag = 7
        C = self.COEFFS_SOIL[imt]
        return C['sigma0'] + C['magfactor'] * mag
799303
def _compute_imt1100(self, imt, sites, rup, dists):
        
        vs30_1100 = np.zeros_like(sites.vs30) + 1100
        vs30_star, _ = self._compute_vs30_star_factor(imt, vs30_1100)
        C = self.COEFFS[imt]
        mean = (self._compute_base_term(C, rup, dists) +
                self._compute_faulting_style_term(C, rup) +
                self._compute_hanging_wall_term(C, dists, rup) +
                self._compute_top_of_rupture_depth_term(C, rup) +
                self._compute_large_distance_term(C, dists, rup) +
                self._compute_soil_depth_term(C, imt, sites.z1pt0, vs30_1100) +
                # this is the site response term in case of vs30=1100
                ((C['a10'] + C['b'] * self.CONSTS['n']) *
                np.log(vs30_star / C['VLIN'])))

        return mean
1271674
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]

        imean = (self._get_magnitude_term(C, rup.mag) +
                 self._get_distance_term(C, dists.rrup, sites.backarc) +
                 self._get_site_term(C, sites.vs30) +
                 self._get_scaling_term(C, dists.rrup))
        # Convert mean from cm/s and cm/s/s and from common logarithm to
        # natural logarithm
        if imt.name in "SA PGA":
            mean = np.log((10.0 ** (imean - 2.0)) / g)
        else:
            mean = np.log((10.0 ** (imean)))
        stddevs = self._get_stddevs(C, len(dists.rrup), stddev_types)
        return mean, stddevs
1094826
def get_step_f(step_f, lR2, lS2):
    
    mu, tau = 10, 2
    if lR2 > mu*lS2:
        return step_f * tau
    elif lS2 > mu*lR2:
        return step_f / tau
    return step_f
951229
def _get_sa_at_1180(self, C, imt, sites, rup, dists):
        
        # reference vs30 = 1180 m/s
        vs30_1180 = np.ones_like(sites.vs30) * 1180.
        # reference shaking intensity = 0
        ref_iml = np.zeros_like(sites.vs30)
        # fake Z1.0 - Since negative it will be replaced by the default Z1.0
        # for the corresponding region
        fake_z1pt0 = np.ones_like(sites.vs30) * -1
        return (self._get_basic_term(C, rup, dists) +
                self._get_faulting_style_term(C, rup) +
                self._get_site_response_term(C, imt, vs30_1180, ref_iml) +
                self._get_hanging_wall_term(C, dists, rup) +
                self._get_top_of_rupture_depth_term(C, imt, rup) +
                self._get_soil_depth_term(C, fake_z1pt0, vs30_1180) +
                self._get_regional_term(C, imt, vs30_1180, dists.rrup)
                )
254440
def _get_anelastic_attenuation_term(self, C, rrup):
        
        f_atn = np.zeros(len(rrup))
        idx = rrup >= 80.0
        f_atn[idx] = (C["c20"] + C["Dc20"]) * (rrup[idx] - 80.0)
        return f_atn
696295
def mean_velocity(self, tol=0.05, Lx=None):
        '''
        
        '''
        dx = None
        dt = None
        p = None
        ind = None
        t_end = None

        if self.area == 0:
            return dict(dx=dx, dt=dt, p=p, ind=ind, t_end=t_end)

        x = self.x_position(Lx=Lx)

        # find the first and last valid indices
        ind_start = mlab.find(x.mask==False)[0]
        ind_last = mlab.find(x.mask==False)[-1]

        # if the original x value is within tol % of the final x value, include
        # all samples
        if x[ind_start] > (1 - tol) * x[ind_last] or x[ind_last] < 0:
            ind_stop = ind_last
        else: # otherwise, stop when x reaches (1 - tol) % of it's final value
            ind_stop = mlab.find(x > (1 - tol) * x[ind_last])[0]

        ind = [ind_start, ind_stop]

        # if we have at least 2 valid samples
        if len(ind) >=2:
            dx = np.diff(x[ind])[0]
            dt = np.diff(self.time[ind])[0] # ms

            # suppress polyfit warnings
            with warnings.catch_warnings():
                warnings.simplefilter("ignore")
                # fit a line to the data
                p = np.polyfit(self.time[ind[0]:ind[1]], x[ind[0]:ind[1]], 1)

            # find time when the the line intercepts x[ind_last]
            ind_stop = mlab.find(self.time > \
                                 (x[ind_last] - p[1]) / p[0])
            if len(ind_stop):
                t_end = self.time[ind_stop[0]]
            else:
                t_end = self.time[-1]
        return dict(dx=dx, dt=dt, p=p, ind=ind, t_end=t_end)
318642
def _compute_mean(self, C, S, mag, rrup, rvol, hypo_depth, CN, CR, f4HW,
                      delta_C, delta_D):
        

        # Stage 1: compute PGA_ABCD and PGA'_ABCD which are then used in
        # equation 6
        # Equation 1 PGA unprimed version
        lnSA_AB = self._compute_mean_on_rock(C, mag, rrup, rvol, hypo_depth,
                                             CN, CR, f4HW)

        # Equation 4 PGA unprimed version
        lnSA_ABCD = lnSA_AB + S *\
            self._compute_nonlinear_soil_term(C, lnSA_AB, delta_C, delta_D)

        return lnSA_ABCD
1365366
def _compute_site_term(self, C, vs30):
        
        # for rock values the site term is zero
        site_term = np.zeros_like(vs30)

        # hard soil
        site_term[(vs30 >= 360) & (vs30 < 800)] = C['aB']

        # medium soil
        site_term[(vs30 >= 180) & (vs30 < 360)] = C['aC']

        # soft soil
        site_term[vs30 < 180] = C['aD']

        return site_term
1343512
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]
        mean = (
            self._get_magnitude_term(C, rup.mag) +
            self._get_distance_term(C, rup.mag, dists.rrup) +
            self._get_style_of_faulting_term(C, rup.rake) +
            self._get_site_response_term(C, sites.vs30))

        stddevs = self._get_stddevs(C, stddev_types, len(sites.vs30))
        return mean, stddevs
252611
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
                   for stddev_type in stddev_types)

        # GMPE differentiates strike-slip, reverse and normal ruptures,
        # but combines normal and strike-slip into one category. See page 180.
        is_reverse = (45 <= rup.rake <= 135)

        stddevs = [numpy.zeros_like(sites.vs30) for _ in stddev_types]
        means = numpy.zeros_like(sites.vs30)

        [rocks_i] = (sites.vs30 > self.ROCK_VS30).nonzero()
        if len(rocks_i):
            rrup = dists.rrup.take(rocks_i)
            mean_rock = self._get_mean_rock(rup.mag, rup.rake, rrup,
                                            is_reverse, imt)
            means.put(rocks_i, mean_rock)
            for stddev_arr in stddevs:
                stddev_rock = self._get_stddev_rock(rup.mag, imt)
                stddev_arr.put(rocks_i, stddev_rock)

        [soils_i] = (sites.vs30 <= self.ROCK_VS30).nonzero()
        if len(soils_i):
            rrup = dists.rrup.take(soils_i)
            mean_soil = self._get_mean_deep_soil(rup.mag, rup.rake, rrup,
                                                 is_reverse, imt)
            means.put(soils_i, mean_soil)
            for stddev_arr in stddevs:
                stddev_soil = self._get_stddev_deep_soil(rup.mag, imt)
                stddev_arr.put(soils_i, stddev_soil)

        return means, stddevs
425888
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        # extracting dictionary of coefficients specific to required
        # intensity measure type.
        C = self.COEFFS[imt]
        C_PGA = self.COEFFS[PGA()]
        imt_per = 0 if imt.name == 'PGV' else imt.period
        pga_rock = self._get_pga_on_rock(C_PGA, rup, dists)
        mean = (self._get_magnitude_scaling_term(C, rup) +
                self._get_path_scaling(C, dists, rup.mag) +
                self._get_site_scaling(C, pga_rock, sites, imt_per, dists.rjb))
        stddevs = self._get_stddevs(C, rup, dists, sites, stddev_types)
        return mean, stddevs
1140865
def _compute_mean_on_rock(self, C, mag, rrup, rvol, hypo_depth, CN, CR,
                              f4HW):
        

        # Define subduction flag (page 23)
        # SI=1 for subduction interface, 0 otherwise
        # DS=1 for subduction intraslab, 0 otherwise
        SI = 0
        DS = 1

        lnSA_AB = (
            # line 1 and 2 of equation 2
            C['c11'] + (C['c12y'] + (C['c15'] - C['c17']) * C['c19y']) *
            (mag - 6) +
            # line 3
            C['c13y'] * (10 - mag) ** 3 +
            # line 4
            C['c17'] * np.log(rrup + C['c18y'] * np.exp(C['c19y'] * mag)) +
            # line 5
            C['c20'] * hypo_depth + C['c24'] * SI +
            # line 6
            C['c46'] * rvol * (1 - DS)
        )

        return lnSA_AB
1146156
def _compute_mean_on_rock(self, C, mag, rrup, rvol, hypo_depth, CN, CR,
                              f4HW):
        

        lnSA_AB = (
            # line 1 of equation 1
            C['c1'] + C['c4as'] * (mag - 6) +
            # line 2
            C['c3as'] * (8.5 - mag) ** 2 +
            # line 3
            C['c5'] * rrup +
            # line 3 and 4
            (C['c8'] + C['c6as'] * (mag - 6)) *
            np.log((rrup ** 2 + C['c10as'] ** 2) ** 0.5) +
            # line 5
            C['c46'] * rvol +
            # line 6
            C['c32'] * CN + C['c33as'] * CR + f4HW
        )

        return lnSA_AB
656512
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]

        # clip rhypo at 10 (this is the minimum distance used in
        # deriving the equation), see page 22, this avoids singularity
        # in mean value equation
        rhypo = dists.rhypo.copy()
        rhypo[rhypo < 10] = 10

        # convert magnitude from Mblg to Mw
        mag = rup.mag * 0.98 - 0.39 if rup.mag <= 5.5 else \
              2.715 - 0.277 * rup.mag + 0.127 * rup.mag * rup.mag

        # functional form as explained in 'Youngs_fit_to_AB95lookup.doc'
        f1 = np.minimum(np.log(rhypo), np.log(70.))
        f2 = np.maximum(np.log(rhypo / 130.), 0)
        mean = (
            C['c1'] + C['c2'] * mag + C['c3'] * mag ** 2 +
            (C['c4'] + C['c5'] * mag) * f1 +
            (C['c6'] + C['c7'] * mag) * f2 +
            C['c8'] * rhypo
        )

        stddevs = self._get_stddevs(stddev_types,  dists.rhypo.shape[0])

        return mean, stddevs