Returns the effective precession spin argument for the larger mass.
331279
def primary_xi(mass1, mass2, spin1x, spin1y, spin2x, spin2y):
    
    spinx = primary_spin(mass1, mass2, spin1x, spin2x)
    spiny = primary_spin(mass1, mass2, spin1y, spin2y)
    return chi_perp_from_spinx_spiny(spinx, spiny)
959573
def mass_3d(self, r, rho0, Ra, Rs):
        
        m_3d = 4 * np.pi * rho0 * Ra ** 2 * Rs ** 2 / (Rs ** 2 - Ra ** 2) * (Rs * np.arctan(r / Rs) - Ra * np.arctan(r / Ra))
        return m_3d
960641
def pointing_vector(self, s, time):
        
        t = time - self._t0
        longitude = self._longitude + self._dlongdt * t

        # define the basis vectors in the spin (primed) coordinates in terms of
        # the Roche coordinates.
        # ez' = s
        # ex' =  (ex - s(s.ex)) /|i - s(s.ex)|
        # ey' = s x ex'
        ex = np.array([1., 0., 0.])
        ezp = s
        exp = (ex - s*np.dot(s,ex))
        eyp = np.cross(s, exp)

        return np.sin(self._colat)*np.cos(longitude)*exp +\
                  np.sin(self._colat)*np.sin(longitude)*eyp +\
                  np.cos(self._colat)*ezp
745740
private static int spinsFor(Node pred, boolean haveData) {
        if (MP && pred != null) {
            if (pred.isData != haveData)      // phase change
                return FRONT_SPINS + CHAINED_SPINS;
            if (pred.isMatched())             // probably at front
                return FRONT_SPINS;
            if (pred.waiter == null)          // pred apparently spinning
                return CHAINED_SPINS;
        }
        return 0;
    }
1088788
def mass_2d(self, r, rho0, gamma):
        
        alpha = np.sqrt(np.pi) * special.gamma(1. / 2 * (-1 + gamma)) / special.gamma(gamma / 2.) * r ** (2 - gamma)/(3 - gamma) *np.pi * 2 * rho0
        mass_2d = alpha*r
        return mass_2d
21049
def d_hkl(self, miller_index: Vector3Like) -> float:
        

        gstar = self.reciprocal_lattice_crystallographic.metric_tensor
        hkl = np.array(miller_index)
        return 1 / ((dot(dot(hkl, gstar), hkl.T)) ** (1 / 2))
132513
def mass_3d(self, r, rho0, gamma):
        
        mass_3d = 4 * np.pi * rho0 /(-gamma + 3) * r ** (-gamma + 3)
        return mass_3d
990283
def moment_inertia(self):
        
        tensor = inertia.sphere_inertia(mass=self.volume,
                                        radius=self.primitive.radius)
        return tensor
883031
def L(self,*args,**kwargs):
        
        #Make sure you are not using physical coordinates
        old_physical= kwargs.get('use_physical',None)
        kwargs['use_physical']= False
        Omega= kwargs.pop('Omega',None)
        t= kwargs.pop('t',None)
        thiso= self(*args,**kwargs)
        if not len(thiso.shape) == 2: thiso= thiso.reshape((thiso.shape[0],1))
        if len(thiso[:,0]) < 3:
            raise AttributeError("'linearOrbit has no angular momentum")
        elif len(thiso[:,0]) == 3 or len(thiso[:,0]) == 4:
            if Omega is None:
                out= thiso[0,:]*thiso[2,:]
            else:
                out= thiso[0,:]*(thiso[2,:]-Omega*t*thiso[0,:])
        elif len(thiso[:,0]) == 5:
            raise AttributeError("You must track the azimuth to get the angular momentum of a 3D Orbit")
        else: #len(thiso[:,0]) == 6
            vx= self.vx(*args,**kwargs)
            vy= self.vy(*args,**kwargs)
            vz= self.vz(*args,**kwargs)
            x= self.x(*args,**kwargs)
            y= self.y(*args,**kwargs)
            z= self.z(*args,**kwargs)
            out= nu.zeros((len(x),3))
            out[:,0]= y*vz-z*vy
            out[:,1]= z*vx-x*vz
            out[:,2]= x*vy-y*vx
        if not old_physical is None:
            kwargs['use_physical']= old_physical
        else:
            kwargs.pop('use_physical')
        return out
