Returns if an object's 'ra' and 'decl' 
    is above the horizon at a specific latitude, 
    given the MC's right ascension.
175906
def isAboveHorizon(ra, decl, mcRA, lat):
    
    # This function checks if the equatorial distance from 
    # the object to the MC is within its diurnal semi-arc.
    
    dArc, _ = dnarcs(decl, lat)
    dist = abs(angle.closestdistance(mcRA, ra))
    return dist <= dArc/2.0 + 0.0003
1251204
def densm(alt, d0, xm, tz, mn3, zn3, tn3, tgn3, mn2, zn2, tn2, tgn2):
    '''

'''
    xs = [0.0]*10
    ys = [0.0]*10
    y2out = [0.0]*10
    rgas = 831.4
    #rgas = 831.44621    #maybe make this a global constant?
    densm_tmp=d0
    if (alt>zn2[0]): # pragma: no cover
        if(xm==0.0):
            return tz[0]
        else:
            return d0

    #/* STRATOSPHERE/MESOSPHERE TEMPERATURE */
    if (alt>zn2[mn2-1]):
        z=alt
    else:
        z=zn2[mn2-1];
    mn=mn2;
    z1=zn2[0];
    z2=zn2[mn-1];
    t1=tn2[0];
    t2=tn2[mn-1];
    zg = zeta(z, z1);
    zgdif = zeta(z2, z1);

    #/* set up spline nodes */
    for k in range(mn):
        xs[k]=zeta(zn2[k],z1)/zgdif;
        ys[k]=1.0 / tn2[k];
    yd1=-tgn2[0] / (t1*t1) * zgdif;
    yd2=-tgn2[1] / (t2*t2) * zgdif * (pow(((re[0]+z2)/(re[0]+z1)),2.0));

    #/* calculate spline coefficients */
    spline (xs, ys, mn, yd1, yd2, y2out);   #No need to change this
    x = zg/zgdif;
    y = [0.0]
    splint (xs, ys, y2out, mn, x, y);

    #/* temperature at altitude */
    tz[0] = 1.0 / y[0];
    if (xm!=0.0):
        #/* calaculate stratosphere / mesospehere density */
        glb = gsurf[0] / (pow((1.0 + z1/re[0]),2.0));
        gamm = xm * glb * zgdif / rgas;

        #/* Integrate temperature profile */
        yi = [0.0]
        splini(xs, ys, y2out, mn, x, yi);
        expl=gamm*yi[0];
        if (expl>50.0): # pragma: no cover
            expl=50.0

        #/* Density at altitude */
        densm_tmp = densm_tmp * (t1 / tz[0]) * exp(-expl);

    if (alt>zn3[0]):
        if (xm==0.0):
            return tz[0]
        else:
            return densm_tmp

    #/* troposhere / stratosphere temperature */
    z = alt;
    mn = mn3;
    z1=zn3[0];
    z2=zn3[mn-1];
    t1=tn3[0];
    t2=tn3[mn-1];
    zg=zeta(z,z1);
    zgdif=zeta(z2,z1);



    #/* set up spline nodes */
    for k in range(mn):
        xs[k] = zeta(zn3[k],z1) / zgdif;
        ys[k] = 1.0 / tn3[k];
    
    yd1=-tgn3[0] / (t1*t1) * zgdif;
    yd2=-tgn3[1] / (t2*t2) * zgdif * (pow(((re[0]+z2)/(re[0]+z1)),2.0));

    #/* calculate spline coefficients */
    spline (xs, ys, mn, yd1, yd2, y2out);
    x = zg/zgdif;
    y = [0.0]
    splint (xs, ys, y2out, mn, x, y);

    #/* temperature at altitude */
    tz[0] = 1.0 / y[0];
    if (xm!=0.0):
        #/* calaculate tropospheric / stratosphere density */
        glb = gsurf[0] / (pow((1.0 + z1/re[0]),2.0));
        gamm = xm * glb * zgdif / rgas;

        #/* Integrate temperature profile */
        yi = [0.0]
        splini(xs, ys, y2out, mn, x, yi);
        expl=gamm*yi[0];
        if (expl>50.0): # pragma: no cover
            expl=50.0;

        #/* Density at altitude */
        densm_tmp = densm_tmp * (t1 / tz[0]) * exp(-expl);
    
    if (xm==0.0):
        return tz[0];
    else:
        return densm_tmp;
