Python units.d方法代码示例

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def test_zeroFAP(self):
        """
        Test case where false alarm probability is 0 and missed det prob is 0
        All observations above SNR limit should be detected
        """

        SNRin = np.array([1.0,2.0,3.0,4.0,5.0,5.1,6.0])
        expected_MDresult = [True,True,True,True,False,False,False] 

        for mod in self.allmods:
            with RedirectStreams(stdout=self.dev_null):
                obj = mod(FAP=0.0,MDP=0.0,**self.specs)

            FA, MD = obj.det_occur(SNRin,self.TL.OpticalSystem.observingModes[0],self.TL,0,1*u.d)     
            for x,y in zip(MD,expected_MDresult):
                self.assertEqual( x,y )
            self.assertEqual( FA, False)

    #another limiting case 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def test_oneFAP(self):
        """
        Test case where false alarm probability is 1 and missed det prob is 0
        """

        SNRin = np.array([1.0,2.0,3.0,4.0,5.0,5.1,6.0])
        expected_MDresult = [True,True,True,True,False,False,False] 

        for mod in self.allmods:
            with RedirectStreams(stdout=self.dev_null):
                obj = mod(FAP=1.0,MDP=0.0,**self.specs)

            FA, MD = obj.det_occur(SNRin,self.TL.OpticalSystem.observingModes[0],self.TL,0,1*u.d)     
            for x,y in zip(MD,expected_MDresult):
                self.assertEqual( x,y )
            self.assertEqual( FA, True) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def test_advancetToStartOfNextOB(self):
        r""" Test advancetToStartOfNextOB method
                Strategy is to call the method once and ensure it advances the Observing Block
        """  
        tk = self.fixture()

        tk.OBstartTimes = [0,10,20,30]*u.d
        tk.OBendTimes = [5,15,25,35]*u.d
        tk.OBnumber = 0

        # 1) Set current Time to End of OB
        tk.currentTimeNorm = tk.OBendTimes[0]
        tk.currentTimeAbs = tk.missionStart + tk.currentTimeNorm
        tmpOBnumber = tk.OBnumber
        tmpcurrentTimeNorm = tk.currentTimeNorm.copy()
        tmpcurrentTimeAbs = tk.currentTimeAbs.copy()
        tmpexoplanetObsTime = tk.exoplanetObsTime.copy()
        tmpOBstartTimes = tk.OBstartTimes.copy()
        tmpOBendTimes = tk.OBendTimes.copy()
        tk.advancetToStartOfNextOB()
        self.assertTrue(tmpOBnumber + 1 == tk.OBnumber)
        self.assertTrue(tk.currentTimeNorm == tk.OBstartTimes[tk.OBnumber])
        self.assertTrue(tk.currentTimeAbs == tk.OBstartTimes[tk.OBnumber] + tk.missionStart) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def test_reset_sim(self):
        r"""Test reset_sim method.

        Approach: Ensures the simulation resets completely
        """

        with RedirectStreams(stdout=self.dev_null):
            sim = self.fixture(SimpleScript)
            sim.run_sim()
        
        self.assertGreater(len(sim.DRM), 0)
        self.assertGreater(sim.TimeKeeping.currentTimeNorm,0.0*u.d)

        sim.reset_sim()
        self.assertEqual(sim.TimeKeeping.currentTimeNorm,0.0*u.d)
        self.assertEqual(len(sim.DRM), 0) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def tests_comp_per_intTime(self):

        comp = self.fixture.comp_per_intTime(np.ones(self.nStars)*u.d,self,np.arange(self.nStars),0/(u.arcsec**2),0/(u.arcsec**2),0*u.arcsec,{})
        
        self.assertTrue(len(comp),self.nStars)

        with self.assertRaises(AssertionError):
            comp = self.fixture.comp_per_intTime(np.ones(self.nStars-1)*u.d,self,np.arange(self.nStars),0/(u.arcsec**2),0/(u.arcsec**2),0*u.arcsec,{})

        with self.assertRaises(AssertionError):
            comp = self.fixture.comp_per_intTime(np.ones(self.nStars)*u.d,self,np.arange(self.nStars),np.zeros(2)/(u.arcsec**2),0/(u.arcsec**2),0*u.arcsec,{})

        with self.assertRaises(AssertionError):
            comp = self.fixture.comp_per_intTime(np.ones(self.nStars)*u.d,self,np.arange(self.nStars),0/(u.arcsec**2),np.zeros(2)/(u.arcsec**2),0*u.arcsec,{})

        with self.assertRaises(AssertionError):
            comp = self.fixture.comp_per_intTime(np.ones(self.nStars-1)*u.d,self,np.arange(self.nStars),0/(u.arcsec**2),0/(u.arcsec**2),np.zeros(2)*u.arcsec,{}) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def revisitFilter(self, sInds, tmpCurrentTimeNorm):
        """Helper method for Overloading Revisit Filtering

        Args:
            sInds - indices of stars still in observation list
            tmpCurrentTimeNorm (MJD) - the simulation time after overhead was added in MJD form
        Returns:
            sInds - indices of stars still in observation list
        """
        tovisit = np.zeros(self.TargetList.nStars, dtype=bool)#tovisit is a boolean array containing the 
        if len(sInds) > 0:#so long as there is at least 1 star left in sInds
            tovisit[sInds] = ((self.starVisits[sInds] == min(self.starVisits[sInds])) \
                    & (self.starVisits[sInds] < self.nVisitsMax))# Checks that no star has exceeded the number of revisits
            if self.starRevisit.size != 0:#There is at least one revisit planned in starRevisit
                dt_rev = self.starRevisit[:,1]*u.day - tmpCurrentTimeNorm#absolute temporal spacing between revisit and now.

                #return indices of all revisits within a threshold dt_max of revisit day and indices of all revisits with no detections past the revisit time
                ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds)]
                tovisit[ind_rev2] = (self.starVisits[ind_rev2] < self.nVisitsMax)
            sInds = np.where(tovisit)[0]

        return sInds 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def __init__(self, coeffs=[1,1,2,1], revisit_wait=91.25, **specs):
        
        SurveySimulation.__init__(self, **specs)
        
        #verify that coefficients input is iterable 4x1
        if not(isinstance(coeffs,(list,tuple,np.ndarray))) or (len(coeffs) != 4):
            raise TypeError("coeffs must be a 4 element iterable")


        #Add to outspec
        self._outspec['coeffs'] = coeffs
        self._outspec['revisit_wait'] = revisit_wait
        
        # normalize coefficients
        coeffs = np.array(coeffs)
        coeffs = coeffs/np.linalg.norm(coeffs)
        
        self.coeffs = coeffs

        self.revisit_wait = revisit_wait*u.d
        self.no_dets = np.ones(self.TargetList.nStars, dtype=bool) 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def revisitFilter(self, sInds, tmpCurrentTimeNorm):
        """Helper method for Overloading Revisit Filtering

        Args:
            sInds - indices of stars still in observation list
            tmpCurrentTimeNorm (MJD) - the simulation time after overhead was added in MJD form
        Returns:
            sInds - indices of stars still in observation list
        """
        tovisit = np.zeros(self.TargetList.nStars, dtype=bool)#tovisit is a boolean array containing the 
        if len(sInds) > 0:#so long as there is at least 1 star left in sInds
            tovisit[sInds] = ((self.starVisits[sInds] == min(self.starVisits[sInds])) \
                    & (self.starVisits[sInds] < self.nVisitsMax))# Checks that no star has exceeded the number of revisits
            if self.starRevisit.size != 0:#There is at least one revisit planned in starRevisit
                dt_rev = self.starRevisit[:,1]*u.day - tmpCurrentTimeNorm#absolute temporal spacing between revisit and now.

                #return indices of all revisits within a threshold dt_max of revisit day and indices of all revisits with no detections past the revisit time
                # ind_rev = [int(x) for x in self.starRevisit[np.abs(dt_rev) < self.dt_max, 0] if (x in sInds and self.no_dets[int(x)] == False)]
                # ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds and self.no_dets[int(x)] == True)]
                # tovisit[ind_rev] = (self.starVisits[ind_rev] < self.nVisitsMax)#IF duplicates exist in ind_rev, the second occurence takes priority
                ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds)]
                tovisit[ind_rev2] = (self.starVisits[ind_rev2] < self.nVisitsMax)
            sInds = np.where(tovisit)[0]

        return sInds 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def revisitFilter(self, sInds, tmpCurrentTimeNorm):
        """Helper method for Overloading Revisit Filtering

        Args:
            sInds - indices of stars still in observation list
            tmpCurrentTimeNorm (MJD) - the simulation time after overhead was added in MJD form
        Returns:
            sInds - indices of stars still in observation list
        """
        tovisit = np.zeros(self.TargetList.nStars, dtype=bool)#tovisit is a boolean array containing the 
        if len(sInds) > 0:#so long as there is at least 1 star left in sInds
            tovisit[sInds] = ((self.starVisits[sInds] == min(self.starVisits[sInds])) \
                    & (self.starVisits[sInds] < self.nVisitsMax))# Checks that no star has exceeded the number of revisits
            if self.starRevisit.size != 0:#There is at least one revisit planned in starRevisit
                dt_rev = self.starRevisit[:,1]*u.day - tmpCurrentTimeNorm#absolute temporal spacing between revisit and now.

                ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds)]
                tovisit[ind_rev2] = (self.starVisits[ind_rev2] < self.nVisitsMax)
            sInds = np.where(tovisit)[0]
        return sInds 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def sacrificeStarCbyT(self, sInds, t_dets, fZ, fEZ, WA, overheadTime):
        """Sacrifice the worst performing CbyT star
        Args:
            sInds[nStars] - indicies of stars in the list
            t_dets[nStars] - time to observe each star (in days)
            fZ[nStars] - zodiacal light for each target
            fEZ - 0 
            WA - inner working angle of the instrument
            overheadTime - overheadTime added to each observation
        Return:
            sInds[nStars] - indicies of stars in the list
            t_dets[nStars] - time to observe each star (in days)
            sacrificedStarTime - time to distribute in days       
        """
        CbyT = self.Completeness.comp_per_intTime(t_dets*u.d, self.TargetList, sInds, self.valfZmin[sInds], fEZ, WA, self.mode, self.Cb[sInds], self.Csp[sInds])/t_dets#takes 5 seconds to do 1 time for all stars

        sacrificeIndex = np.argmin(CbyT)#finds index of star to sacrifice

        #Need index of sacrificed star by this point
        sacrificedStarTime = t_dets[sacrificeIndex] + overheadTime#saves time being sacrificed
        sInds = np.delete(sInds,sacrificeIndex)
        t_dets = np.delete(t_dets,sacrificeIndex)
        return sInds, t_dets, sacrificedStarTime 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def objfun(self,t,sInds,fZ):
        """
        Objective Function for SLSQP minimization. Purpose is to maximize summed completeness

        Args:
            t (ndarray):
                Integration times in days. NB: NOT an astropy quantity.
            sInds (ndarray):
                Target star indices (of same size as t)
            fZ (astropy Quantity):
                Surface brightness of local zodiacal light in units of 1/arcsec2
                Same size as t

        """
        good = t*u.d >= 0.1*u.s # inds that were not downselected by initial MIP

        comp = self.Completeness.comp_per_intTime(t[good]*u.d, self.TargetList, sInds[good], fZ[good], 
                self.ZodiacalLight.fEZ0, self.WAint[sInds][good], self.detmode)
        #self.vprint(-comp.sum()) # for verifying SLSQP output
        return -comp.sum() 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def objfun_deriv(self,t,sInds,fZ):
        """
        Jacobian of objective Function for SLSQP minimization. 

        Args:
            t (astropy Quantity):
                Integration times in days. NB: NOT an astropy quantity.
            sInds (ndarray):
                Target star indices (of same size as t)
            fZ (astropy Quantity):
                Surface brightness of local zodiacal light in units of 1/arcsec2
                Same size as t

        """
        good = t*u.d >= 0.1*u.s # inds that were not downselected by initial MIP

        tmp = self.Completeness.dcomp_dt(t[good]*u.d, self.TargetList, sInds[good], fZ[good], 
                self.ZodiacalLight.fEZ0, self.WAint[sInds][good], self.detmode).to("1/d").value

        jac = np.zeros(len(t))
        jac[good] = tmp
        return -jac 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def revisitFilter(self, sInds, tmpCurrentTimeNorm):
        """Helper method for Overloading Revisit Filtering
        Args:
            sInds - indices of stars still in observation list
            tmpCurrentTimeNorm (MJD) - the simulation time after overhead was added in MJD form
        Returns:
            sInds - indices of stars still in observation list
        """
        tovisit = np.zeros(self.TargetList.nStars, dtype=bool)#tovisit is a boolean array containing the 
        if len(sInds) > 0:#so long as there is at least 1 star left in sInds
            tovisit[sInds] = ((self.starVisits[sInds] == min(self.starVisits[sInds])) \
                    & (self.starVisits[sInds] < self.nVisitsMax))# Checks that no star has exceeded the number of revisits
            if self.starRevisit.size != 0:#There is at least one revisit planned in starRevisit
                dt_rev = self.starRevisit[:,1]*u.day - tmpCurrentTimeNorm#absolute temporal spacing between revisit and now.

                #return indices of all revisits within a threshold dt_max of revisit day and indices of all revisits with no detections past the revisit time
                # ind_rev = [int(x) for x in self.starRevisit[np.abs(dt_rev) < self.dt_max, 0] if (x in sInds and self.no_dets[int(x)] == False)]
                # ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds and self.no_dets[int(x)] == True)]
                # tovisit[ind_rev] = (self.starVisits[ind_rev] < self.nVisitsMax)#IF duplicates exist in ind_rev, the second occurence takes priority
                ind_rev2 = [int(x) for x in self.starRevisit[dt_rev < 0*u.d, 0] if (x in sInds)]
                tovisit[ind_rev2] = (self.starVisits[ind_rev2] < self.nVisitsMax)
            sInds = np.where(tovisit)[0]

        return sInds 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def semi_major_axis(P, Mtotal):
    """Semi-major axis

    Kepler's third law

    Args:
        P (float): Orbital period [days]
        Mtotal (float): Mass [Msun]

    Returns:
        float or array: semi-major axis in AU
    """

    # convert inputs to array so they work with units
    P = np.array(P)
    Mtotal = np.array(Mtotal)

    Mtotal = Mtotal*c.M_sun.value
    P = (P * u.d).to(u.second).value
    G = c.G.value
    
    a = ((P**2)*G*Mtotal/(4*(np.pi)**2))**(1/3.)
    a = a/c.au.value

    return a 

# 需要导入模块: from astropy import units [as 别名]
# 或者: from astropy.units import d [as 别名]
def test_timedelta_setitem():
    t = TimeDelta([1, 2, 3] * u.d, format='jd')

    t[0] = 0.5
    assert allclose_jd(t.value, [0.5, 2, 3])

    t[1:] = 4.5
    assert allclose_jd(t.value, [0.5, 4.5, 4.5])

    t[:] = 86400 * u.s
    assert allclose_jd(t.value, [1, 1, 1])

    t[1] = TimeDelta(2, format='jd')
    assert allclose_jd(t.value, [1, 2, 1])

    with pytest.raises(ValueError) as err:
        t[1] = 1 * u.m
    assert 'cannot convert value to a compatible TimeDelta' in str(err.value)