938075
def spin_sy(self):
        """"""
        return conversions.secondary_spin(self.mass1, self.mass2, self.spin1y,
                                        self.spin2y)
267653
def spa_tmplt(**kwds):
    
    # Pull out the input arguments
    f_lower = kwds['f_lower']
    delta_f = kwds['delta_f']
    distance = kwds['distance']
    mass1 = kwds['mass1']
    mass2 = kwds['mass2']
    s1z = kwds['spin1z']
    s2z = kwds['spin2z']
    phase_order = int(kwds['phase_order'])
    #amplitude_order = int(kwds['amplitude_order'])
    spin_order = int(kwds['spin_order'])

    if 'out' in kwds:
        out = kwds['out']
    else:
        out = None

    amp_factor = spa_amplitude_factor(mass1=mass1, mass2=mass2) / distance

    lal_pars = lal.CreateDict()
    if phase_order != -1:
        lalsimulation.SimInspiralWaveformParamsInsertPNPhaseOrder(
            lal_pars, phase_order)

    if spin_order != -1:
        lalsimulation.SimInspiralWaveformParamsInsertPNSpinOrder(
            lal_pars, spin_order)

    #Calculate the PN terms
    phasing = lalsimulation.SimInspiralTaylorF2AlignedPhasing(
                                    float(mass1), float(mass2),
                                    float(s1z), float(s2z),
                                    lal_pars)

    pfaN = phasing.v[0]
    pfa2 = phasing.v[2] / pfaN
    pfa3 = phasing.v[3] / pfaN
    pfa4 = phasing.v[4] / pfaN
    pfa5 = phasing.v[5] / pfaN
    pfa6 = (phasing.v[6] - phasing.vlogv[6] * log(4)) / pfaN
    pfa7 = phasing.v[7] / pfaN

    pfl5 = phasing.vlogv[5] / pfaN
    pfl6 = phasing.vlogv[6] / pfaN

    piM = lal.PI * (mass1 + mass2) * lal.MTSUN_SI

    kmin = int(f_lower / float(delta_f))

    vISCO = 1. / sqrt(6.)
    fISCO = vISCO * vISCO * vISCO / piM
    kmax = int(fISCO / delta_f)
    f_max = ceilpow2(fISCO)
    n = int(f_max / delta_f) + 1

    if not out:
        htilde = FrequencySeries(zeros(n, dtype=numpy.complex64), delta_f=delta_f, copy=False)
    else:
        if type(out) is not Array:
            raise TypeError("Output must be an instance of Array")
        if len(out) < kmax:
            kmax = len(out)
        if out.dtype != complex64:
            raise TypeError("Output array is the wrong dtype")
        htilde = FrequencySeries(out, delta_f=delta_f, copy=False)

    spa_tmplt_engine(htilde[kmin:kmax], kmin, phase_order, delta_f, piM, pfaN,
                     pfa2, pfa3, pfa4, pfa5, pfl5,
                     pfa6, pfl6, pfa7, amp_factor)
    return htilde
