signature:
def hms_to_decimal ( hours , minutes , seconds , returndeg = True ) :
    '''Converts from HH , MM , SS to a decimal value .
 Parameters
 hours : int
 The HH part of a RA coordinate .
 minutes : int
 The MM part of a RA coordinate .
 seconds : float
 The SS . sss part of a RA coordinate .
 returndeg : bool
 If this is True , then will return decimal degrees as the output .
 If this is False , then will return decimal HOURS as the output .
 Decimal hours are sometimes used in FITS headers .
 Returns
 float
 The right ascension value in either decimal degrees or decimal hours
 depending on ` returndeg ` .'''

implementation:
if hours > 24 :
    return None
else :
    dec_hours = fabs ( hours ) + fabs ( minutes ) / 60.0 + fabs ( seconds ) / 3600.0
    if returndeg :
        dec_deg = dec_hours * 15.0
        if dec_deg < 0 :
            dec_deg = dec_deg + 360.0
        dec_deg = dec_deg % 360.0
        return dec_deg
    else :
        return dec_hours


signature:
def dms_to_decimal ( sign , degrees , minutes , seconds ) :
    '''Converts from DD : MM : SS to a decimal value .
 Parameters
 sign : { ' + ' , ' - ' , ' ' }
 The sign part of a Dec coordinate .
 degrees : int
 The DD part of a Dec coordinate .
 minutes : int
 The MM part of a Dec coordinate .
 seconds : float
 The SS . sss part of a Dec coordinate .
 Returns
 float
 The declination value in decimal degrees .'''

implementation:
dec_deg = fabs ( degrees ) + fabs ( minutes ) / 60.0 + fabs ( seconds ) / 3600.0
if sign == '-' :
    return - dec_deg
else :
    return dec_deg


signature:
def great_circle_dist ( ra1 , dec1 , ra2 , dec2 ) :
    '''Calculates the great circle angular distance between two coords .
 This calculates the great circle angular distance in arcseconds between two
 coordinates ( ra1 , dec1 ) and ( ra2 , dec2 ) . This is basically a clone of GCIRC
 from the IDL Astrolib .
 Parameters
 ra1 , dec1 : float or array - like
 The first coordinate ' s right ascension and declination value ( s ) in
 decimal degrees .
 ra2 , dec2 : float or array - like
 The second coordinate ' s right ascension and declination value ( s ) in
 decimal degrees .
 Returns
 float or array - like
 Great circle distance between the two coordinates in arseconds .
 Notes
 If ( ` ra1 ` , ` dec1 ` ) is scalar and ( ` ra2 ` , ` dec2 ` ) is scalar : the result is a
 float distance in arcseconds .
 If ( ` ra1 ` , ` dec1 ` ) is scalar and ( ` ra2 ` , ` dec2 ` ) is array - like : the result
 is an np . array with distance in arcseconds between ( ` ra1 ` , ` dec1 ` ) and each
 element of ( ` ra2 ` , ` dec2 ` ) .
 If ( ` ra1 ` , ` dec1 ` ) is array - like and ( ` ra2 ` , ` dec2 ` ) is scalar : the result
 is an np . array with distance in arcseconds between ( ` ra2 ` , ` dec2 ` ) and each
 element of ( ` ra1 ` , ` dec1 ` ) .
 If ( ` ra1 ` , ` dec1 ` ) and ( ` ra2 ` , ` dec2 ` ) are both array - like : the result is an
 np . array with the pair - wise distance in arcseconds between each element of
 the two coordinate lists . In this case , if the input array - likes are not the
 same length , then excess elements of the longer one will be ignored .'''

implementation:
# wrap RA if negative or larger than 360.0 deg
in_ra1 = ra1 % 360.0
in_ra1 = in_ra1 + 360.0 * ( in_ra1 < 0.0 )
in_ra2 = ra2 % 360.0
in_ra2 = in_ra2 + 360.0 * ( in_ra1 < 0.0 )
# convert to radians
ra1_rad , dec1_rad = np . deg2rad ( in_ra1 ) , np . deg2rad ( dec1 )
ra2_rad , dec2_rad = np . deg2rad ( in_ra2 ) , np . deg2rad ( dec2 )
del_dec2 = ( dec2_rad - dec1_rad ) / 2.0
del_ra2 = ( ra2_rad - ra1_rad ) / 2.0
sin_dist = np . sqrt ( np . sin ( del_dec2 ) * np . sin ( del_dec2 ) + np . cos ( dec1_rad ) * np . cos ( dec2_rad ) * np . sin ( del_ra2 ) * np . sin ( del_ra2 ) )
dist_rad = 2.0 * np . arcsin ( sin_dist )
# return the distance in arcseconds
return np . rad2deg ( dist_rad ) * 3600.0


signature:
def xmatch_basic ( ra1 , dec1 , ra2 , dec2 , match_radius = 5.0 ) :
    '''Finds the closest object in ( ` ra2 ` , ` dec2 ` ) to scalar coordinate pair
 ( ` ra1 ` , ` dec1 ` ) and returns the distance in arcseconds .
 This is a quick matcher that uses the ` great _ circle _ dist ` function to find
 the closest object in ( ` ra2 ` , ` dec2 ` ) within ` match _ radius ` arcseconds to
 ( ` ra1 ` , ` dec1 ` ) . ( ` ra1 ` , ` dec1 ` ) must be a scalar pair , while
 ( ` ra2 ` , ` dec2 ` ) must be array - likes of the same lengths .
 Parameters
 ra1 , dec1 : float
 Coordinate of the object to find matches to . In decimal degrees .
 ra2 , dec2 : array - like
 The coordinates that will be searched for matches . In decimal degrees .
 match _ radius : float
 The match radius in arcseconds to use for the match .
 Returns
 tuple
 A two element tuple like the following : :
 ( True - > no match found or False - > found a match ,
 minimum distance between target and list in arcseconds )'''

implementation:
min_dist_arcsec = np . min ( great_circle_dist ( ra1 , dec1 , ra2 , dec2 ) )
if ( min_dist_arcsec < match_radius ) :
    return ( True , min_dist_arcsec )
else :
    return ( False , min_dist_arcsec )


signature:
def xmatch_neighbors ( ra1 , dec1 , ra2 , dec2 , match_radius = 60.0 , includeself = False , sortresults = True ) :
    '''Finds the closest objects in ( ` ra2 ` , ` dec2 ` ) to scalar coordinate pair
 ( ` ra1 ` , ` dec1 ` ) and returns the indices of the objects that match .
 This is a quick matcher that uses the ` great _ circle _ dist ` function to find
 the closest object in ( ` ra2 ` , ` dec2 ` ) within ` match _ radius ` arcseconds to
 ( ` ra1 ` , ` dec1 ` ) . ( ` ra1 ` , ` dec1 ` ) must be a scalar pair , while
 ( ` ra2 ` , ` dec2 ` ) must be array - likes of the same lengths .
 Parameters
 ra1 , dec1 : float
 Coordinate of the object to find matches to . In decimal degrees .
 ra2 , dec2 : array - like
 The coordinates that will be searched for matches . In decimal degrees .
 match _ radius : float
 The match radius in arcseconds to use for the match .
 includeself : bool
 If this is True , the object itself will be included in the match
 results .
 sortresults : bool
 If this is True , the match indices will be sorted by distance .
 Returns
 tuple
 A tuple like the following is returned : :
 ( True - > matches found or False - > no matches found ,
 minimum distance between target and list ,
 np . array of indices where list of coordinates is
 closer than ` match _ radius ` arcseconds from the target ,
 np . array of distances in arcseconds )'''

implementation:
dist = great_circle_dist ( ra1 , dec1 , ra2 , dec2 )
if includeself :
    match_dist_ind = np . where ( dist < match_radius )
else : # make sure we match only objects that are not the same as this object
    match_dist_ind = np . where ( ( dist < match_radius ) & ( dist > 0.1 ) )
if len ( match_dist_ind ) > 0 :
    match_dists = dist [ match_dist_ind ]
    dist_sort_ind = np . argsort ( match_dists )
    if sortresults :
        match_dist_ind = ( match_dist_ind [ 0 ] ) [ dist_sort_ind ]
    min_dist = np . min ( match_dists )
    return ( True , min_dist , match_dist_ind , match_dists [ dist_sort_ind ] )
else :
    return ( False , )


signature:
def make_kdtree ( ra , decl ) :
    '''This makes a ` scipy . spatial . CKDTree ` on ( ` ra ` , ` decl ` ) .
 Parameters
 ra , decl : array - like
 The right ascension and declination coordinate pairs in decimal degrees .
 Returns
 ` scipy . spatial . CKDTree `
 The cKDTRee object generated by this function is returned and can be
 used to run various spatial queries .'''

implementation:
# get the xyz unit vectors from ra , decl
# since i had to remind myself :
# https : / / en . wikipedia . org / wiki / Equatorial _ coordinate _ system
cosdecl = np . cos ( np . radians ( decl ) )
sindecl = np . sin ( np . radians ( decl ) )
cosra = np . cos ( np . radians ( ra ) )
sinra = np . sin ( np . radians ( ra ) )
xyz = np . column_stack ( ( cosra * cosdecl , sinra * cosdecl , sindecl ) )
# generate the kdtree
kdt = sps . cKDTree ( xyz , copy_data = True )
return kdt


signature:
def conesearch_kdtree ( kdtree , racenter , declcenter , searchradiusdeg , conesearchworkers = 1 ) :
    '''This does a cone - search around ( ` racenter ` , ` declcenter ` ) in ` kdtree ` .
 Parameters
 kdtree : scipy . spatial . CKDTree
 This is a kdtree object generated by the ` make _ kdtree ` function .
 racenter , declcenter : float or array - like
 This is the center coordinate to run the cone - search around in decimal
 degrees . If this is an np . array , will search for all coordinate pairs in
 the array .
 searchradiusdeg : float
 The search radius to use for the cone - search in decimal degrees .
 conesearchworkers : int
 The number of parallel workers to launch for the cone - search .
 Returns
 list or np . array of lists
 If ( ` racenter ` , ` declcenter ` ) is a single coordinate , this will return a
 list of the indices of the matching objects in the kdtree . If
 ( ` racenter ` , ` declcenter ` ) are array - likes , this will return an object
 array containing lists of matching object indices for each coordinate
 searched .'''

implementation:
cosdecl = np . cos ( np . radians ( declcenter ) )
sindecl = np . sin ( np . radians ( declcenter ) )
cosra = np . cos ( np . radians ( racenter ) )
sinra = np . sin ( np . radians ( racenter ) )
# this is the search distance in xyz unit vectors
xyzdist = 2.0 * np . sin ( np . radians ( searchradiusdeg ) / 2.0 )
# look up the coordinates
kdtindices = kdtree . query_ball_point ( [ cosra * cosdecl , sinra * cosdecl , sindecl ] , xyzdist , n_jobs = conesearchworkers )
return kdtindices


signature:
def xmatch_kdtree ( kdtree , extra , extdecl , xmatchdistdeg , closestonly = True ) :
    '''This cross - matches between ` kdtree ` and ( ` extra ` , ` extdecl ` ) arrays .
 Returns the indices of the kdtree and the indices of extra , extdecl that
 xmatch successfully .
 Parameters
 kdtree : scipy . spatial . CKDTree
 This is a kdtree object generated by the ` make _ kdtree ` function .
 extra , extdecl : array - like
 These are np . arrays of ' external ' coordinates in decimal degrees that
 will be cross - matched against the objects in ` kdtree ` .
 xmatchdistdeg : float
 The match radius to use for the cross - match in decimal degrees .
 closestonly : bool
 If closestonly is True , then this function returns only the closest
 matching indices in ( extra , extdecl ) for each object in kdtree if there
 are any matches . Otherwise , it returns a list of indices in ( extra ,
 extdecl ) for all matches within xmatchdistdeg between kdtree and ( extra ,
 extdecl ) .
 Returns
 tuple of lists
 Returns a tuple of the form : :
 ( list of ` kdtree ` indices matching to external objects ,
 list of all ` extra ` / ` extdecl ` indices that match to each
 element in ` kdtree ` within the specified cross - match distance )'''

implementation:
ext_cosdecl = np . cos ( np . radians ( extdecl ) )
ext_sindecl = np . sin ( np . radians ( extdecl ) )
ext_cosra = np . cos ( np . radians ( extra ) )
ext_sinra = np . sin ( np . radians ( extra ) )
ext_xyz = np . column_stack ( ( ext_cosra * ext_cosdecl , ext_sinra * ext_cosdecl , ext_sindecl ) )
ext_xyzdist = 2.0 * np . sin ( np . radians ( xmatchdistdeg ) / 2.0 )
# get our kdtree
our_kdt = kdtree
# get the external kdtree
ext_kdt = sps . cKDTree ( ext_xyz )
# do a query _ ball _ tree
extkd_matchinds = our_kdt . query_ball_tree ( ext_kdt , ext_xyzdist )
ext_matchinds = [ ]
kdt_matchinds = [ ]
for extind , mind in enumerate ( extkd_matchinds ) :
    if len ( mind ) > 0 : # our object indices
        kdt_matchinds . append ( extind )
        # external object indices
        if closestonly :
            ext_matchinds . append ( mind [ 0 ] )
        else :
            ext_matchinds . append ( mind )