334525
def check_heading(self, m):
        ''''''
        if 'GPS_RAW' in self.status.msgs:
            gps = self.status.msgs['GPS_RAW']
            if gps.v < 3:
                return
            diff = math.fabs(angle_diff(m.heading, gps.hdg))
        elif 'GPS_RAW_INT' in self.status.msgs:
            gps = self.status.msgs['GPS_RAW_INT']
            if gps.vel < 300:
                return
            diff = math.fabs(angle_diff(m.heading, gps.cog / 100.0))
        else:
            return
        self.report('heading', diff < 20, 'heading error %u' % diff)
1259062
def _calc(cls, **kwargs):
        
        zenith = 90.833333 # offical value

        jd = kwargs.get("jd", None)
        lat = kwargs.get("lat", None)
        lon = kwargs.get("lon", None)
        stage = kwargs.get("stage", None)
        if jd is None or stage is None or lat is None or lon is None:
            raise ValueError("Must supply an 'jd', 'lat, 'lon', and 'stage' parameter")

        if stage != SunCycles.RISING and stage != SunCycles.SETTING:
            raise ValueError("'stage' parameter must be %s or %s" % (SunCycles.RISING, SunCycles.SETTING))

        longhr = lon / 15.
        if stage == SunCycles.RISING:
            apx = jd + ( (6 - longhr) / 24 )
        elif stage == SunCycles.SETTING:
            apx = jd + ( (18 - longhr) / 24 )

        sun_mean_anom = ( 0.9856 * apx ) - 3.289 # sun's mean anomaly
        #sun's longitude
        sun_lon = sun_mean_anom + (1.916 * np.sin( np.radians(sun_mean_anom) )) \
                + (0.02 * np.sin( np.radians(2 * sun_mean_anom) )) + 282.634

        if sun_lon > 360:
            sun_lon = sun_lon - 360
        elif sun_lon < 0:
            sun_lon = sun_lon + 360

        right_ascension = np.degrees(np.arctan( 0.91764 * np.tan( np.radians(sun_lon) ) )) # sun's right ascension

        if right_ascension > 360:
            right_ascension = right_ascension - 360
        elif right_ascension < 0:
            right_ascension = right_ascension + 360

        # put sun's right ascension value in the same quadrant as the sun's
        # true longitude
        lQuad = 90. * np.floor(sun_lon / 90.)
        raQuad = 90. * np.floor(right_ascension / 90.)
        right_ascension = right_ascension + ( lQuad - raQuad)
        right_ascension = right_ascension / 15. # Convert to hours

        # Sun's declination
        sinDecl = 0.39782 * np.sin( np.radians(sun_lon) )
        cosDecl = np.cos( np.arcsin( sinDecl ) )

        # Sun's local hour angle
        cosHr = (np.cos( np.radians(zenith) ) - ( sinDecl * np.sin(np.radians(lat)) )) \
                / ( cosDecl * np.cos( np.radians(lat) ) )

        if cosHr > 1: # Sun doesnt rise on this loc on this date
            return -1, -1
        elif cosHr < -1: # Sun doesnt set on this location on this date
            return -1, -1
        elif stage == SunCycles.RISING: # Sunrise
            hr = 360 - np.degrees(np.arccos(cosHr))
        elif stage == SunCycles.SETTING:  # Sunset
            hr = np.degrees(np.arccos(cosHr))

        hr = hr / 15. # Convert angle to hours

        localTime = hr + right_ascension - ( 0.06571 * apx ) - 6.622# local meantime of rise/set

        UTtime = localTime - longhr # adjust to UTC

        if UTtime < 0:
            UTtime = UTtime + 24
        elif UTtime > 24:
            UTtime = UTtime - 24

        hour = np.floor(UTtime)

        minute = (UTtime - hour) * 60
        if minute == 60:
            hour = hour + 1
            minute = 0

        return hour, minute