979174
private void calculateRotationalAxis(Matrix rotation, Atom translation,
			double c) {
		// Calculate magnitude of rotationAxis components, but not signs
		double sum=0;
		double[] rotAx = new double[3];
		for(int i=0;i<3;i++) {
			rotAx[i] = Math.sqrt(rotation.get(i, i)-c);
			sum+=rotAx[i]*rotAx[i];
		}
		for(int i=0;i<3;i++) {
			rotAx[i] /= Math.sqrt(sum);
		}

		// Now determine signs
		double d0 = rotation.get(2,1)-rotation.get(1,2); //=2u[0]*sin(theta)
		double d1 = rotation.get(0,2)-rotation.get(2,0); //=2u[1]*sin(theta)
		double d2 = rotation.get(1,0)-rotation.get(0,1); //=2u[2]*sin(theta)

		double s12 = rotation.get(2,1)+rotation.get(1,2); //=2*u[1]*u[2]*(1-cos(theta))
		double s02 = rotation.get(0,2)+rotation.get(2,0); //=2*u[0]*u[2]*(1-cos(theta))
		double s01 = rotation.get(1,0)+rotation.get(0,1); //=2*u[0]*u[1]*(1-cos(theta))

		//Only use biggest d for the sign directly, for numerical stability around 180deg
		if( Math.abs(d0) < Math.abs(d1) ) { // not d0
			if( Math.abs(d1) < Math.abs(d2) ) { //d2
				if(d2>=0){ //u[2] positive
					if( s02 < 0 ) rotAx[0] = -rotAx[0];
					if( s12 < 0 ) rotAx[1] = -rotAx[1];
				} else { //u[2] negative
					rotAx[2] = -rotAx[2];
					if( s02 >= 0 ) rotAx[0] = -rotAx[0];
					if( s12 >= 0 ) rotAx[1] = -rotAx[1];
				}
			} else { //d1
				if(d1>=0) {//u[1] positive
					if( s01 < 0) rotAx[0] = -rotAx[0];
					if( s12 < 0) rotAx[2] = -rotAx[2];
				} else { //u[1] negative
					rotAx[1] = -rotAx[1];
					if( s01 >= 0) rotAx[0] = -rotAx[0];
					if( s12 >= 0) rotAx[2] = -rotAx[2];
				}
			}
		} else { // not d1
			if( Math.abs(d0) < Math.abs(d2) ) { //d2
				if(d2>=0){ //u[2] positive
					if( s02 < 0 ) rotAx[0] = -rotAx[0];
					if( s12 < 0 ) rotAx[1] = -rotAx[1];
				} else { //u[2] negative
					rotAx[2] = -rotAx[2];
					if( s02 >= 0 ) rotAx[0] = -rotAx[0];
					if( s12 >= 0 ) rotAx[1] = -rotAx[1];
				}
			} else { //d0
				if(d0>=0) { //u[0] positive
					if( s01 < 0 ) rotAx[1] = -rotAx[1];
					if( s02 < 0 ) rotAx[2] = -rotAx[2];
				} else { //u[0] negative
					rotAx[0] = -rotAx[0];
					if( s01 >= 0 ) rotAx[1] = -rotAx[1];
					if( s02 >= 0 ) rotAx[2] = -rotAx[2];
				}
			}
		}

		rotationAxis = new AtomImpl();
		rotationAxis.setCoords(rotAx);

		// Calculate screw = (rotationAxis dot translation)*u
		double dotProduct = Calc.scalarProduct(rotationAxis, translation);

		screwTranslation = Calc.scale(rotationAxis, dotProduct);
		otherTranslation = Calc.subtract(translation, screwTranslation);

		Atom hypot = Calc.vectorProduct(otherTranslation,rotationAxis);
		Calc.scaleEquals(hypot,.5/Math.tan(theta/2.0));

		// Calculate rotation axis position
		rotationPos = Calc.scaleAdd(.5,otherTranslation, hypot);
	}