return kdt_matchinds , ext_matchinds


signature:
def total_proper_motion ( pmra , pmdecl , decl ) :
    '''This calculates the total proper motion of an object .
 Parameters
 pmra : float or array - like
 The proper motion ( s ) in right ascension , measured in mas / yr .
 pmdecl : float or array - like
 The proper motion ( s ) in declination , measured in mas / yr .
 decl : float or array - like
 The declination of the object ( s ) in decimal degrees .
 Returns
 float or array - like
 The total proper motion ( s ) of the object ( s ) in mas / yr .'''

implementation:
pm = np . sqrt ( pmdecl * pmdecl + pmra * pmra * np . cos ( np . radians ( decl ) ) * np . cos ( np . radians ( decl ) ) )
return pm


signature:
def equatorial_to_galactic ( ra , decl , equinox = 'J2000' ) :
    '''This converts from equatorial coords to galactic coords .
 Parameters
 ra : float or array - like
 Right ascension values ( s ) in decimal degrees .
 decl : float or array - like
 Declination value ( s ) in decimal degrees .
 equinox : str
 The equinox that the coordinates are measured at . This must be
 recognizable by Astropy ' s ` SkyCoord ` class .
 Returns
 tuple of ( float , float ) or tuple of ( np . array , np . array )
 The galactic coordinates ( l , b ) for each element of the input
 ( ` ra ` , ` decl ` ) .'''

implementation:
# convert the ra / decl to gl , gb
radecl = SkyCoord ( ra = ra * u . degree , dec = decl * u . degree , equinox = equinox )
gl = radecl . galactic . l . degree
gb = radecl . galactic . b . degree
return gl , gb


signature:
def galactic_to_equatorial ( gl , gb ) :
    '''This converts from galactic coords to equatorial coordinates .
 Parameters
 gl : float or array - like
 Galactic longitude values ( s ) in decimal degrees .
 gb : float or array - like
 Galactic latitude value ( s ) in decimal degrees .
 Returns
 tuple of ( float , float ) or tuple of ( np . array , np . array )
 The equatorial coordinates ( RA , DEC ) for each element of the input
 ( ` gl ` , ` gb ` ) in decimal degrees . These are reported in the ICRS frame .'''

implementation:
gal = SkyCoord ( gl * u . degree , gl * u . degree , frame = 'galactic' )
transformed = gal . transform_to ( 'icrs' )
return transformed . ra . degree , transformed . dec . degree


signature:
def xieta_from_radecl ( inra , indecl , incenterra , incenterdecl , deg = True ) :
    '''This returns the image - plane projected xi - eta coords for inra , indecl .
 Parameters
 inra , indecl : array - like
 The equatorial coordinates to get the xi , eta coordinates for in decimal
 degrees or radians .
 incenterra , incenterdecl : float
 The center coordinate values to use to calculate the plane - projected
 coordinates around .
 deg : bool
 If this is True , the input angles are assumed to be in degrees and the
 output is in degrees as well .
 Returns
 tuple of np . arrays
 This is the ( ` xi ` , ` eta ` ) coordinate pairs corresponding to the
 image - plane projected coordinates for each pair of input equatorial
 coordinates in ( ` inra ` , ` indecl ` ) .'''

implementation:
if deg :
    ra = np . radians ( inra )
    decl = np . radians ( indecl )
    centerra = np . radians ( incenterra )
    centerdecl = np . radians ( incenterdecl )
else :
    ra = inra
    decl = indecl
    centerra = incenterra
    centerdecl = incenterdecl
cdecc = np . cos ( centerdecl )
sdecc = np . sin ( centerdecl )
crac = np . cos ( centerra )
srac = np . sin ( centerra )
uu = np . cos ( decl ) * np . cos ( ra )
vv = np . cos ( decl ) * np . sin ( ra )
ww = np . sin ( decl )
uun = uu * cdecc * crac + vv * cdecc * srac + ww * sdecc
vvn = - uu * srac + vv * crac
wwn = - uu * sdecc * crac - vv * sdecc * srac + ww * cdecc
denom = vvn * vvn + wwn * wwn
aunn = np . zeros_like ( uun )
aunn [ uun >= 1.0 ] = 0.0
aunn [ uun < 1.0 ] = np . arccos ( uun )
xi , eta = np . zeros_like ( aunn ) , np . zeros_like ( aunn )
xi [ ( aunn <= 0.0 ) | ( denom <= 0.0 ) ] = 0.0
eta [ ( aunn <= 0.0 ) | ( denom <= 0.0 ) ] = 0.0
sdenom = np . sqrt ( denom )
xi [ ( aunn > 0.0 ) | ( denom > 0.0 ) ] = aunn * vvn / sdenom
eta [ ( aunn > 0.0 ) | ( denom > 0.0 ) ] = aunn * wwn / sdenom
if deg :
    return np . degrees ( xi ) , np . degrees ( eta )
else :
    return xi , eta


signature:
def generate_transit_lightcurve ( times , mags = None , errs = None , paramdists = { 'transitperiod' : sps . uniform ( loc = 0.1 , scale = 49.9 ) , 'transitdepth' : sps . uniform ( loc = 1.0e-4 , scale = 2.0e-2 ) , 'transitduration' : sps . uniform ( loc = 0.01 , scale = 0.29 ) } , magsarefluxes = False , ) :
    '''This generates fake planet transit light curves .
 Parameters
 times : np . array
 This is an array of time values that will be used as the time base .
 mags , errs : np . array
 These arrays will have the model added to them . If either is
 None , ` np . full _ like ( times , 0.0 ) ` will used as a substitute and the model
 light curve will be centered around 0.0.
 paramdists : dict
 This is a dict containing parameter distributions to use for the
 model params , containing the following keys : :
 { ' transitperiod ' , ' transitdepth ' , ' transitduration ' }
 The values of these keys should all be ' frozen ' scipy . stats distribution
 objects , e . g . :
 https : / / docs . scipy . org / doc / scipy / reference / stats . html # continuous - distributions
 The variability epoch will be automatically chosen from a uniform
 distribution between ` times . min ( ) ` and ` times . max ( ) ` .
 The ingress duration will be automatically chosen from a uniform
 distribution ranging from 0.05 to 0.5 of the transitduration .
 The transitdepth will be flipped automatically as appropriate if
 ` magsarefluxes = True ` .
 magsarefluxes : bool
 If the generated time series is meant to be a flux time - series , set this
 to True to get the correct sign of variability amplitude .
 Returns
 dict
 A dict of the form below is returned : :
 { ' vartype ' : ' planet ' ,
 ' params ' : { ' transitperiod ' : generated value of period ,
 ' transitepoch ' : generated value of epoch ,
 ' transitdepth ' : generated value of transit depth ,
 ' transitduration ' : generated value of transit duration ,
 ' ingressduration ' : generated value of transit ingress
 duration } ,
 ' times ' : the model times ,
 ' mags ' : the model mags ,
 ' errs ' : the model errs ,
 ' varperiod ' : the generated period of variability = = ' transitperiod '
 ' varamplitude ' : the generated amplitude of
 variability = = ' transitdepth ' }'''

implementation:
if mags is None :
    mags = np . full_like ( times , 0.0 )
if errs is None :
    errs = np . full_like ( times , 0.0 )
# choose the epoch
epoch = npr . random ( ) * ( times . max ( ) - times . min ( ) ) + times . min ( )
# choose the period , depth , duration
period = paramdists [ 'transitperiod' ] . rvs ( size = 1 )
depth = paramdists [ 'transitdepth' ] . rvs ( size = 1 )
duration = paramdists [ 'transitduration' ] . rvs ( size = 1 )
# figure out the ingress duration
ingduration = npr . random ( ) * ( 0.5 * duration - 0.05 * duration ) + 0.05 * duration
# fix the transit depth if it needs to be flipped
if magsarefluxes and depth < 0.0 :
    depth = - depth
elif not magsarefluxes and depth > 0.0 :
    depth = - depth
# generate the model
modelmags , phase , ptimes , pmags , perrs = ( transits . trapezoid_transit_func ( [ period , epoch , depth , duration , ingduration ] , times , mags , errs ) )
# resort in original time order
timeind = np . argsort ( ptimes )
mtimes = ptimes [ timeind ]
mmags = modelmags [ timeind ]
merrs = perrs [ timeind ]
# return a dict with everything
modeldict = { 'vartype' : 'planet' , 'params' : { x : np . asscalar ( y ) for x , y in zip ( [ 'transitperiod' , 'transitepoch' , 'transitdepth' , 'transitduration' , 'ingressduration' ] , [ period , epoch , depth , duration , ingduration ] ) } , 'times' : mtimes , 'mags' : mmags , 'errs' : merrs , # these are standard keys that help with later characterization of
# variability as a function period , variability amplitude , object mag ,
# ndet , etc .
'varperiod' : period , 'varamplitude' : depth }
return modeldict


signature:
def generate_eb_lightcurve ( times , mags = None , errs = None , paramdists = { 'period' : sps . uniform ( loc = 0.2 , scale = 99.8 ) , 'pdepth' : sps . uniform ( loc = 1.0e-4 , scale = 0.7 ) , 'pduration' : sps . uniform ( loc = 0.01 , scale = 0.44 ) , 'depthratio' : sps . uniform ( loc = 0.01 , scale = 0.99 ) , 'secphase' : sps . norm ( loc = 0.5 , scale = 0.1 ) } , magsarefluxes = False , ) :
    '''This generates fake EB light curves .
 Parameters
 times : np . array
 This is an array of time values that will be used as the time base .
 mags , errs : np . array
 These arrays will have the model added to them . If either is
 None , ` np . full _ like ( times , 0.0 ) ` will used as a substitute and the model
 light curve will be centered around 0.0.
 paramdists : dict
 This is a dict containing parameter distributions to use for the
 model params , containing the following keys : :
 { ' period ' , ' pdepth ' , ' pduration ' , ' depthratio ' , ' secphase ' }
 The values of these keys should all be ' frozen ' scipy . stats distribution
 objects , e . g . :
 https : / / docs . scipy . org / doc / scipy / reference / stats . html # continuous - distributions
 The variability epoch will be automatically chosen from a uniform
 distribution between ` times . min ( ) ` and ` times . max ( ) ` .
 The ` pdepth ` will be flipped automatically as appropriate if
 ` magsarefluxes = True ` .
 magsarefluxes : bool
 If the generated time series is meant to be a flux time - series , set this
 to True to get the correct sign of variability amplitude .
 Returns
 dict
 A dict of the form below is returned : :
 { ' vartype ' : ' EB ' ,
 ' params ' : { ' period ' : generated value of period ,
 ' epoch ' : generated value of epoch ,
 ' pdepth ' : generated value of priary eclipse depth ,
 ' pduration ' : generated value of prim eclipse duration ,
 ' depthratio ' : generated value of prim / sec eclipse
 depth ratio } ,
 ' times ' : the model times ,
 ' mags ' : the model mags ,
 ' errs ' : the model errs ,
 ' varperiod ' : the generated period of variability = = ' period '
 ' varamplitude ' : the generated amplitude of
 variability = = ' pdepth ' }'''

implementation:
if mags is None :
    mags = np . full_like ( times , 0.0 )
if errs is None :
    errs = np . full_like ( times , 0.0 )
# choose the epoch
epoch = npr . random ( ) * ( times . max ( ) - times . min ( ) ) + times . min ( )
# choose the period , pdepth , duration , depthratio
period = paramdists [ 'period' ] . rvs ( size = 1 )
pdepth = paramdists [ 'pdepth' ] . rvs ( size = 1 )
pduration = paramdists [ 'pduration' ] . rvs ( size = 1 )
depthratio = paramdists [ 'depthratio' ] . rvs ( size = 1 )
secphase = paramdists [ 'secphase' ] . rvs ( size = 1 )
# fix the transit depth if it needs to be flipped
if magsarefluxes and pdepth < 0.0 :
    pdepth = - pdepth
elif not magsarefluxes and pdepth > 0.0 :
    pdepth = - pdepth
# generate the model
modelmags , phase , ptimes , pmags , perrs = ( eclipses . invgauss_eclipses_func ( [ period , epoch , pdepth , pduration , depthratio , secphase ] , times , mags , errs ) )
# resort in original time order
timeind = np . argsort ( ptimes )
mtimes = ptimes [ timeind ]
mmags = modelmags [ timeind ]
merrs = perrs [ timeind ]
# return a dict with everything
modeldict = { 'vartype' : 'EB' , 'params' : { x : np . asscalar ( y ) for x , y in zip ( [ 'period' , 'epoch' , 'pdepth' , 'pduration' , 'depthratio' ] , [ period , epoch , pdepth , pduration , depthratio ] ) } , 'times' : mtimes , 'mags' : mmags , 'errs' : merrs , 'varperiod' : period , 'varamplitude' : pdepth , }
return modeldict