467071
public static Geometry computeShadow(Geometry geometry, Geometry sunPosition, double height) {
        if (sunPosition != null) {
            if (sunPosition instanceof Point) {
                return computeShadow(geometry, sunPosition.getCoordinate().x, sunPosition.getCoordinate().y, height, true);
            }
            throw new IllegalArgumentException("The sun position must be stored in a point with \n"
                    + "x = sun azimuth in radians (direction along the horizon, measured from north to\n"
                    + "east and y = sun altitude above the horizon in radians.");
        }
        return null;
    }
1100779
def location(self, relative_alt=False):
        ''''''
        self.wait_gps_fix()
        # wait for another VFR_HUD, to ensure we have correct altitude
        self.recv_match(type='VFR_HUD', blocking=True)
        self.recv_match(type='GLOBAL_POSITION_INT', blocking=True)
        if relative_alt:
            alt = self.messages['GLOBAL_POSITION_INT'].relative_alt*0.001
        else:
            alt = self.messages['VFR_HUD'].alt
        return location(self.messages['GPS_RAW_INT'].lat*1.0e-7,
                        self.messages['GPS_RAW_INT'].lon*1.0e-7,
                        alt,
                        self.messages['VFR_HUD'].heading)
46669
def galactic_xyz(self):
        """"""
        vector = _GALACTIC.dot(self.position.au)
        return Distance(vector)
780751
def is_circumpolar(self, var):
        
        xcoord = self.get_x(var)
        return xcoord is not None and xcoord.ndim == 2
883937
def sunRelation(obj, sun):
    
    if obj.id == const.SUN:
        return None
    dist = abs(angle.closestdistance(sun.lon, obj.lon))
    if dist < 0.2833: return CAZIMI
    elif dist < 8.0: return COMBUST
    elif dist < 16.0: return UNDER_SUN
    else:
        return None
310434
def coord(self):
        """"""
        if self._coordinate is None:
            self._coordinate = SkyCoord(self.origin.ra, self.origin.dec + 45 * units.arcsec)
        return self._coordinate
62895
def right_ascension_at_time(self,t):
        ''''''
        α0 = self.right_ascension_at_epoch
        t0 = self.epoch
        ω = self.rotational_speed
        return (α0 + ω * (t - t0)) % (2*π)
686570
def latitude(self, value):
        """"""
        self._latitude = math.radians(float(value))
        assert -self.PI / 2 <= self._latitude <= self.PI / 2, \
            "latitude value should be between -90..90."
121947
def coneSearch(self, center, radius=3*u.arcmin, magnitudelimit=25):
        '''
        
        '''


        self.magnitudelimit = magnitudelimit

        # run the query
        self.speak('querying GALEX, centered on {} with radius {}'.format(center, radius, magnitudelimit))

        coordinatetosearch = '{0.ra.deg} {0.dec.deg}'.format(center)
        table = astroquery.mast.Catalogs.query_region(coordinates=center, radius=radius, catalog='GALEX')



        # the gaia DR2 epoch is 2015.5
        epoch = 2005#???

        # create skycoord objects
        self.coordinates = coord.SkyCoord(  ra=table['ra'].data*u.deg,
                                        dec=table['dec'].data*u.deg,
                                        obstime=Time(epoch, format='decimalyear'))

        self.magnitudes = dict(NUV=table['nuv_mag'].data, FUV=table['fuv_mag'].data)
        self.magnitude = self.magnitudes['NUV']
98201
def earth_accel_df(IMU,ATT):
    ''''''
    r = rotation_df(ATT)
    accel = Vector3(IMU.AccX, IMU.AccY, IMU.AccZ)
    return r * accel
81111
function(obj) {
    if (obj.normal) { return obj.intersects(this); }
    return (!this.isParallelTo(obj) && this.distanceFrom(obj) <= Sylvester.precision);
  }
88859
def get_b(self):
        

        try:

            return self.b.value

        except AttributeError:

            # Transform from L,B to R.A., Dec

            return self.sky_coord.transform_to('galactic').b.value
5517
def orientality(obj, sun):
    
    dist = angle.distance(sun.lon, obj.lon)
    return OCCIDENTAL if dist < 180 else ORIENTAL