1141379
def tau0_from_mass1_mass2(mass1, mass2, f_lower):
    
    """
    mtotal = mass1 + mass2
    eta = eta_from_mass1_mass2(mass1, mass2)
    return tau0_from_mtotal_eta(mtotal, eta, f_lower)
867241
def _dEndfr(self):
        
        Mc = self._chirp_mass()
        return (np.pi**(2./3.)*Mc**(5./3.)/(3.*(1.+self.z)**(1./3.)
                * (self.freqs_orb/(1.+self.z))**(1./3.))*(2./self.n)**(2./3.)
                * self._g_func()/self._f_func())
333383
def spin2_polar(self):
        """"""
        return coordinates.cartesian_to_spherical_polar(
                                     self.spin2x, self.spin2y, self.spin2z)
1184815
def imrphenomc_tmplt(**kwds):
    
    # Pull out the input arguments
    f_min = float128(kwds['f_lower'])
    f_max = float128(kwds['f_final'])
    delta_f = float128(kwds['delta_f'])
    distance = float128(kwds['distance'])
    mass1 = float128(kwds['mass1'])
    mass2 = float128(kwds['mass2'])
    spin1z = float128(kwds['spin1z'])
    spin2z = float128(kwds['spin2z'])

    if 'out' in kwds:
        out = kwds['out']
    else:
        out = None

    # Calculate binary parameters
    M = mass1 + mass2
    eta = mass1 * mass2 / (M * M)
    Xi = (mass1 * spin1z / M) + (mass2 * spin2z / M)
    Xisum = 2.*Xi
    Xiprod = Xi*Xi
    Xi2 = Xi*Xi

    m_sec = M * lal.MTSUN_SI;
    piM = lal.PI * m_sec;

    ## The units of distance given as input is taken to pe Mpc. Converting to SI
    distance *= (1.0e6 * lal.PC_SI / (2. * sqrt(5. / (64.*lal.PI)) * M * lal.MRSUN_SI * M * lal.MTSUN_SI))

    # Check if the value of f_max is correctly given, else replace with the fCut
    # used in the PhenomB code in lalsimulation. The various coefficients come
    # from Eq.(4.18) of http://arxiv.org/pdf/0710.2335 and
    # Table I of http://arxiv.org/pdf/0712.0343
    if not f_max:
        f_max = (1.7086 * eta * eta - 0.26592 * eta + 0.28236) / piM

    # Transform the eta, chi to Lambda parameters, using Eq 5.14, Table II of Main
    # paper.
    z101 = -2.417e-03
    z102 = -1.093e-03
    z111 = -1.917e-02
    z110 = 7.267e-02
    z120 = -2.504e-01

    z201 = 5.962e-01
    z202 = -5.600e-02
    z211 = 1.520e-01
    z210 = -2.970e+00
    z220 = 1.312e+01

    z301 = -3.283e+01
    z302 = 8.859e+00
    z311 = 2.931e+01
    z310 = 7.954e+01
    z320 = -4.349e+02

    z401 = 1.619e+02
    z402 = -4.702e+01
    z411 = -1.751e+02
    z410 = -3.225e+02
    z420 = 1.587e+03

    z501 = -6.320e+02
    z502 = 2.463e+02
    z511 = 1.048e+03
    z510 = 3.355e+02
    z520 = -5.115e+03

    z601 = -4.809e+01
    z602 = -3.643e+02
    z611 = -5.215e+02
    z610 = 1.870e+03
    z620 = 7.354e+02

    z701 = 4.149e+00
    z702 = -4.070e+00
    z711 = -8.752e+01
    z710 = -4.897e+01
    z720 = 6.665e+02

    z801 = -5.472e-02
    z802 = 2.094e-02
    z811 = 3.554e-01
    z810 = 1.151e-01
    z820 = 9.640e-01

    z901 = -1.235e+00
    z902 = 3.423e-01
    z911 = 6.062e+00
    z910 = 5.949e+00
    z920 = -1.069e+01

    eta2 = eta*eta
    Xi2 = Xiprod

    # Calculate alphas, gamma, deltas from Table II and Eq 5.14 of Main paper
    a1 = z101 * Xi + z102 * Xi2 + z111 * eta * Xi + z110 * eta + z120 * eta2
    a2 = z201 * Xi + z202 * Xi2 + z211 * eta * Xi + z210 * eta + z220 * eta2
    a3 = z301 * Xi + z302 * Xi2 + z311 * eta * Xi + z310 * eta + z320 * eta2
    a4 = z401 * Xi + z402 * Xi2 + z411 * eta * Xi + z410 * eta + z420 * eta2
    a5 = z501 * Xi + z502 * Xi2 + z511 * eta * Xi + z510 * eta + z520 * eta2
    a6 = z601 * Xi + z602 * Xi2 + z611 * eta * Xi + z610 * eta + z620 * eta2

    g1 = z701 * Xi + z702 * Xi2 + z711 * eta * Xi + z710 * eta + z720 * eta2

    del1 = z801 * Xi + z802 * Xi2 + z811 * eta * Xi + z810 * eta + z820 * eta2
    del2 = z901 * Xi + z902 * Xi2 + z911 * eta * Xi + z910 * eta + z920 * eta2

    # Get the spin of the final BH
    afin = FinalSpin( Xi, eta )
    Q = Qa( abs(afin) )

    # Get the fRD
    frd = fRD( abs(afin), M)
    Mfrd = frd * m_sec

    # Define the frequencies where SPA->PM->RD
    f1 = 0.1 * frd
    Mf1 = m_sec * f1
    f2 = frd
    Mf2 = m_sec * f2
    d1 = 0.005
    d2 = 0.005
    f0 = 0.98 * frd
    Mf0 = m_sec * f0
    d0 = 0.015

    # Now use this frequency for calculation of betas
    # calculate beta1 and beta2, that appear in Eq 5.7 in the main paper.
    b2 = ((-5./3.)* a1 * pow(Mfrd,(-8./3.)) - a2/(Mfrd*Mfrd) - \
      (a3/3.)*pow(Mfrd,(-4./3.)) + (2./3.)* a5 * pow(Mfrd,(-1./3.)) + a6)/eta

    psiPMrd = (a1 * pow(Mfrd,(-5./3.)) + a2/Mfrd + a3 * pow(Mfrd,(-1./3.)) + \
      a4 + a5 * pow(Mfrd,(2./3.)) + a6 * Mfrd)/eta
    b1 = psiPMrd - (b2 * Mfrd)

    ### Calculate the PN coefficients, Eq A3 - A5 of main paper ###
    pfaN = 3.0/(128.0 * eta)
    pfa2 = (3715./756.) + (55.*eta/9.0)
    pfa3 = -16.0*lal.PI + (113./3.)*Xi - 38.*eta*Xisum/3.
    pfa4 = (152.93365/5.08032) - 50.*Xi2 + eta*(271.45/5.04 + 1.25*Xiprod) + \
        3085.*eta2/72.
    pfa5 = lal.PI*(386.45/7.56 - 65.*eta/9.) - \
        Xi*(735.505/2.268 + 130.*eta/9.) + Xisum*(1285.0*eta/8.1 + 170.*eta2/9.) - \
        10.*Xi2*Xi/3. + 10.*eta*Xi*Xiprod
    pfa6 = 11583.231236531/4.694215680 - 640.0*lal.PI*lal.PI/3. - \
        6848.0*lal.GAMMA/21. - 684.8*log(64.)/6.3 + \
        eta*(2255.*lal.PI*lal.PI/12. - 15737.765635/3.048192) + \
        76.055*eta2/1.728 - (127.825*eta2*eta/1.296) + \
        2920.*lal.PI*Xi/3. - (175. - 1490.*eta)*Xi2/3. - \
        (1120.*lal.PI/3. - 1085.*Xi/3.)*eta*Xisum + \
        (269.45*eta/3.36 - 2365.*eta2/6.)*Xiprod

    pfa6log = -6848./63.

    pfa7 = lal.PI*(770.96675/2.54016 + 378.515*eta/1.512 - 740.45*eta2/7.56) - \
        Xi*(20373.952415/3.048192 + 1509.35*eta/2.24 - 5786.95*eta2/4.32) + \
        Xisum*(4862.041225*eta/1.524096 + 1189.775*eta2/1.008 - 717.05*eta2*eta/2.16 - 830.*eta*Xi2/3. + 35.*eta2*Xiprod/3.) - \
        560.*lal.PI*Xi2 + 20.*lal.PI*eta*Xiprod + \
        Xi2*Xi*(945.55/1.68 - 85.*eta) + Xi*Xiprod*(396.65*eta/1.68 + 255.*eta2)


    xdotaN = 64.*eta/5.
    xdota2 = -7.43/3.36 - 11.*eta/4.
    xdota3 = 4.*lal.PI - 11.3*Xi/1.2 + 19.*eta*Xisum/6.
    xdota4 = 3.4103/1.8144 + 5*Xi2 + eta*(13.661/2.016 - Xiprod/8.) + 5.9*eta2/1.8
    xdota5 = -lal.PI*(41.59/6.72 + 189.*eta/8.) - Xi*(31.571/1.008 - 116.5*eta/2.4) + \
          Xisum*(21.863*eta/1.008 - 79.*eta2/6.) - 3*Xi*Xi2/4. + \
          9.*eta*Xi*Xiprod/4.
    xdota6 = 164.47322263/1.39708800 - 17.12*lal.GAMMA/1.05 + \
          16.*lal.PI*lal.PI/3 - 8.56*log(16.)/1.05 + \
          eta*(45.1*lal.PI*lal.PI/4.8 - 561.98689/2.17728) + \
          5.41*eta2/8.96 - 5.605*eta*eta2/2.592 - 80.*lal.PI*Xi/3. + \
          eta*Xisum*(20.*lal.PI/3. - 113.5*Xi/3.6) + \
          Xi2*(64.153/1.008 - 45.7*eta/3.6) - \
          Xiprod*(7.87*eta/1.44 - 30.37*eta2/1.44)

    xdota6log = -856./105.

    xdota7 = -lal.PI*(4.415/4.032 - 358.675*eta/6.048 - 91.495*eta2/1.512) - \
          Xi*(252.9407/2.7216 - 845.827*eta/6.048 + 415.51*eta2/8.64) + \
          Xisum*(158.0239*eta/5.4432 - 451.597*eta2/6.048 + 20.45*eta2*eta/4.32 + 107.*eta*Xi2/6. - 5.*eta2*Xiprod/24.) + \
          12.*lal.PI*Xi2 - Xi2*Xi*(150.5/2.4 + eta/8.) + \
          Xi*Xiprod*(10.1*eta/2.4 + 3.*eta2/8.)


    AN = 8.*eta*sqrt(lal.PI/5.)
    A2 = (-107. + 55.*eta)/42.
    A3 = 2.*lal.PI - 4.*Xi/3. + 2.*eta*Xisum/3.
    A4 = -2.173/1.512 - eta*(10.69/2.16 - 2.*Xiprod) + 2.047*eta2/1.512
    A5 = -10.7*lal.PI/2.1 + eta*(3.4*lal.PI/2.1)

    A5imag = -24.*eta

    A6 = 270.27409/6.46800 - 8.56*lal.GAMMA/1.05 + \
      2.*lal.PI*lal.PI/3. + \
      eta*(4.1*lal.PI*lal.PI/9.6 - 27.8185/3.3264) - \
      20.261*eta2/2.772 + 11.4635*eta*eta2/9.9792 - \
      4.28*log(16.)/1.05

    A6log = -428./105.

    A6imag = 4.28*lal.PI/1.05

    ### Define other parameters needed by waveform generation ###
    kmin = int(f_min / delta_f)
    kmax = int(f_max / delta_f)
    n = kmax + 1;

    if not out:
        htilde = FrequencySeries(zeros(n,dtype=numpy.complex128), delta_f=delta_f, copy=False)
    else:
        if type(out) is not Array:
            raise TypeError("Output must be an instance of Array")
        if len(out) < kmax:
            raise TypeError("Output array is too small")
        if out.dtype != complex64:
            raise TypeError("Output array is the wrong dtype")
        htilde = FrequencySeries(out, delta_f=delta_f, copy=False)

    phenomC_kernel(htilde.data[kmin:kmax], kmin, delta_f, eta, Xi, distance,
                                       m_sec,  piM,  Mfrd,
                                       pfaN,  pfa2,  pfa3,  pfa4, pfa5,  pfa6,  pfa6log,  pfa7,
                                       a1,  a2,  a3,  a4, a5,  a6,  b1,  b2,
                                       Mf1,  Mf2,  Mf0, d1,  d2,  d0,
                                       xdota2,  xdota3,  xdota4, xdota5,  xdota6,  xdota6log,
                                       xdota7,  xdotaN,  AN, A2,  A3,  A4,  A5,
                                       A5imag,  A6,  A6log,  A6imag,
                                       g1,  del1,  del2,  Q )
    hp = htilde
    hc = htilde * 1j
    return hp, hc
522280
def secondary_spin(mass1, mass2, spin1, spin2):
    """"""
    mass1, mass2, spin1, spin2, input_is_array = ensurearray(
        mass1, mass2, spin1, spin2)
    ss = copy.copy(spin2)
    mask = mass1 < mass2
    ss[mask] = spin1[mask]
    return formatreturn(ss, input_is_array)