signature:
def generate_flare_lightcurve ( times , mags = None , errs = None , paramdists = { # flare peak amplitude from 0.01 mag to 1.0 mag above median . this
# is tuned for redder bands , flares are much stronger in bluer
# bands , so tune appropriately for your situation .
'amplitude' : sps . uniform ( loc = 0.01 , scale = 0.99 ) , # up to 5 flares per LC and at least 1
'nflares' : [ 1 , 5 ] , # 10 minutes to 1 hour for rise stdev
'risestdev' : sps . uniform ( loc = 0.007 , scale = 0.04 ) , # 1 hour to 4 hours for decay time constant
'decayconst' : sps . uniform ( loc = 0.04 , scale = 0.163 ) } , magsarefluxes = False , ) :
    '''This generates fake flare light curves .
 Parameters
 times : np . array
 This is an array of time values that will be used as the time base .
 mags , errs : np . array
 These arrays will have the model added to them . If either is
 None , ` np . full _ like ( times , 0.0 ) ` will used as a substitute and the model
 light curve will be centered around 0.0.
 paramdists : dict
 This is a dict containing parameter distributions to use for the
 model params , containing the following keys : :
 { ' amplitude ' , ' nflares ' , ' risestdev ' , ' decayconst ' }
 The values of these keys should all be ' frozen ' scipy . stats distribution
 objects , e . g . :
 https : / / docs . scipy . org / doc / scipy / reference / stats . html # continuous - distributions
 The ` flare _ peak _ time ` for each flare will be generated automatically
 between ` times . min ( ) ` and ` times . max ( ) ` using a uniform distribution .
 The ` amplitude ` will be flipped automatically as appropriate if
 ` magsarefluxes = True ` .
 magsarefluxes : bool
 If the generated time series is meant to be a flux time - series , set this
 to True to get the correct sign of variability amplitude .
 Returns
 dict
 A dict of the form below is returned : :
 { ' vartype ' : ' flare ' ,
 ' params ' : { ' amplitude ' : generated value of flare amplitudes ,
 ' nflares ' : generated value of number of flares ,
 ' risestdev ' : generated value of stdev of rise time ,
 ' decayconst ' : generated value of decay constant ,
 ' peaktime ' : generated value of flare peak time } ,
 ' times ' : the model times ,
 ' mags ' : the model mags ,
 ' errs ' : the model errs ,
 ' varamplitude ' : the generated amplitude of
 variability = = ' amplitude ' }'''

implementation:
if mags is None :
    mags = np . full_like ( times , 0.0 )
if errs is None :
    errs = np . full_like ( times , 0.0 )
nflares = npr . randint ( paramdists [ 'nflares' ] [ 0 ] , high = paramdists [ 'nflares' ] [ 1 ] )
# generate random flare peak times based on the number of flares
flarepeaktimes = ( npr . random ( size = nflares ) * ( times . max ( ) - times . min ( ) ) + times . min ( ) )
# now add the flares to the time - series
params = { 'nflares' : nflares }
for flareind , peaktime in zip ( range ( nflares ) , flarepeaktimes ) : # choose the amplitude , rise stdev and decay time constant
    amp = paramdists [ 'amplitude' ] . rvs ( size = 1 )
    risestdev = paramdists [ 'risestdev' ] . rvs ( size = 1 )
    decayconst = paramdists [ 'decayconst' ] . rvs ( size = 1 )
    # fix the transit depth if it needs to be flipped
    if magsarefluxes and amp < 0.0 :
        amp = - amp
    elif not magsarefluxes and amp > 0.0 :
        amp = - amp
    # add this flare to the light curve
    modelmags , ptimes , pmags , perrs = ( flares . flare_model ( [ amp , peaktime , risestdev , decayconst ] , times , mags , errs ) )
    # update the mags
    mags = modelmags
    # add the flare params to the modeldict
    params [ flareind ] = { 'peaktime' : peaktime , 'amplitude' : amp , 'risestdev' : risestdev , 'decayconst' : decayconst }
# done with all flares
# return a dict with everything
modeldict = { 'vartype' : 'flare' , 'params' : params , 'times' : times , 'mags' : mags , 'errs' : errs , 'varperiod' : None , # FIXME : this is complicated because we can have multiple flares
# figure out a good way to handle this upstream
'varamplitude' : [ params [ x ] [ 'amplitude' ] for x in range ( params [ 'nflares' ] ) ] , }
return modeldict


signature:
def generate_sinusoidal_lightcurve ( times , mags = None , errs = None , paramdists = { 'period' : sps . uniform ( loc = 0.04 , scale = 500.0 ) , 'fourierorder' : [ 2 , 10 ] , 'amplitude' : sps . uniform ( loc = 0.1 , scale = 0.9 ) , 'phioffset' : 0.0 , } , magsarefluxes = False ) :
    '''This generates fake sinusoidal light curves .
 This can be used for a variety of sinusoidal variables , e . g . RRab , RRc ,
 Cepheids , Miras , etc . The functions that generate these model LCs below
 implement the following table : :
 # # FOURIER PARAMS FOR SINUSOIDAL VARIABLES
 # type fourier period [ days ]
 # order dist limits dist
 # RRab 8 to 10 uniform 0.45 - - 0.80 uniform
 # RRc 3 to 6 uniform 0.10 - - 0.40 uniform
 # HADS 7 to 9 uniform 0.04 - - 0.10 uniform
 # rotator 2 to 5 uniform 0.80 - - 120.0 uniform
 # LPV 2 to 5 uniform 250 - - 500.0 uniform
 FIXME : for better model LCs , figure out how scipy . signal . butter works and
 low - pass filter using scipy . signal . filtfilt .
 Parameters
 times : np . array
 This is an array of time values that will be used as the time base .
 mags , errs : np . array
 These arrays will have the model added to them . If either is
 None , ` np . full _ like ( times , 0.0 ) ` will used as a substitute and the model
 light curve will be centered around 0.0.
 paramdists : dict
 This is a dict containing parameter distributions to use for the
 model params , containing the following keys : :
 { ' period ' , ' fourierorder ' , ' amplitude ' , ' phioffset ' }
 The values of these keys should all be ' frozen ' scipy . stats distribution
 objects , e . g . :
 https : / / docs . scipy . org / doc / scipy / reference / stats . html # continuous - distributions
 The variability epoch will be automatically chosen from a uniform
 distribution between ` times . min ( ) ` and ` times . max ( ) ` .
 The ` amplitude ` will be flipped automatically as appropriate if
 ` magsarefluxes = True ` .
 magsarefluxes : bool
 If the generated time series is meant to be a flux time - series , set this
 to True to get the correct sign of variability amplitude .
 Returns
 dict
 A dict of the form below is returned : :
 { ' vartype ' : ' sinusoidal ' ,
 ' params ' : { ' period ' : generated value of period ,
 ' epoch ' : generated value of epoch ,
 ' amplitude ' : generated value of amplitude ,
 ' fourierorder ' : generated value of fourier order ,
 ' fourieramps ' : generated values of fourier amplitudes ,
 ' fourierphases ' : generated values of fourier phases } ,
 ' times ' : the model times ,
 ' mags ' : the model mags ,
 ' errs ' : the model errs ,
 ' varperiod ' : the generated period of variability = = ' period '
 ' varamplitude ' : the generated amplitude of
 variability = = ' amplitude ' }'''

implementation:
if mags is None :
    mags = np . full_like ( times , 0.0 )
if errs is None :
    errs = np . full_like ( times , 0.0 )
# choose the epoch
epoch = npr . random ( ) * ( times . max ( ) - times . min ( ) ) + times . min ( )
# choose the period , fourierorder , and amplitude
period = paramdists [ 'period' ] . rvs ( size = 1 )
fourierorder = npr . randint ( paramdists [ 'fourierorder' ] [ 0 ] , high = paramdists [ 'fourierorder' ] [ 1 ] )
amplitude = paramdists [ 'amplitude' ] . rvs ( size = 1 )
# fix the amplitude if it needs to be flipped
if magsarefluxes and amplitude < 0.0 :
    amplitude = - amplitude
elif not magsarefluxes and amplitude > 0.0 :
    amplitude = - amplitude
# generate the amplitudes and phases of the Fourier components
ampcomps = [ abs ( amplitude / 2.0 ) / float ( x ) for x in range ( 1 , fourierorder + 1 ) ]
phacomps = [ paramdists [ 'phioffset' ] * float ( x ) for x in range ( 1 , fourierorder + 1 ) ]
# now that we have our amp and pha components , generate the light curve
modelmags , phase , ptimes , pmags , perrs = sinusoidal . sine_series_sum ( [ period , epoch , ampcomps , phacomps ] , times , mags , errs )
# resort in original time order
timeind = np . argsort ( ptimes )
mtimes = ptimes [ timeind ]
mmags = modelmags [ timeind ]
merrs = perrs [ timeind ]
mphase = phase [ timeind ]
# return a dict with everything
modeldict = { 'vartype' : 'sinusoidal' , 'params' : { x : y for x , y in zip ( [ 'period' , 'epoch' , 'amplitude' , 'fourierorder' , 'fourieramps' , 'fourierphases' ] , [ period , epoch , amplitude , fourierorder , ampcomps , phacomps ] ) } , 'times' : mtimes , 'mags' : mmags , 'errs' : merrs , 'phase' : mphase , # these are standard keys that help with later characterization of
# variability as a function period , variability amplitude , object mag ,
# ndet , etc .
'varperiod' : period , 'varamplitude' : amplitude }
return modeldict


signature:
def generate_rrab_lightcurve ( times , mags = None , errs = None , paramdists = { 'period' : sps . uniform ( loc = 0.45 , scale = 0.35 ) , 'fourierorder' : [ 8 , 11 ] , 'amplitude' : sps . uniform ( loc = 0.4 , scale = 0.5 ) , 'phioffset' : np . pi , } , magsarefluxes = False ) :
    '''This generates fake RRab light curves .
 Parameters
 times : np . array
 This is an array of time values that will be used as the time base .
 mags , errs : np . array
 These arrays will have the model added to them . If either is
 None , ` np . full _ like ( times , 0.0 ) ` will used as a substitute and the model
 light curve will be centered around 0.0.
 paramdists : dict
 This is a dict containing parameter distributions to use for the
 model params , containing the following keys : :
 { ' period ' , ' fourierorder ' , ' amplitude ' }
 The values of these keys should all be ' frozen ' scipy . stats distribution
 objects , e . g . :
 https : / / docs . scipy . org / doc / scipy / reference / stats . html # continuous - distributions
 The variability epoch will be automatically chosen from a uniform
 distribution between ` times . min ( ) ` and ` times . max ( ) ` .
 The ` amplitude ` will be flipped automatically as appropriate if
 ` magsarefluxes = True ` .
 magsarefluxes : bool
 If the generated time series is meant to be a flux time - series , set this
 to True to get the correct sign of variability amplitude .
 Returns
 dict
 A dict of the form below is returned : :
 { ' vartype ' : ' RRab ' ,
 ' params ' : { ' period ' : generated value of period ,
 ' epoch ' : generated value of epoch ,
 ' amplitude ' : generated value of amplitude ,
 ' fourierorder ' : generated value of fourier order ,
 ' fourieramps ' : generated values of fourier amplitudes ,
 ' fourierphases ' : generated values of fourier phases } ,
 ' times ' : the model times ,
 ' mags ' : the model mags ,
 ' errs ' : the model errs ,
 ' varperiod ' : the generated period of variability = = ' period '
 ' varamplitude ' : the generated amplitude of
 variability = = ' amplitude ' }'''

implementation:
modeldict = generate_sinusoidal_lightcurve ( times , mags = mags , errs = errs , paramdists = paramdists , magsarefluxes = magsarefluxes )
modeldict [ 'vartype' ] = 'RRab'
return modeldict


signature:
def make_fakelc ( lcfile , outdir , magrms = None , randomizemags = True , randomizecoords = False , lcformat = 'hat-sql' , lcformatdir = None , timecols = None , magcols = None , errcols = None ) :
    '''This preprocesses an input real LC and sets it up to be a fake LC .
 Parameters
 lcfile : str
 This is an input light curve file that will be used to copy over the
 time - base . This will be used to generate the time - base for fake light
 curves to provide a realistic simulation of the observing window
 function .
 outdir : str
 The output directory where the the fake light curve will be written .
 magrms : dict
 This is a dict containing the SDSS r mag - RMS ( SDSS rmag - MAD preferably )
 relation based on all light curves that the input lcfile is from . This
 will be used to generate the median mag and noise corresponding to the
 magnitude chosen for this fake LC .
 randomizemags : bool
 If this is True , then a random mag between the first and last magbin in
 magrms will be chosen as the median mag for this light curve . This
 choice will be weighted by the mag bin probability obtained from the
 magrms kwarg . Otherwise , the median mag will be taken from the input
 lcfile ' s lcdict [ ' objectinfo ' ] [ ' sdssr ' ] key or a transformed SDSS r mag
 generated from the input lcfile ' s lcdict [ ' objectinfo ' ] [ ' jmag ' ] ,
 [ ' hmag ' ] , and [ ' kmag ' ] keys . The magrms relation for each magcol will be
 used to generate Gaussian noise at the correct level for the magbin this
 light curve ' s median mag falls into .
 randomizecoords : bool
 If this is True , will randomize the RA , DEC of the output fake object
 and not copy over the RA / DEC from the real input object .
 lcformat : str
 This is the ` formatkey ` associated with your input real light curve
 format , which you previously passed in to the ` lcproc . register _ lcformat `
 function . This will be used to look up how to find and read the light
 curve specified in ` lcfile ` .
 lcformatdir : str or None
 If this is provided , gives the path to a directory when you ' ve stored
 your lcformat description JSONs , other than the usual directories lcproc
 knows to search for them in . Use this along with ` lcformat ` to specify
 an LC format JSON file that ' s not currently registered with lcproc .
 timecols : list of str or None
 The timecol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curve .
 magcols : list of str or None
 The magcol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curve .
 errcols : list of str or None
 The errcol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curve .
 Returns
 tuple
 A tuple of the following form is returned : :
 ( fakelc _ fpath ,
 fakelc _ lcdict [ ' columns ' ] ,
 fakelc _ lcdict [ ' objectinfo ' ] ,
 fakelc _ lcdict [ ' moments ' ] )'''

