Standard deviation on the magnitude for the WC1994 area relation.
182948
def get_std_dev_mag(self, rake):
        
        assert rake is None or -180 <= rake <= 180
        if rake is None:
            # their "All" case
            return 0.24
        elif (-45 <= rake <= 45) or (rake >= 135) or (rake <= -135):
            # strike slip
            return 0.23
        elif rake > 0:
            # thrust/reverse
            return 0.25
        else:
            # normal
            return 0.25
936233
def _compute_mean(self, C, rup, dists, sites, imt):
        
        mean = (self._get_magnitude_scaling_term(C, rup.mag) +
                self._get_distance_scaling_term(C, rup.mag, dists.rrup) +
                self._get_site_amplification_term(C, sites.vs30))
        # convert from cm/s**2 to g for SA and from m/s**2 to g for PGA (PGV
        # is already in cm/s) and also convert from base 10 to base e.
        if imt.name == "PGA":
            mean = np.log((10 ** mean) * ((2 * np.pi / 0.01) ** 2) *
                          1e-2 / g)
        elif imt.name == "SA":
            mean = np.log((10 ** mean) * ((2 * np.pi / imt.period) ** 2) *
                          1e-2 / g)
        else:
            mean = np.log(10 ** mean)

        return mean
149274
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        # Check that the GSIM supports the standard deviations requested
        assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
                   for stddev_type in stddev_types)
        coe = self.COEFFS[imt]
        mean = (self._get_magnitude_scaling(coe, rup.mag) +
                self._get_distance_scaling(coe, dists.rhypo) +
                self._get_azimuth_correction(coe, dists.azimuth))
        # Convert to g
        if imt.name in "SA PGA":
            mean = np.log(np.exp(mean) / (100.0 * g))
        # Compute std
        stddevs = self._compute_std(coe, stddev_types, dists.azimuth.shape)
        return mean, stddevs
860907
def _compute_magnitude_distance_term(self, C, rjb, mag):
        
        rval = np.sqrt(rjb ** 2. + C["h"] ** 2.)
        return (C["b4"] + C["b5"] * (mag - 4.5)) * np.log(rval)
517790
def get_magnitude_scaling_term(self, C, rup):
        
        if rup.mag <= self.CONSTANTS["m_c"]:
            return C["ccr"] * rup.mag
        else:
            return (C["ccr"] * self.CONSTANTS["m_c"]) +\
                (C["dcr"] * (rup.mag - self.CONSTANTS["m_c"]))
939576
def _get_magnitude_scaling(self, C, mag):
        
        if mag < self.CONSTANTS["mh"]:
            return C["e1"] + C["b1"] * (mag - self.CONSTANTS["mref"]) +\
                C["b2"] * ((mag - self.CONSTANTS["mref"]) ** 2.)
        else:
            d_m = self.CONSTANTS["mh"] - self.CONSTANTS["mref"]
            return C["e1"] + C["b3"] * (mag - self.CONSTANTS["mh"]) +\
                (C["b1"] * d_m) + C["b2"] * (d_m ** 2.)
669550
def _get_stddevs(self, C, stddev_types, num_sites):
        
        assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES
                   for stddev_type in stddev_types)

        if self.DEFINED_FOR_INTENSITY_MEASURE_COMPONENT == 'Random horizontal':
            # Using equation 8 of the paper,
            # corrected as indicated in the erratum
            Sr = np.sqrt(C['SlZ']**2 + (C['S3']/np.sqrt(2))**2)
        else:
            Sr = C['SlZ']

        stddevs = [np.zeros(num_sites) + Sr for _ in stddev_types]
        return stddevs
927101
def average_over_area(q, x, y):
    """"""
    area = (np.max(x) - np.min(x))*(np.max(y) - np.min(y))
    integral = np.trapz(np.trapz(q, y, axis=0), x)
    return integral/area
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
1063341
def mjd2gmst(mjd):
    """"""
    
    D2PI= 6.2831853071795864769252867665590057683943387987502
    DS2R= 7.2722052166430399038487115353692196393452995355905e-5
    
    tu = (mjd-51544.5)/36525.0;
    
    st = fmod(mjd,1.0) * D2PI + ( 24110.54841 + ( 8640184.812866 + ( 0.093104 - 6.2e-6 * tu ) * tu ) * tu ) * DS2R 
    
    w = fmod(st,D2PI)
    if w >= 0.0 :
        return w
    else:
        return w+D2PI
116696
def surface_area(self):
        
        m = self.lattice.matrix
        return np.linalg.norm(np.cross(m[0], m[1]))
1221958
def dim_dtau(self, pars):
        
        """
        self._set_parameters(pars)
        # term 1
        num1 = -self.m * (self.w ** self.c) * self.c\
            * (self.tau ** (self.c - 1)) * np.sin(self.ang)
        term1 = self.sigmai * num1 / self.denom

        # term 2
        num2a = -self.m * self.otc * np.sin(self.ang)
        num2b = 2 * (self.w ** 2.0) * self.c * (self.tau ** (self.c - 1)) *\
            np.cos(self.ang)
        num2c = 2 * self.c * (self.w ** (self.c * 2)) *\
            (self.tau ** (2 * self.c - 1))
        term2 = self.sigma0 * num2a * (num2b + num2c) / (self.denom ** 2)

        result = term1 + term2

        return result
576122
def _get_magnitude_term(self, C, mag):
        
        f_mag = C["c0"] + C["c1"] * mag
        if (mag > 4.5) and (mag <= 5.5):
            return f_mag + (C["c2"] * (mag - 4.5))
        elif (mag > 5.5) and (mag <= 6.5):
            return f_mag + (C["c2"] * (mag - 4.5)) + (C["c3"] * (mag - 5.5))
        elif mag > 6.5:
            return f_mag + (C["c2"] * (mag - 4.5)) + (C["c3"] * (mag - 5.5)) +\
                (C["c4"] * (mag - 6.5))
        else:
            return f_mag
277715
def get_magnitude_squared(self):
        """"""
        return self.x*self.x + self.y*self.y
695855
def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types):
        
        C = self.COEFFS[imt]
        mean = (np.log(self.get_magnitude_term(C, rup) +
                       self.get_distance_term(C, dists.rrup)) +
                self.get_site_amplification(C, sites))
        stddevs = self.get_stddevs(C, sites.vs30.shape, rup.mag, stddev_types)
        return mean, stddevs
1126258
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]

        mean = (self.get_magnitude_scaling(C, rup.mag) +
                self.get_distance_term(C, rup, dists.rjb, imt) +
                self.get_site_amplification(C, sites))
        # GMPE originally in cm/s/s - convert to g
        if imt.name in "PGA SA":
            mean -= np.log(100.0 * g)
        stddevs = self.get_stddevs(C, dists.rjb.shape, stddev_types, sites)
        if self.sigma_mu_epsilon:
            # Apply the epistemic uncertainty factor (sigma_mu) multiplied by
            # the number of standard deviations
            sigma_mu = self.get_sigma_mu_adjustment(C, imt, rup, dists)
            # Cap sigma_mu at 0.5 ln units
            sigma_mu[sigma_mu > 0.5] = 0.5
            # Sigma mu should not be less than the standard deviation of the
            # fault-to-fault variability
            sigma_mu[sigma_mu < C["tau_fault"]] = C["tau_fault"]
            mean += (self.sigma_mu_epsilon * sigma_mu)
        return mean, stddevs
409516
def mean_size_ISO(self, k, r):
        '''        
        
        '''
        p = k + r
        q = r
        return self.mean_size(p=p, q=q)