implementation:
try :
    formatinfo = get_lcformat ( lcformat , use_lcformat_dir = lcformatdir )
    if formatinfo :
        ( fileglob , readerfunc , dtimecols , dmagcols , derrcols , magsarefluxes , normfunc ) = formatinfo
    else :
        LOGERROR ( "can't figure out the light curve format" )
        return None
except Exception as e :
    LOGEXCEPTION ( "can't figure out the light curve format" )
    return None
# override the default timecols , magcols , and errcols
# using the ones provided to the function
if timecols is None :
    timecols = dtimecols
if magcols is None :
    magcols = dmagcols
if errcols is None :
    errcols = derrcols
# read in the light curve
lcdict = readerfunc ( lcfile )
if isinstance ( lcdict , tuple ) and isinstance ( lcdict [ 0 ] , dict ) :
    lcdict = lcdict [ 0 ]
# set up the fakelcdict with a randomly assigned objectid
fakeobjectid = sha512 ( npr . bytes ( 12 ) ) . hexdigest ( ) [ - 8 : ]
fakelcdict = { 'objectid' : fakeobjectid , 'objectinfo' : { 'objectid' : fakeobjectid } , 'columns' : [ ] , 'moments' : { } , 'origformat' : lcformat , }
# now , get the actual mag of this object and other info and use that to
# populate the corresponding entries of the fakelcdict objectinfo
if ( 'objectinfo' in lcdict and isinstance ( lcdict [ 'objectinfo' ] , dict ) ) :
    objectinfo = lcdict [ 'objectinfo' ]
    # get the RA
    if ( not randomizecoords and 'ra' in objectinfo and objectinfo [ 'ra' ] is not None and np . isfinite ( objectinfo [ 'ra' ] ) ) :
        fakelcdict [ 'objectinfo' ] [ 'ra' ] = objectinfo [ 'ra' ]
    else : # if there ' s no RA available , we ' ll assign a random one between 0
    # and 360.0
        LOGWARNING ( '%s: assigning a random right ascension' % lcfile )
        fakelcdict [ 'objectinfo' ] [ 'ra' ] = npr . random ( ) * 360.0
    # get the DEC
    if ( not randomizecoords and 'decl' in objectinfo and objectinfo [ 'decl' ] is not None and np . isfinite ( objectinfo [ 'decl' ] ) ) :
        fakelcdict [ 'objectinfo' ] [ 'decl' ] = objectinfo [ 'decl' ]
    else : # if there ' s no DECL available , we ' ll assign a random one between
    # -90.0 and + 90.0
        LOGWARNING ( ' %s: assigning a random declination' % lcfile )
        fakelcdict [ 'objectinfo' ] [ 'decl' ] = npr . random ( ) * 180.0 - 90.0
    # get the SDSS r mag for this object
    # this will be used for getting the eventual mag - RMS relation later
    if ( ( not randomizemags ) and 'sdssr' in objectinfo and objectinfo [ 'sdssr' ] is not None and np . isfinite ( objectinfo [ 'sdssr' ] ) ) :
        fakelcdict [ 'objectinfo' ] [ 'sdssr' ] = objectinfo [ 'sdssr' ]
    # if the SDSS r is unavailable , but we have J , H , K : use those to get
    # the SDSS r by using transformations
    elif ( ( not randomizemags ) and ( 'jmag' in objectinfo and objectinfo [ 'jmag' ] is not None and np . isfinite ( objectinfo [ 'jmag' ] ) ) and ( 'hmag' in objectinfo and objectinfo [ 'hmag' ] is not None and np . isfinite ( objectinfo [ 'hmag' ] ) ) and ( 'kmag' in objectinfo and objectinfo [ 'kmag' ] is not None and np . isfinite ( objectinfo [ 'kmag' ] ) ) ) :
        LOGWARNING ( 'used JHK mags to generate an SDSS r mag for %s' % lcfile )
        fakelcdict [ 'objectinfo' ] [ 'sdssr' ] = jhk_to_sdssr ( objectinfo [ 'jmag' ] , objectinfo [ 'hmag' ] , objectinfo [ 'kmag' ] )
    # if there are no mags available or we ' re specically told to randomize
    # them , generate a random mag between 8 and 16.0
    elif randomizemags and magrms :
        LOGWARNING ( ' %s: assigning a random mag weighted by mag ' 'bin probabilities' % lcfile )
        magbins = magrms [ magcols [ 0 ] ] [ 'binned_sdssr_median' ]
        binprobs = magrms [ magcols [ 0 ] ] [ 'magbin_probabilities' ]
        # this is the center of the magbin chosen
        magbincenter = npr . choice ( magbins , size = 1 , p = binprobs )
        # in this magbin , choose between center and - + 0.25 mag
        chosenmag = ( npr . random ( ) * ( ( magbincenter + 0.25 ) - ( magbincenter - 0.25 ) ) + ( magbincenter - 0.25 ) )
        fakelcdict [ 'objectinfo' ] [ 'sdssr' ] = np . asscalar ( chosenmag )
    # if there are no mags available at all , generate a random mag
    # between 8 and 16.0
    else :
        LOGWARNING ( ' %s: assigning a random mag from ' 'uniform distribution between 8.0 and 16.0' % lcfile )
        fakelcdict [ 'objectinfo' ] [ 'sdssr' ] = npr . random ( ) * 8.0 + 8.0
# if there ' s no info available , generate fake info
else :
    LOGWARNING ( 'no object information found in %s, ' 'generating random ra, decl, sdssr' % lcfile )
    fakelcdict [ 'objectinfo' ] [ 'ra' ] = npr . random ( ) * 360.0
    fakelcdict [ 'objectinfo' ] [ 'decl' ] = npr . random ( ) * 180.0 - 90.0
    fakelcdict [ 'objectinfo' ] [ 'sdssr' ] = npr . random ( ) * 8.0 + 8.0
# NOW FILL IN THE TIMES , MAGS , ERRS
# get the time columns
for tcind , tcol in enumerate ( timecols ) :
    if '.' in tcol :
        tcolget = tcol . split ( '.' )
    else :
        tcolget = [ tcol ]
    if tcol not in fakelcdict :
        fakelcdict [ tcol ] = _dict_get ( lcdict , tcolget )
        fakelcdict [ 'columns' ] . append ( tcol )
        # update the ndet with the first time column ' s size . it ' s possible
        # that different time columns have different lengths , but that would
        # be weird and we won ' t deal with it for now
        if tcind == 0 :
            fakelcdict [ 'objectinfo' ] [ 'ndet' ] = fakelcdict [ tcol ] . size
# get the mag columns
for mcol in magcols :
    if '.' in mcol :
        mcolget = mcol . split ( '.' )
    else :
        mcolget = [ mcol ]
    # put the mcol in only once
    if mcol not in fakelcdict :
        measuredmags = _dict_get ( lcdict , mcolget )
        measuredmags = measuredmags [ np . isfinite ( measuredmags ) ]
        # if we ' re randomizing , get the mags from the interpolated mag - RMS
        # relation
        if ( randomizemags and magrms and mcol in magrms and 'interpolated_magmad' in magrms [ mcol ] and magrms [ mcol ] [ 'interpolated_magmad' ] is not None ) :
            interpfunc = magrms [ mcol ] [ 'interpolated_magmad' ]
            lcmad = interpfunc ( fakelcdict [ 'objectinfo' ] [ 'sdssr' ] )
            fakelcdict [ 'moments' ] [ mcol ] = { 'median' : fakelcdict [ 'objectinfo' ] [ 'sdssr' ] , 'mad' : lcmad }
        # if we ' re not randomizing , get the median and MAD from the light
        # curve itself
        else : # we require at least 10 finite measurements
            if measuredmags . size > 9 :
                measuredmedian = np . median ( measuredmags )
                measuredmad = np . median ( np . abs ( measuredmags - measuredmedian ) )
                fakelcdict [ 'moments' ] [ mcol ] = { 'median' : measuredmedian , 'mad' : measuredmad }
            # if there aren ' t enough measurements in this LC , try to get the
            # median and RMS from the interpolated mag - RMS relation first
            else :
                if ( magrms and mcol in magrms and 'interpolated_magmad' in magrms [ mcol ] and magrms [ mcol ] [ 'interpolated_magmad' ] is not None ) :
                    LOGWARNING ( 'input LC %s does not have enough ' 'finite measurements, ' 'generating mag moments from ' 'fakelc sdssr and the mag-RMS relation' % lcfile )
                    interpfunc = magrms [ mcol ] [ 'interpolated_magmad' ]
                    lcmad = interpfunc ( fakelcdict [ 'objectinfo' ] [ 'sdssr' ] )
                    fakelcdict [ 'moments' ] [ mcol ] = { 'median' : fakelcdict [ 'objectinfo' ] [ 'sdssr' ] , 'mad' : lcmad }
                # if we don ' t have the mag - RMS relation either , then we
                # can ' t do anything for this light curve , generate a random
                # MAD between 5e - 4 and 0.1
                else :
                    LOGWARNING ( 'input LC %s does not have enough ' 'finite measurements and ' 'no mag-RMS relation provided ' 'assigning a random MAD between 5.0e-4 and 0.1' % lcfile )
                    fakelcdict [ 'moments' ] [ mcol ] = { 'median' : fakelcdict [ 'objectinfo' ] [ 'sdssr' ] , 'mad' : npr . random ( ) * ( 0.1 - 5.0e-4 ) + 5.0e-4 }
        # the magnitude column is set to all zeros initially . this will be
        # filled in by the add _ fakelc _ variability function below
        fakelcdict [ mcol ] = np . full_like ( _dict_get ( lcdict , mcolget ) , 0.0 )
        fakelcdict [ 'columns' ] . append ( mcol )
# get the err columns
for mcol , ecol in zip ( magcols , errcols ) :
    if '.' in ecol :
        ecolget = ecol . split ( '.' )
    else :
        ecolget = [ ecol ]
    if ecol not in fakelcdict :
        measurederrs = _dict_get ( lcdict , ecolget )
        measurederrs = measurederrs [ np . isfinite ( measurederrs ) ]
        # if we ' re randomizing , get the errs from the interpolated mag - RMS
        # relation
        if ( randomizemags and magrms and mcol in magrms and 'interpolated_magmad' in magrms [ mcol ] and magrms [ mcol ] [ 'interpolated_magmad' ] is not None ) :
            interpfunc = magrms [ mcol ] [ 'interpolated_magmad' ]
            lcmad = interpfunc ( fakelcdict [ 'objectinfo' ] [ 'sdssr' ] )
            # the median of the errs = lcmad
            # the mad of the errs is 0.1 x lcmad
            fakelcdict [ 'moments' ] [ ecol ] = { 'median' : lcmad , 'mad' : 0.1 * lcmad }
        else : # we require at least 10 finite measurements
        # we ' ll calculate the median and MAD of the errs to use later on
            if measurederrs . size > 9 :
                measuredmedian = np . median ( measurederrs )
                measuredmad = np . median ( np . abs ( measurederrs - measuredmedian ) )
                fakelcdict [ 'moments' ] [ ecol ] = { 'median' : measuredmedian , 'mad' : measuredmad }
            else :
                if ( magrms and mcol in magrms and 'interpolated_magmad' in magrms [ mcol ] and magrms [ mcol ] [ 'interpolated_magmad' ] is not None ) :
                    LOGWARNING ( 'input LC %s does not have enough ' 'finite measurements, ' 'generating err moments from ' 'the mag-RMS relation' % lcfile )
                    interpfunc = magrms [ mcol ] [ 'interpolated_magmad' ]
                    lcmad = interpfunc ( fakelcdict [ 'objectinfo' ] [ 'sdssr' ] )
                    fakelcdict [ 'moments' ] [ ecol ] = { 'median' : lcmad , 'mad' : 0.1 * lcmad }
                # if we don ' t have the mag - RMS relation either , then we
                # can ' t do anything for this light curve , generate a random
                # MAD between 5e - 4 and 0.1
                else :
                    LOGWARNING ( 'input LC %s does not have ' 'enough finite measurements and ' 'no mag-RMS relation provided, ' 'generating errs randomly' % lcfile )
                    fakelcdict [ 'moments' ] [ ecol ] = { 'median' : npr . random ( ) * ( 0.01 - 5.0e-4 ) + 5.0e-4 , 'mad' : npr . random ( ) * ( 0.01 - 5.0e-4 ) + 5.0e-4 }
        # the errors column is set to all zeros initially . this will be
        # filled in by the add _ fakelc _ variability function below .
        fakelcdict [ ecol ] = np . full_like ( _dict_get ( lcdict , ecolget ) , 0.0 )
        fakelcdict [ 'columns' ] . append ( ecol )
# add the timecols , magcols , errcols to the lcdict
fakelcdict [ 'timecols' ] = timecols
fakelcdict [ 'magcols' ] = magcols
fakelcdict [ 'errcols' ] = errcols
# generate an output file name
fakelcfname = '%s-fakelc.pkl' % fakelcdict [ 'objectid' ]
fakelcfpath = os . path . abspath ( os . path . join ( outdir , fakelcfname ) )
# write this out to the output directory
with open ( fakelcfpath , 'wb' ) as outfd :
    pickle . dump ( fakelcdict , outfd , protocol = pickle . HIGHEST_PROTOCOL )
# return the fakelc path , its columns , info , and moments so we can put them
# into a collection DB later on
LOGINFO ( 'real LC %s -> fake LC %s OK' % ( lcfile , fakelcfpath ) )
return ( fakelcfpath , fakelcdict [ 'columns' ] , fakelcdict [ 'objectinfo' ] , fakelcdict [ 'moments' ] )


signature:
def collection_worker ( task ) :
    '''This wraps ` process _ fakelc ` for ` make _ fakelc _ collection ` below .
 Parameters
 task : tuple
 This is of the form : :
 task [ 0 ] = lcfile
 task [ 1 ] = outdir
 task [ 2 ] = magrms
 task [ 3 ] = dict with keys : { ' lcformat ' , ' timecols ' , ' magcols ' ,
 ' errcols ' , ' randomizeinfo ' }
 Returns
 tuple
 This returns a tuple of the form : :
 ( fakelc _ fpath ,
 fakelc _ lcdict [ ' columns ' ] ,
 fakelc _ lcdict [ ' objectinfo ' ] ,
 fakelc _ lcdict [ ' moments ' ] )'''

implementation:
lcfile , outdir , kwargs = task
try :
    fakelcresults = make_fakelc ( lcfile , outdir , ** kwargs )
    return fakelcresults
except Exception as e :
    LOGEXCEPTION ( 'could not process %s into a fakelc' % lcfile )
    return None


signature:
def make_fakelc_collection ( lclist , simbasedir , magrmsfrom , magrms_interpolate = 'quadratic' , magrms_fillvalue = 'extrapolate' , maxlcs = 25000 , maxvars = 2000 , randomizemags = True , randomizecoords = False , vartypes = ( 'EB' , 'RRab' , 'RRc' , 'cepheid' , 'rotator' , 'flare' , 'HADS' , 'planet' , 'LPV' ) , lcformat = 'hat-sql' , lcformatdir = None , timecols = None , magcols = None , errcols = None ) :
    '''This prepares light curves for the recovery sim .
 Collects light curves from ` lclist ` using a uniform sampling among
 them . Copies them to the ` simbasedir ` , zeroes out their mags and errs but
 keeps their time bases , also keeps their RMS and median mags for later
 use . Calculates the mag - rms relation for the entire collection and writes
 that to the ` simbasedir ` as well .
 The purpose of this function is to copy over the time base and mag - rms
 relation of an existing light curve collection to use it as the basis for a
 variability recovery simulation .
 This returns a pickle written to the ` simbasedir ` that contains all the
 information for the chosen ensemble of fake light curves and writes all
 generated light curves to the ` simbasedir / lightcurves ` directory . Run the
 ` add _ variability _ to _ fakelc _ collection ` function after this function to add
 variability of the specified type to these generated light curves .
 Parameters
 lclist : list of str
 This is a list of existing project light curves . This can be generated
 from : py : func : ` astrobase . lcproc . catalogs . make _ lclist ` or similar .
 simbasedir : str
 This is the directory to where the fake light curves and their
 information will be copied to .
 magrmsfrom : str or dict
 This is used to generate magnitudes and RMSes for the objects in the
 output collection of fake light curves . This arg is either a string
 pointing to an existing pickle file that must contain a dict or a dict
 variable that MUST have the following key - vals at a minimum : :
 { ' < magcol1 _ name > ' : {
 ' binned _ sdssr _ median ' : array of median mags for each magbin
 ' binned _ lcmad _ median ' : array of LC MAD values per magbin
 ' < magcol2 _ name > ' : {
 ' binned _ sdssr _ median ' : array of median mags for each magbin
 ' binned _ lcmad _ median ' : array of LC MAD values per magbin
 where ` magcol1 _ name ` , etc . are the same as the ` magcols ` listed in the
 magcols kwarg ( or the default magcols for the specified
 lcformat ) . Examples of the magrmsfrom dict ( or pickle ) required can be
 generated by the
 : py : func : ` astrobase . lcproc . varthreshold . variability _ threshold ` function .
 magrms _ interpolate , magrms _ fillvalue : str
 These are arguments that will be passed directly to the
 scipy . interpolate . interp1d function to generate interpolating functions
 for the mag - RMS relation . See :
 https : / / docs . scipy . org / doc / scipy / reference / generated / scipy . interpolate . interp1d . html
 for details .
 maxlcs : int
 This is the total number of light curves to choose from ` lclist ` and
 generate as fake LCs .
 maxvars : int
 This is the total number of fake light curves that will be marked as
 variable .
 vartypes : list of str
 This is a list of variable types to put into the collection . The
 vartypes for each fake variable star will be chosen uniformly from this
 list .
 lcformat : str
 This is the ` formatkey ` associated with your input real light curves '
 format , which you previously passed in to the ` lcproc . register _ lcformat `
 function . This will be used to look up how to find and read the light
 curves specified in ` lclist ` .
 lcformatdir : str or None
 If this is provided , gives the path to a directory when you ' ve stored
 your lcformat description JSONs , other than the usual directories lcproc
 knows to search for them in . Use this along with ` lcformat ` to specify
 an LC format JSON file that ' s not currently registered with lcproc .
 timecols : list of str or None
 The timecol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curves .
 magcols : list of str or None
 The magcol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curves .
 errcols : list of str or None
 The errcol keys to use from the input lcdict in generating the fake
 light curve . Fake LCs will be generated for each each
 timecol / magcol / errcol combination in the input light curves .
 Returns
 str
 Returns the string file name of a pickle containing all of the
 information for the fake LC collection that has been generated .'''

implementation:
try :
    formatinfo = get_lcformat ( lcformat , use_lcformat_dir = lcformatdir )
    if formatinfo :
        ( fileglob , readerfunc , dtimecols , dmagcols , derrcols , magsarefluxes , normfunc ) = formatinfo
    else :
        LOGERROR ( "can't figure out the light curve format" )
        return None
except Exception as e :
    LOGEXCEPTION ( "can't figure out the light curve format" )
    return None
# override the default timecols , magcols , and errcols
# using the ones provided to the function
if timecols is None :
    timecols = dtimecols
if magcols is None :
    magcols = dmagcols
if errcols is None :
    errcols = derrcols
if not isinstance ( lclist , np . ndarray ) :
    lclist = np . array ( lclist )
chosenlcs = npr . choice ( lclist , maxlcs , replace = False )
fakelcdir = os . path . join ( simbasedir , 'lightcurves' )
if not os . path . exists ( fakelcdir ) :
    os . makedirs ( fakelcdir )
# get the magrms relation needed from the pickle or input dict
if isinstance ( magrmsfrom , str ) and os . path . exists ( magrmsfrom ) :
    with open ( magrmsfrom , 'rb' ) as infd :
        xmagrms = pickle . load ( infd )
elif isinstance ( magrmsfrom , dict ) :
    xmagrms = magrmsfrom
magrms = { }
# get the required items from the magrms dict . interpolate the mag - rms
# relation for the magcol so the make _ fake _ lc function can use it directly .
for magcol in magcols :
    if ( magcol in xmagrms and 'binned_sdssr_median' in xmagrms [ magcol ] and 'binned_lcmad_median' in xmagrms [ magcol ] ) :
        magrms [ magcol ] = { 'binned_sdssr_median' : np . array ( xmagrms [ magcol ] [ 'binned_sdssr_median' ] ) , 'binned_lcmad_median' : np . array ( xmagrms [ magcol ] [ 'binned_lcmad_median' ] ) , }
        # interpolate the mag - MAD relation
        interpolated_magmad = spi . interp1d ( xmagrms [ magcol ] [ 'binned_sdssr_median' ] , xmagrms [ magcol ] [ 'binned_lcmad_median' ] , kind = magrms_interpolate , fill_value = magrms_fillvalue , )
        # save the magrms
        magrms [ magcol ] [ 'interpolated_magmad' ] = interpolated_magmad
        # generate the probability distribution in magbins . this is needed
        # to correctly sample the objects in this population
        magbins = np . array ( xmagrms [ magcol ] [ 'binned_sdssr_median' ] )
        bincounts = np . array ( xmagrms [ magcol ] [ 'binned_count' ] )
        binprobs = bincounts / np . sum ( bincounts )
        # save the bin probabilities as well
        magrms [ magcol ] [ 'magbin_probabilities' ] = binprobs
    else :
        LOGWARNING ( 'input magrms dict does not have ' 'required info for magcol: %s' % magcol )
        magrms [ magcol ] = { 'binned_sdssr_median' : None , 'binned_lcmad_median' : None , 'interpolated_magmad' : None , 'magbin_probabilities' : None , }
tasks = [ ( x , fakelcdir , { 'lcformat' : lcformat , 'timecols' : timecols , 'magcols' : magcols , 'errcols' : errcols , 'magrms' : magrms , 'randomizemags' : randomizemags , 'randomizecoords' : randomizecoords } ) for x in chosenlcs ]
# we can ' t parallelize because it messes up the random number generation ,
# causing all the IDs to clash . FIXME : figure out a way around this
# ( probably initial a seed in each worker process ? )
fakeresults = [ collection_worker ( task ) for task in tasks ]
fakedb = { 'simbasedir' : simbasedir , 'lcformat' : lcformat , 'timecols' : timecols , 'magcols' : magcols , 'errcols' : errcols , 'magsarefluxes' : magsarefluxes }
fobjects , fpaths = [ ] , [ ]
fras , fdecls , fndets = [ ] , [ ] , [ ]
fmags , fmagmads = [ ] , [ ]
ferrmeds , ferrmads = [ ] , [ ]
totalvars = 0
# these are the indices for the variable objects chosen randomly
isvariableind = npr . randint ( 0 , high = len ( fakeresults ) , size = maxvars )
isvariable = np . full ( len ( fakeresults ) , False , dtype = np . bool )
isvariable [ isvariableind ] = True
fakedb [ 'isvariable' ] = isvariable
LOGINFO ( 'added %s variable stars' % maxvars )
# these are the variable types for each variable object
vartypeind = npr . randint ( 0 , high = len ( vartypes ) , size = maxvars )
vartypearr = np . array ( [ vartypes [ x ] for x in vartypeind ] )
fakedb [ 'vartype' ] = vartypearr
for vt in sorted ( vartypes ) :
    LOGINFO ( '%s: %s stars' % ( vt , vartypearr [ vartypearr == vt ] . size ) )
# now go through the collection and get the mag / rms and err / rms for each
# star . these will be used later to add noise to light curves
LOGINFO ( 'collecting info...' )
for fr in fakeresults :
    if fr is not None :
        fpath , fcols , finfo , fmoments = fr
        fobjects . append ( finfo [ 'objectid' ] )
        fpaths . append ( fpath )
        fras . append ( finfo [ 'ra' ] )
        fdecls . append ( finfo [ 'decl' ] )
        fndets . append ( finfo [ 'ndet' ] )
        fmags . append ( finfo [ 'sdssr' ] )
        # this is per magcol
        fmagmads . append ( [ fmoments [ x ] [ 'mad' ] for x in magcols ] )
        # these are per errcol
        ferrmeds . append ( [ fmoments [ x ] [ 'median' ] for x in errcols ] )
        ferrmads . append ( [ fmoments [ x ] [ 'mad' ] for x in errcols ] )
# convert to nparrays
fobjects = np . array ( fobjects )
fpaths = np . array ( fpaths )
fras = np . array ( fras )
fdecls = np . array ( fdecls )
fndets = np . array ( fndets )
fmags = np . array ( fmags )
fmagmads = np . array ( fmagmads )
ferrmeds = np . array ( ferrmeds )
ferrmads = np . array ( ferrmads )
# put these in the fakedb
fakedb [ 'objectid' ] = fobjects
fakedb [ 'lcfpath' ] = fpaths
fakedb [ 'ra' ] = fras
fakedb [ 'decl' ] = fdecls
fakedb [ 'ndet' ] = fndets
fakedb [ 'sdssr' ] = fmags
fakedb [ 'mad' ] = fmagmads
fakedb [ 'errmedian' ] = ferrmeds
fakedb [ 'errmad' ] = ferrmads
# get the mag - RMS curve for this light curve collection for each magcol
fakedb [ 'magrms' ] = magrms
# finally , write the collection DB to a pickle in simbasedir
dboutfname = os . path . join ( simbasedir , 'fakelcs-info.pkl' )
with open ( dboutfname , 'wb' ) as outfd :
    pickle . dump ( fakedb , outfd )
LOGINFO ( 'wrote %s fake LCs to: %s' % ( len ( fakeresults ) , simbasedir ) )
LOGINFO ( 'fake LC info written to: %s' % dboutfname )
return dboutfname


signature:
def add_fakelc_variability ( fakelcfile , vartype , override_paramdists = None , magsarefluxes = False , overwrite = False ) :
    '''This adds variability of the specified type to the fake LC .
 The procedure is ( for each ` magcol ` ) :
 - read the fakelcfile , get the stored moments and vartype info
 - add the periodic variability specified in vartype and varparamdists . if
 ` vartype = = None ` , then do nothing in this step . If ` override _ vartype ` is
 not None , override stored vartype with specified vartype . If
 ` override _ varparamdists ` provided , override with specified
 ` varparamdists ` . NOTE : the varparamdists must make sense for the vartype ,
 otherwise , weird stuff will happen .
 - add the median mag level stored in ` fakelcfile ` to the time series
 - add Gaussian noise to the light curve as specified in ` fakelcfile `
 - add a varinfo key and dict to the lcdict with ` varperiod ` , ` varepoch ` ,
 ` varparams `
 - write back to fake LC pickle
 - return the ` varinfo ` dict to the caller
 Parameters
 fakelcfile : str
 The name of the fake LC file to process .
 vartype : str
 The type of variability to add to this fake LC file .
 override _ paramdists : dict
 A parameter distribution dict as in the ` generate _ XX _ lightcurve `
 functions above . If provided , will override the distribution stored in
 the input fake LC file itself .
 magsarefluxes : bool
 Sets if the variability amplitude is in fluxes and not magnitudes .
 overwite : bool
 This overwrites the input fake LC file with a new variable LC even if
 it ' s been processed before .
 Returns
 dict
 A dict of the following form is returned : :
 { ' objectid ' : lcdict [ ' objectid ' ] ,
 ' lcfname ' : fakelcfile ,
 ' actual _ vartype ' : vartype ,
 ' actual _ varparams ' : lcdict [ ' actual _ varparams ' ] }'''

implementation:
# read in the fakelcfile
lcdict = _read_pklc ( fakelcfile )
# make sure to bail out if this light curve already has fake variability
# added
if ( 'actual_vartype' in lcdict and 'actual_varparams' in lcdict and not overwrite ) :
    LOGERROR ( '%s has existing variability type: %s ' 'and params: %s and overwrite = False, ' 'skipping this file...' % ( fakelcfile , lcdict [ 'actual_vartype' ] , repr ( lcdict [ 'actual_varparams' ] ) ) )
    return None
# get the times , mags , errs from this LC
timecols , magcols , errcols = ( lcdict [ 'timecols' ] , lcdict [ 'magcols' ] , lcdict [ 'errcols' ] )
# get the correct function to apply variability
if vartype in VARTYPE_LCGEN_MAP :
    vargenfunc = VARTYPE_LCGEN_MAP [ vartype ]
elif vartype is None :
    vargenfunc = None
else :
    LOGERROR ( 'unknown variability type: %s, choose from: %s' % ( vartype , repr ( list ( VARTYPE_LCGEN_MAP . keys ( ) ) ) ) )
    return None
# 1 . generate the variability , including the overrides if provided we do
# this outside the loop below to get the period , etc . distributions once
# only per object . NOTE : in doing so , we ' re assuming that the difference
# between magcols is just additive and the timebases for each magcol are the
# same ; this is not strictly correct
if vargenfunc is not None :
    if ( override_paramdists is not None and isinstance ( override_paramdists , dict ) ) :
        variablelc = vargenfunc ( lcdict [ timecols [ 0 ] ] , paramdists = override_paramdists , magsarefluxes = magsarefluxes )
    else :
        variablelc = vargenfunc ( lcdict [ timecols [ 0 ] ] , magsarefluxes = magsarefluxes )
# for nonvariables , don ' t execute vargenfunc , but return a similar dict
# so we can add the required noise to it
else :
    variablelc = { 'vartype' : None , 'params' : None , 'times' : lcdict [ timecols [ 0 ] ] , 'mags' : np . full_like ( lcdict [ timecols [ 0 ] ] , 0.0 ) , 'errs' : np . full_like ( lcdict [ timecols [ 0 ] ] , 0.0 ) }
# now iterate over the time , mag , err columns
for tcol , mcol , ecol in zip ( timecols , magcols , errcols ) :
    times , mags , errs = lcdict [ tcol ] , lcdict [ mcol ] , lcdict [ ecol ]
    # 2 . get the moments for this magcol
    mag_median = lcdict [ 'moments' ] [ mcol ] [ 'median' ]
    mag_mad = lcdict [ 'moments' ] [ mcol ] [ 'mad' ]
    # add up to 5 mmag of extra RMS for systematics and red - noise
    mag_rms = mag_mad * 1.483
    err_median = lcdict [ 'moments' ] [ ecol ] [ 'median' ]
    err_mad = lcdict [ 'moments' ] [ ecol ] [ 'mad' ]
    err_rms = err_mad * 1.483
    # 3 . add the median level + gaussian noise
    magnoise = npr . normal ( size = variablelc [ 'mags' ] . size ) * mag_rms
    errnoise = npr . normal ( size = variablelc [ 'errs' ] . size ) * err_rms
    finalmags = mag_median + ( variablelc [ 'mags' ] + magnoise )
    finalerrs = err_median + ( variablelc [ 'errs' ] + errnoise )
    # 4 . update these tcol , mcol , ecol values in the lcdict
    lcdict [ mcol ] = finalmags
    lcdict [ ecol ] = finalerrs
# all done with updating mags and errs
# 5 . update the light curve with the variability info
lcdict [ 'actual_vartype' ] = variablelc [ 'vartype' ]
lcdict [ 'actual_varparams' ] = variablelc [ 'params' ]
# these standard keys are set to help out later with characterizing recovery
# rates by magnitude , period , amplitude , ndet , etc .
if vartype is not None :
    lcdict [ 'actual_varperiod' ] = variablelc [ 'varperiod' ]
    lcdict [ 'actual_varamplitude' ] = variablelc [ 'varamplitude' ]
else :
    lcdict [ 'actual_varperiod' ] = np . nan
    lcdict [ 'actual_varamplitude' ] = np . nan
# 6 . write back , making sure to do it safely
tempoutf = '%s.%s' % ( fakelcfile , md5 ( npr . bytes ( 4 ) ) . hexdigest ( ) [ - 8 : ] )
with open ( tempoutf , 'wb' ) as outfd :
    pickle . dump ( lcdict , outfd , pickle . HIGHEST_PROTOCOL )
if os . path . exists ( tempoutf ) :
    shutil . copy ( tempoutf , fakelcfile )
    os . remove ( tempoutf )
else :
    LOGEXCEPTION ( 'could not write output light curve file to dir: %s' % os . path . dirname ( tempoutf ) )
    # fail here
    raise
LOGINFO ( 'object: %s, vartype: %s -> %s OK' % ( lcdict [ 'objectid' ] , vartype , fakelcfile ) )
return { 'objectid' : lcdict [ 'objectid' ] , 'lcfname' : fakelcfile , 'actual_vartype' : vartype , 'actual_varparams' : lcdict [ 'actual_varparams' ] }


signature:
def add_variability_to_fakelc_collection ( simbasedir , override_paramdists = None , overwrite_existingvar = False ) :
    '''This adds variability and noise to all fake LCs in ` simbasedir ` .
 If an object is marked as variable in the ` fakelcs - info ` . pkl file in
 ` simbasedir ` , a variable signal will be added to its light curve based on
 its selected type , default period and amplitude distribution , the
 appropriate params , etc . the epochs for each variable object will be chosen
 uniformly from its time - range ( and may not necessarily fall on a actual
 observed time ) . Nonvariable objects will only have noise added as determined
 by their params , but no variable signal will be added .
 Parameters
 simbasedir : str
 The directory containing the fake LCs to process .
 override _ paramdists : dict
 This can be used to override the stored variable parameters in each fake
 LC . It should be a dict of the following form : :
 { ' < vartype1 > ' : { ' < param1 > : a scipy . stats distribution function or
 the np . random . randint function ,
 ' < paramN > : a scipy . stats distribution function
 or the np . random . randint function }
 for any vartype in VARTYPE _ LCGEN _ MAP . These are used to override the
 default parameter distributions for each variable type .
 overwrite _ existingvar : bool
 If this is True , then will overwrite any existing variability in the
 input fake LCs in ` simbasedir ` .
 Returns
 dict
 This returns a dict containing the fake LC filenames as keys and
 variability info for each as values .'''

implementation:
# open the fakelcs - info . pkl
infof = os . path . join ( simbasedir , 'fakelcs-info.pkl' )
with open ( infof , 'rb' ) as infd :
    lcinfo = pickle . load ( infd )
lclist = lcinfo [ 'lcfpath' ]
varflag = lcinfo [ 'isvariable' ]
vartypes = lcinfo [ 'vartype' ]
vartind = 0
varinfo = { }
# go through all the LCs and add the required type of variability
for lc , varf , _lcind in zip ( lclist , varflag , range ( len ( lclist ) ) ) : # if this object is variable , add variability
    if varf :
        thisvartype = vartypes [ vartind ]
        if ( override_paramdists and isinstance ( override_paramdists , dict ) and thisvartype in override_paramdists and isinstance ( override_paramdists [ thisvartype ] , dict ) ) :
            thisoverride_paramdists = override_paramdists [ thisvartype ]
        else :
            thisoverride_paramdists = None
        varlc = add_fakelc_variability ( lc , thisvartype , override_paramdists = thisoverride_paramdists , overwrite = overwrite_existingvar )
        varinfo [ varlc [ 'objectid' ] ] = { 'params' : varlc [ 'actual_varparams' ] , 'vartype' : varlc [ 'actual_vartype' ] }
        # update vartind
        vartind = vartind + 1
    else :
        varlc = add_fakelc_variability ( lc , None , overwrite = overwrite_existingvar )
        varinfo [ varlc [ 'objectid' ] ] = { 'params' : varlc [ 'actual_varparams' ] , 'vartype' : varlc [ 'actual_vartype' ] }
# done with all objects
# write the varinfo back to the dict and fakelcs - info . pkl
lcinfo [ 'varinfo' ] = varinfo
tempoutf = '%s.%s' % ( infof , md5 ( npr . bytes ( 4 ) ) . hexdigest ( ) [ - 8 : ] )
with open ( tempoutf , 'wb' ) as outfd :
    pickle . dump ( lcinfo , outfd , pickle . HIGHEST_PROTOCOL )
if os . path . exists ( tempoutf ) :
    shutil . copy ( tempoutf , infof )
    os . remove ( tempoutf )
else :
    LOGEXCEPTION ( 'could not write output light curve file to dir: %s' % os . path . dirname ( tempoutf ) )
    # fail here
    raise
return lcinfo


signature:
def add_flare_model ( flareparams , times , mags , errs ) :
    '''This adds a flare model function to the input magnitude / flux time - series .
 Parameters
 flareparams : list of float
 This defines the flare model : :
 [ amplitude ,
 flare _ peak _ time ,
 rise _ gaussian _ stdev ,
 decay _ time _ constant ]
 where :
 ` amplitude ` : the maximum flare amplitude in mags or flux . If flux , then
 amplitude should be positive . If mags , amplitude should be negative .
 ` flare _ peak _ time ` : time at which the flare maximum happens .
 ` rise _ gaussian _ stdev ` : the stdev of the gaussian describing the rise of
 the flare .
 ` decay _ time _ constant ` : the time constant of the exponential fall of the
 flare .
 times , mags , errs : np . array
 The input time - series of measurements and associated errors for which
 the model will be generated . The times will be used to generate
 model mags .
 magsarefluxes : bool
 Sets the correct direction of the flare amplitude ( + ve ) for fluxes if
 True and for mags ( - ve ) if False .
 Returns
 dict
 A dict of the form below is returned : :
 { ' times ' : the original times array
 ' mags ' : the original mags + the flare model mags evaluated at times ,
 ' errs ' : the original errs array ,
 ' flareparams ' : the input list of flare params }'''

implementation:
modelmags , ftimes , fmags , ferrs = flares . flare_model ( flareparams , times , mags , errs )
return { 'times' : times , 'mags' : mags + modelmags , 'errs' : errs , 'flareparams' : flareparams }


signature:
def simple_flare_find ( times , mags , errs , smoothbinsize = 97 , flare_minsigma = 4.0 , flare_maxcadencediff = 1 , flare_mincadencepoints = 3 , magsarefluxes = False , savgol_polyorder = 2 , ** savgol_kwargs ) :
    '''This finds flares in time series using the method in Walkowicz + 2011.
 FIXME : finish this .
 Parameters
 times , mags , errs : np . array
 The input time - series to find flares in .
 smoothbinsize : int
 The number of consecutive light curve points to smooth over in the time
 series using a Savitsky - Golay filter . The smoothed light curve is then
 subtracted from the actual light curve to remove trends that potentially
 last ` smoothbinsize ` light curve points . The default value is chosen as
 ~ 6.5 hours ( 97 x 4 minute cadence for HATNet / HATSouth ) .
 flare _ minsigma : float
 The minimum sigma above the median LC level to designate points as
 belonging to possible flares .
 flare _ maxcadencediff : int
 The maximum number of light curve points apart each possible flare event
 measurement is allowed to be . If this is 1 , then we ' ll look for
 consecutive measurements .
 flare _ mincadencepoints : int
 The minimum number of light curve points ( each ` flare _ maxcadencediff `
 points apart ) required that are at least ` flare _ minsigma ` above the
 median light curve level to call an event a flare .
 magsarefluxes : bool
 If True , indicates that mags is actually an array of fluxes .
 savgol _ polyorder : int
 The polynomial order of the function used by the Savitsky - Golay filter .
 savgol _ kwargs : extra kwargs
 Any remaining keyword arguments are passed directly to the
 ` savgol _ filter ` function from ` scipy . signal ` .
 Returns
 ( nflares , flare _ indices ) : tuple
 Returns the total number of flares found and their time - indices ( start ,
 end ) as tuples .'''

implementation:
# if no errs are given , assume 0.1 % errors
if errs is None :
    errs = 0.001 * mags
# get rid of nans first
finiteind = np . isfinite ( times ) & np . isfinite ( mags ) & np . isfinite ( errs )
ftimes = times [ finiteind ]
fmags = mags [ finiteind ]
ferrs = errs [ finiteind ]
# now get the smoothed mag series using the filter
# kwargs are provided to the savgol _ filter function
smoothed = savgol_filter ( fmags , smoothbinsize , savgol_polyorder , ** savgol_kwargs )
subtracted = fmags - smoothed
# calculate some stats
# the series _ median is ~ zero after subtraction
series_mad = np . median ( np . abs ( subtracted ) )
series_stdev = 1.483 * series_mad
# find extreme positive deviations
if magsarefluxes :
    extind = np . where ( subtracted > ( flare_minsigma * series_stdev ) )
else :
    extind = np . where ( subtracted < ( - flare_minsigma * series_stdev ) )
# see if there are any extrema
if extind and extind [ 0 ] :
    extrema_indices = extind [ 0 ]
    flaregroups = [ ]
    # find the deviations within the requested flaremaxcadencediff
    for ind , extrema_index in enumerate ( extrema_indices ) : # FIXME : finish this
        pass


signature:
def _smooth_acf ( acf , windowfwhm = 7 , windowsize = 21 ) :
    '''This returns a smoothed version of the ACF .
 Convolves the ACF with a Gaussian of given ` windowsize ` and ` windowfwhm ` .
 Parameters
 acf : np . array
 The auto - correlation function array to smooth .
 windowfwhm : int
 The smoothing window Gaussian kernel ' s FWHM .
 windowsize : int
 The number of input points to apply the smoothing over .
 Returns
 np . array
 Smoothed version of the input ACF array .'''

implementation:
convkernel = Gaussian1DKernel ( windowfwhm , x_size = windowsize )
smoothed = convolve ( acf , convkernel , boundary = 'extend' )
return smoothed


signature:
def _smooth_acf_savgol ( acf , windowsize = 21 , polyorder = 2 ) :
    '''This returns a smoothed version of the ACF .
 This version uses the Savitsky - Golay smoothing filter .
 Parameters
 acf : np . array
 The auto - correlation function array to smooth .
 windowsize : int
 The number of input points to apply the smoothing over .
 polyorder : int
 The order of the polynomial to use in the Savitsky - Golay filter .
 Returns
 np . array
 Smoothed version of the input ACF array .'''

implementation:
smoothed = savgol_filter ( acf , windowsize , polyorder )
return smoothed


signature:
def _get_acf_peakheights ( lags , acf , npeaks = 20 , searchinterval = 1 ) :
    '''This calculates the relative peak heights for first npeaks in ACF .
 Usually , the first peak or the second peak ( if its peak height > first peak )
 corresponds to the correct lag . When we know the correct lag , the period is
 then : :
 bestperiod = time [ lags = = bestlag ] - time [ 0]
 Parameters
 lags : np . array
 An array of lags that the ACF is calculated at .
 acf : np . array
 The array containing the ACF values .
 npeaks : int
 THe maximum number of peaks to consider when finding peak heights .
 searchinterval : int
 From ` scipy . signal . argrelmax ` : " How many points on each side to use for
 the comparison to consider comparator ( n , n + x ) to be True . " This
 effectively sets how many points on each of the current peak will be
 used to check if the current peak is the local maximum .
 Returns
 dict
 This returns a dict of the following form : :
 { ' maxinds ' : the indices of the lag array where maxes are ,
 ' maxacfs ' : the ACF values at each max ,
 ' maxlags ' : the lag values at each max ,
 ' mininds ' : the indices of the lag array where mins are ,
 ' minacfs ' : the ACF values at each min ,
 ' minlags ' : the lag values at each min ,
 ' relpeakheights ' : the relative peak heights of each rel . ACF peak ,
 ' relpeaklags ' : the lags at each rel . ACF peak found ,
 ' peakindices ' : the indices of arrays where each rel . ACF peak is ,
 ' bestlag ' : the lag value with the largest rel . ACF peak height ,
 ' bestpeakheight ' : the largest rel . ACF peak height ,
 ' bestpeakindex ' : the largest rel . ACF peak ' s number in all peaks }'''

implementation:
maxinds = argrelmax ( acf , order = searchinterval ) [ 0 ]
maxacfs = acf [ maxinds ]
maxlags = lags [ maxinds ]
mininds = argrelmin ( acf , order = searchinterval ) [ 0 ]
minacfs = acf [ mininds ]
minlags = lags [ mininds ]
relpeakheights = npzeros ( npeaks )
relpeaklags = npzeros ( npeaks , dtype = npint64 )
peakindices = npzeros ( npeaks , dtype = npint64 )
for peakind , mxi in enumerate ( maxinds [ : npeaks ] ) : # check if there are no mins to the left
# throw away this peak because it ' s probably spurious
# ( FIXME : is this OK ? )
    if npall ( mxi < mininds ) :
        continue
    leftminind = mininds [ mininds < mxi ] [ - 1 ]
    # the last index to the left
    rightminind = mininds [ mininds > mxi ] [ 0 ]
    # the first index to the right
    relpeakheights [ peakind ] = ( acf [ mxi ] - ( acf [ leftminind ] + acf [ rightminind ] ) / 2.0 )
    relpeaklags [ peakind ] = lags [ mxi ]
    peakindices [ peakind ] = peakind
# figure out the bestperiod if possible
if relpeakheights [ 0 ] > relpeakheights [ 1 ] :
    bestlag = relpeaklags [ 0 ]
    bestpeakheight = relpeakheights [ 0 ]
    bestpeakindex = peakindices [ 0 ]
else :
    bestlag = relpeaklags [ 1 ]
    bestpeakheight = relpeakheights [ 1 ]
    bestpeakindex = peakindices [ 1 ]
return { 'maxinds' : maxinds , 'maxacfs' : maxacfs , 'maxlags' : maxlags , 'mininds' : mininds , 'minacfs' : minacfs , 'minlags' : minlags , 'relpeakheights' : relpeakheights , 'relpeaklags' : relpeaklags , 'peakindices' : peakindices , 'bestlag' : bestlag , 'bestpeakheight' : bestpeakheight , 'bestpeakindex' : bestpeakindex }


signature:
def plot_acf_results ( acfp , outfile , maxlags = 5000 , yrange = ( - 0.4 , 0.4 ) ) :
    '''This plots the unsmoothed / smoothed ACF vs lag .
 Parameters
 acfp : dict
 This is the dict returned from ` macf _ period _ find ` below .
 outfile : str
 The output file the plot will be written to .
 maxlags : int
 The maximum number of lags to include in the plot .
 yrange : sequence of two floats
 The y - range of the ACF vs . lag plot to use .'''

implementation:
import matplotlib . pyplot as plt
lags = acfp [ 'acfresults' ] [ 'lags' ] [ : maxlags ]
smoothedacf = acfp [ 'acf' ] [ : maxlags ]
unsmoothedacf = acfp [ 'acfresults' ] [ 'acf' ] [ : maxlags ]
acfparams = acfp [ 'kwargs' ] [ 'smoothfunckwargs' ] . copy ( )
acfparams . update ( { 'peakinterval' : int ( acfp [ 'kwargs' ] [ 'smoothacf' ] / 2.0 ) } )
# plot the ACFs
fig , ax1 = plt . subplots ( )
# this is lags vs acf
ax1 . plot ( lags , unsmoothedacf , label = 'unsmoothed ACF' , color = '#1f77b4' )
ax1 . plot ( lags , smoothedacf , label = 'smoothed ACF' , color = '#ff7f0e' )
ax1 . set_xlim ( ( 0 , maxlags ) )
ax1 . set_xlabel ( 'lags' )
# overplot the identified peaks
acfmaxinds = acfp [ 'acfpeaks' ] [ 'maxinds' ]
for i , maxind in enumerate ( acfmaxinds ) :
    if i == 0 :
        ax1 . axvline ( maxind , linewidth = 2.0 , color = 'red' , ymin = 0.2 , ymax = 0.3 , label = 'identified ACF peaks' )
    else :
        ax1 . axvline ( maxind , linewidth = 2.0 , color = 'red' , ymin = 0.2 , ymax = 0.3 )
plt . ylabel ( 'ACF' )
plt . ylim ( yrange )
ax1 . legend ( )
plt . title ( '%s' % repr ( acfparams ) )
plt . tight_layout ( )
plt . savefig ( outfile )
plt . close ( 'all' )
return outfile


signature:
def macf_period_find ( times , mags , errs , fillgaps = 0.0 , filterwindow = 11 , forcetimebin = None , maxlags = None , maxacfpeaks = 10 , smoothacf = 21 , # set for Kepler - type LCs , see details below
smoothfunc = _smooth_acf_savgol , smoothfunckwargs = None , magsarefluxes = False , sigclip = 3.0 , verbose = True , periodepsilon = 0.1 , # doesn ' t do anything , for consistent external API
nworkers = None , # doesn ' t do anything , for consistent external API
startp = None , # doesn ' t do anything , for consistent external API
endp = None , # doesn ' t do anything , for consistent external API
autofreq = None , # doesn ' t do anything , for consistent external API
stepsize = None , # doesn ' t do anything , for consistent external API
) :
    '''This finds periods using the McQuillan + ( 2013a , 2014 ) ACF method .
 The kwargs from ` periodepsilon ` to ` stepsize ` don ' t do anything but are used
 to present a consistent API for all periodbase period - finders to an outside
 driver ( e . g . the one in the checkplotserver ) .
 Parameters
 times , mags , errs : np . array
 The input magnitude / flux time - series to run the period - finding for .
 fillgaps : ' noiselevel ' or float
 This sets what to use to fill in gaps in the time series . If this is
 ' noiselevel ' , will smooth the light curve using a point window size of
 ` filterwindow ` ( this should be an odd integer ) , subtract the smoothed LC
 from the actual LC and estimate the RMS . This RMS will be used to fill
 in the gaps . Other useful values here are 0.0 , and npnan .
 filterwindow : int
 The light curve ' s smoothing filter window size to use if
 ` fillgaps = ' noiselevel ` ' .
 forcetimebin : None or float
 This is used to force a particular cadence in the light curve other than
 the automatically determined cadence . This effectively rebins the light
 curve to this cadence . This should be in the same time units as ` times ` .
 maxlags : None or int
 This is the maximum number of lags to calculate . If None , will calculate
 all lags .
 maxacfpeaks : int
 This is the maximum number of ACF peaks to use when finding the highest
 peak and obtaining a fit period .
 smoothacf : int
 This is the number of points to use as the window size when smoothing
 the ACF with the ` smoothfunc ` . This should be an odd integer value . If
 this is None , will not smooth the ACF , but this will probably lead to
 finding spurious peaks in a generally noisy ACF .
 For Kepler , a value between 21 and 51 seems to work fine . For ground
 based data , much larger values may be necessary : between 1001 and 2001
 seem to work best for the HAT surveys . This is dependent on cadence , RMS
 of the light curve , the periods of the objects you ' re looking for , and
 finally , any correlated noise in the light curve . Make a plot of the
 smoothed / unsmoothed ACF vs . lag using the result dict of this function
 and the ` plot _ acf _ results ` function above to see the identified ACF
 peaks and what kind of smoothing might be needed .
 The value of ` smoothacf ` will also be used to figure out the interval to
 use when searching for local peaks in the ACF : this interval is 1/2 of
 the ` smoothacf ` value .
 smoothfunc : Python function
 This is the function that will be used to smooth the ACF . This should
 take at least one kwarg : ' windowsize ' . Other kwargs can be passed in
 using a dict provided in ` smoothfunckwargs ` . By default , this uses a
 Savitsky - Golay filter , a Gaussian filter is also provided but not
 used . Another good option would be an actual low - pass filter ( generated
 using scipy . signal ? ) to remove all high frequency noise from the ACF .
 smoothfunckwargs : dict or None
 The dict of optional kwargs to pass in to the ` smoothfunc ` .
 magsarefluxes : bool
 If your input measurements in ` mags ` are actually fluxes instead of
 mags , set this is True .
 sigclip : float or int or sequence of two floats / ints or None
 If a single float or int , a symmetric sigma - clip will be performed using
 the number provided as the sigma - multiplier to cut out from the input
 time - series .
 If a list of two ints / floats is provided , the function will perform an
 ' asymmetric ' sigma - clip . The first element in this list is the sigma
 value to use for fainter flux / mag values ; the second element in this
 list is the sigma value to use for brighter flux / mag values . For
 example , ` sigclip = [ 10 . , 3 . ] ` , will sigclip out greater than 10 - sigma
 dimmings and greater than 3 - sigma brightenings . Here the meaning of
 " dimming " and " brightening " is set by * physics * ( not the magnitude
 system ) , which is why the ` magsarefluxes ` kwarg must be correctly set .
 If ` sigclip ` is None , no sigma - clipping will be performed , and the
 time - series ( with non - finite elems removed ) will be passed through to
 the output .
 verbose : bool
 If True , will indicate progress and report errors .
 Returns
 dict
 Returns a dict with results . dict [ ' bestperiod ' ] is the estimated best
 period and dict [ ' fitperiodrms ' ] is its estimated error . Other
 interesting things in the output include :
 - dict [ ' acfresults ' ] : all results from calculating the ACF . in
 particular , the unsmoothed ACF might be of interest :
 dict [ ' acfresults ' ] [ ' acf ' ] and dict [ ' acfresults ' ] [ ' lags ' ] .
 - dict [ ' lags ' ] and dict [ ' acf ' ] contain the ACF after smoothing was
 applied .
 - dict [ ' periods ' ] and dict [ ' lspvals ' ] can be used to construct a
 pseudo - periodogram .
 - dict [ ' naivebestperiod ' ] is obtained by multiplying the lag at the
 highest ACF peak with the cadence . This is usually close to the fit
 period ( dict [ ' fitbestperiod ' ] ) , which is calculated by doing a fit to
 the lags vs . peak index relation as in McQuillan + 2014.'''

implementation:
# get the ACF
acfres = autocorr_magseries ( times , mags , errs , maxlags = maxlags , fillgaps = fillgaps , forcetimebin = forcetimebin , sigclip = sigclip , magsarefluxes = magsarefluxes , filterwindow = filterwindow , verbose = verbose )
xlags = acfres [ 'lags' ]
# smooth the ACF if requested
if smoothacf and isinstance ( smoothacf , int ) and smoothacf > 0 :
    if smoothfunckwargs is None :
        sfkwargs = { 'windowsize' : smoothacf }
    else :
        sfkwargs = smoothfunckwargs . copy ( )
        sfkwargs . update ( { 'windowsize' : smoothacf } )
    xacf = smoothfunc ( acfres [ 'acf' ] , ** sfkwargs )
else :
    xacf = acfres [ 'acf' ]
# get the relative peak heights and fit best lag
peakres = _get_acf_peakheights ( xlags , xacf , npeaks = maxacfpeaks , searchinterval = int ( smoothacf / 2 ) )
# this is the best period ' s best ACF peak height
bestlspval = peakres [ 'bestpeakheight' ]
try : # get the fit best lag from a linear fit to the peak index vs time ( peak
# lag ) function as in McQillian + ( 2014)
    fity = npconcatenate ( ( [ 0.0 , peakres [ 'bestlag' ] ] , peakres [ 'relpeaklags' ] [ peakres [ 'relpeaklags' ] > peakres [ 'bestlag' ] ] ) )
    fity = fity * acfres [ 'cadence' ]
    fitx = nparange ( fity . size )
    fitcoeffs , fitcovar = nppolyfit ( fitx , fity , 1 , cov = True )
    # fit best period is the gradient of fit
    fitbestperiod = fitcoeffs [ 0 ]
    bestperiodrms = npsqrt ( fitcovar [ 0 , 0 ] )
# from the covariance matrix
except Exception as e :
    LOGWARNING ( 'linear fit to time at each peak lag ' 'value vs. peak number failed, ' 'naively calculated ACF period may not be accurate' )
    fitcoeffs = nparray ( [ npnan , npnan ] )
    fitcovar = nparray ( [ [ npnan , npnan ] , [ npnan , npnan ] ] )
    fitbestperiod = npnan
    bestperiodrms = npnan
    raise
# calculate the naive best period using delta _ tau = lag * cadence
naivebestperiod = peakres [ 'bestlag' ] * acfres [ 'cadence' ]
if fitbestperiod < naivebestperiod :
    LOGWARNING ( 'fit bestperiod = %.5f may be an alias, ' 'naively calculated bestperiod is = %.5f' % ( fitbestperiod , naivebestperiod ) )
if npisfinite ( fitbestperiod ) :
    bestperiod = fitbestperiod
else :
    bestperiod = naivebestperiod
return { 'bestperiod' : bestperiod , 'bestlspval' : bestlspval , 'nbestpeaks' : maxacfpeaks , # for compliance with the common pfmethod API
'nbestperiods' : npconcatenate ( [ [ fitbestperiod ] , peakres [ 'relpeaklags' ] [ 1 : maxacfpeaks ] * acfres [ 'cadence' ] ] ) , 'nbestlspvals' : peakres [ 'maxacfs' ] [ : maxacfpeaks ] , 'lspvals' : xacf , 'periods' : xlags * acfres [ 'cadence' ] , 'acf' : xacf , 'lags' : xlags , 'method' : 'acf' , 'naivebestperiod' : naivebestperiod , 'fitbestperiod' : fitbestperiod , 'fitperiodrms' : bestperiodrms , 'periodfitcoeffs' : fitcoeffs , 'periodfitcovar' : fitcovar , 'kwargs' : { 'maxlags' : maxlags , 'maxacfpeaks' : maxacfpeaks , 'fillgaps' : fillgaps , 'filterwindow' : filterwindow , 'smoothacf' : smoothacf , 'smoothfunckwargs' : sfkwargs , 'magsarefluxes' : magsarefluxes , 'sigclip' : sigclip } , 'acfresults' : acfres , 'acfpeaks' : peakres }


signature:
def _autocorr_func1 ( mags , lag , maglen , magmed , magstd ) :
    '''Calculates the autocorr of mag series for specific lag .
 This version of the function is taken from : Kim et al . ( ` 2011
 < https : / / dx . doi . org / 10.1088/0004-637X / 735/2/68 > ` _ )
 Parameters
 mags : np . array
 This is the magnitudes array . MUST NOT have any nans .
 lag : float
 The specific lag value to calculate the auto - correlation for . This MUST
 be less than total number of observations in ` mags ` .
 maglen : int
 The number of elements in the ` mags ` array .
 magmed : float
 The median of the ` mags ` array .
 magstd : float
 The standard deviation of the ` mags ` array .
 Returns
 float
 The auto - correlation at this specific ` lag ` value .'''

implementation:
lagindex = nparange ( 1 , maglen - lag )
products = ( mags [ lagindex ] - magmed ) * ( mags [ lagindex + lag ] - magmed )
acorr = ( 1.0 / ( ( maglen - lag ) * magstd ) ) * npsum ( products )
return acorr


signature:
def _autocorr_func2 ( mags , lag , maglen , magmed , magstd ) :
    '''This is an alternative function to calculate the autocorrelation .
 This version is from ( first definition ) :
 https : / / en . wikipedia . org / wiki / Correlogram # Estimation _ of _ autocorrelations
 Parameters
 mags : np . array
 This is the magnitudes array . MUST NOT have any nans .
 lag : float
 The specific lag value to calculate the auto - correlation for . This MUST
 be less than total number of observations in ` mags ` .
 maglen : int
 The number of elements in the ` mags ` array .
 magmed : float
 The median of the ` mags ` array .
 magstd : float
 The standard deviation of the ` mags ` array .
 Returns
 float
 The auto - correlation at this specific ` lag ` value .'''

implementation:
lagindex = nparange ( 0 , maglen - lag )
products = ( mags [ lagindex ] - magmed ) * ( mags [ lagindex + lag ] - magmed )
autocovarfunc = npsum ( products ) / lagindex . size
varfunc = npsum ( ( mags [ lagindex ] - magmed ) * ( mags [ lagindex ] - magmed ) ) / mags . size
acorr = autocovarfunc / varfunc
return acorr


signature:
def _autocorr_func3 ( mags , lag , maglen , magmed , magstd ) :
    '''This is yet another alternative to calculate the autocorrelation .
 Taken from : ` Bayesian Methods for Hackers by Cameron Pilon < http : / / nbviewer . jupyter . org / github / CamDavidsonPilon / Probabilistic - Programming - and - Bayesian - Methods - for - Hackers / blob / master / Chapter3 _ MCMC / Chapter3 . ipynb # Autocorrelation > ` _
 ( This should be the fastest method to calculate ACFs . )
 Parameters
 mags : np . array
 This is the magnitudes array . MUST NOT have any nans .
 lag : float
 The specific lag value to calculate the auto - correlation for . This MUST
 be less than total number of observations in ` mags ` .
 maglen : int
 The number of elements in the ` mags ` array .
 magmed : float
 The median of the ` mags ` array .
 magstd : float
 The standard deviation of the ` mags ` array .
 Returns
 float
 The auto - correlation at this specific ` lag ` value .'''

implementation:
# from http : / / tinyurl . com / afz57c4
result = npcorrelate ( mags , mags , mode = 'full' )
result = result / npmax ( result )
return result [ int ( result . size / 2 ) : ]


signature:
def autocorr_magseries ( times , mags , errs , maxlags = 1000 , func = _autocorr_func3 , fillgaps = 0.0 , filterwindow = 11 , forcetimebin = None , sigclip = 3.0 , magsarefluxes = False , verbose = True ) :
    '''This calculates the ACF of a light curve .
 This will pre - process the light curve to fill in all the gaps and normalize
 everything to zero . If ` fillgaps = ' noiselevel ' ` , fills the gaps with the
 noise level obtained via the procedure above . If ` fillgaps = ' nan ' ` , fills
 the gaps with ` np . nan ` .
 Parameters
 times , mags , errs : np . array
 The measurement time - series and associated errors .
 maxlags : int
 The maximum number of lags to calculate .
 func : Python function
 This is a function to calculate the lags .
 fillgaps : ' noiselevel ' or float
 This sets what to use to fill in gaps in the time series . If this is
 ' noiselevel ' , will smooth the light curve using a point window size of
 ` filterwindow ` ( this should be an odd integer ) , subtract the smoothed LC
 from the actual LC and estimate the RMS . This RMS will be used to fill
 in the gaps . Other useful values here are 0.0 , and npnan .
 filterwindow : int
 The light curve ' s smoothing filter window size to use if
 ` fillgaps = ' noiselevel ` ' .
 forcetimebin : None or float
 This is used to force a particular cadence in the light curve other than
 the automatically determined cadence . This effectively rebins the light
 curve to this cadence . This should be in the same time units as ` times ` .
 sigclip : float or int or sequence of two floats / ints or None
 If a single float or int , a symmetric sigma - clip will be performed using
 the number provided as the sigma - multiplier to cut out from the input
 time - series .
 If a list of two ints / floats is provided , the function will perform an
 ' asymmetric ' sigma - clip . The first element in this list is the sigma
 value to use for fainter flux / mag values ; the second element in this
 list is the sigma value to use for brighter flux / mag values . For
 example , ` sigclip = [ 10 . , 3 . ] ` , will sigclip out greater than 10 - sigma
 dimmings and greater than 3 - sigma brightenings . Here the meaning of
 " dimming " and " brightening " is set by * physics * ( not the magnitude
 system ) , which is why the ` magsarefluxes ` kwarg must be correctly set .
 If ` sigclip ` is None , no sigma - clipping will be performed , and the
 time - series ( with non - finite elems removed ) will be passed through to
 the output .
 magsarefluxes : bool
 If your input measurements in ` mags ` are actually fluxes instead of
 mags , set this is True .
 verbose : bool
 If True , will indicate progress and report errors .
 Returns
 dict
 A dict of the following form is returned : :
 { ' itimes ' : the interpolated time values after gap - filling ,
 ' imags ' : the interpolated mag / flux values after gap - filling ,
 ' ierrs ' : the interpolated mag / flux values after gap - filling ,
 ' cadence ' : the cadence of the output mag / flux time - series ,
 ' minitime ' : the minimum value of the interpolated times array ,
 ' lags ' : the lags used to calculate the auto - correlation function ,
 ' acf ' : the value of the ACF at each lag used }'''

implementation:
# get the gap - filled timeseries
interpolated = fill_magseries_gaps ( times , mags , errs , fillgaps = fillgaps , forcetimebin = forcetimebin , sigclip = sigclip , magsarefluxes = magsarefluxes , filterwindow = filterwindow , verbose = verbose )
if not interpolated :
    print ( 'failed to interpolate light curve to minimum cadence!' )
    return None
itimes , imags = interpolated [ 'itimes' ] , interpolated [ 'imags' ] ,
# calculate the lags up to maxlags
if maxlags :
    lags = nparange ( 0 , maxlags )
else :
    lags = nparange ( itimes . size )
series_stdev = 1.483 * npmedian ( npabs ( imags ) )
if func != _autocorr_func3 : # get the autocorrelation as a function of the lag of the mag series
    autocorr = nparray ( [ func ( imags , x , imags . size , 0.0 , series_stdev ) for x in lags ] )
# this doesn ' t need a lags array
else :
    autocorr = _autocorr_func3 ( imags , lags [ 0 ] , imags . size , 0.0 , series_stdev )
    # return only the maximum number of lags
    if maxlags is not None :
        autocorr = autocorr [ : maxlags ]
interpolated . update ( { 'minitime' : itimes . min ( ) , 'lags' : lags , 'acf' : autocorr } )
return interpolated