290425 commit

# This file is machine-generated - editing it directly is not advised

julia_version = "1.11.4"

julia_version = "1.11.5"

manifest_format = "2.0"

project_hash = "4dfdfa82db68047c3e73c2b212e113ca7569973a"

project_hash = "d562d9224bce463ab6e7468b8369592932cdf2c8"

[[deps.AbstractFFTs]]

deps = ["LinearAlgebra"]

uuid = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f"

version = "1.11.0"

[[deps.BitFlags]]

git-tree-sha1 = "0691e34b3bb8be9307330f88d1a3c3f25466c24d"

uuid = "d1d4a3ce-64b1-5f1a-9ba4-7e7e69966f35"

version = "0.1.9"

[[deps.Bonito]]

deps = ["Base64", "CodecZlib", "Colors", "Dates", "Deno_jll", "HTTP", "Hyperscript", "LinearAlgebra", "Markdown", "MsgPack", "Observables", "RelocatableFolders", "SHA", "Sockets", "Tables", "ThreadPools", "URIs", "UUIDs", "WidgetsBase"]

git-tree-sha1 = "91302336aa0c70d6eed21d6fbae6823168903b89"

uuid = "824d6782-a2ef-11e9-3a09-e5662e0c26f8"

version = "3.2.2"

[[deps.Bzip2_jll]]

deps = ["Artifacts", "JLLWrappers", "Libdl"]

git-tree-sha1 = "1b96ea4a01afe0ea4090c5c8039690672dd13f2e"

    [deps.ChainRulesCore.extensions]

    ChainRulesCoreSparseArraysExt = "SparseArrays"

[[deps.CodecZlib]]

deps = ["TranscodingStreams", "Zlib_jll"]

git-tree-sha1 = "962834c22b66e32aa10f7611c08c8ca4e20749a9"

uuid = "944b1d66-785c-5afd-91f1-9de20f533193"

version = "0.7.8"

[[deps.ColorBrewer]]

deps = ["Colors", "JSON"]

git-tree-sha1 = "e771a63cc8b539eca78c85b0cabd9233d6c8f06f"

uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae"

version = "1.1.1+0"

[[deps.ConcurrentUtilities]]

deps = ["Serialization", "Sockets"]

git-tree-sha1 = "d9d26935a0bcffc87d2613ce14c527c99fc543fd"

uuid = "f0e56b4a-5159-44fe-b623-3e5288b988bb"

version = "2.5.0"

[[deps.ConstructionBase]]

git-tree-sha1 = "76219f1ed5771adbb096743bff43fb5fdd4c1157"

uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9"

uuid = "927a84f5-c5f4-47a5-9785-b46e178433df"

version = "1.6.4"

[[deps.Deno_jll]]

deps = ["Artifacts", "JLLWrappers", "Libdl"]

git-tree-sha1 = "82adce9b61f01f036934b6fb2b100636df941b37"

uuid = "04572ae6-984a-583e-9378-9577a1c2574d"

version = "2.2.6+0"

[[deps.Distributed]]

deps = ["Random", "Serialization", "Sockets"]

uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b"

uuid = "429591f6-91af-11e9-00e2-59fbe8cec110"

version = "2.2.8"

[[deps.ExceptionUnwrapping]]

deps = ["Test"]

git-tree-sha1 = "d36f682e590a83d63d1c7dbd287573764682d12a"

uuid = "460bff9d-24e4-43bc-9d9f-a8973cb893f4"

version = "0.1.11"

[[deps.Expat_jll]]

deps = ["Artifacts", "JLLWrappers", "Libdl"]

git-tree-sha1 = "d55dffd9ae73ff72f1c0482454dcf2ec6c6c4a63"

git-tree-sha1 = "b66970a70db13f45b7e57fbda1736e1cf72174ea"

uuid = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549"

version = "1.17.0"

weakdeps = ["HTTP"]

    [deps.FileIO.extensions]

    HTTPExt = "HTTP"

    [deps.FileIO.weakdeps]

    HTTP = "cd3eb016-35fb-5094-929b-558a96fad6f3"

[[deps.FilePaths]]

deps = ["FilePathsBase", "MacroTools", "Reexport", "Requires"]

git-tree-sha1 = "919d9412dbf53a2e6fe74af62a73ceed0bce0629"

uuid = "42e2da0e-8278-4e71-bc24-59509adca0fe"

version = "1.0.2"

[[deps.HTTP]]

deps = ["Base64", "CodecZlib", "ConcurrentUtilities", "Dates", "ExceptionUnwrapping", "Logging", "LoggingExtras", "MbedTLS", "NetworkOptions", "OpenSSL", "PrecompileTools", "Random", "SimpleBufferStream", "Sockets", "URIs", "UUIDs"]

git-tree-sha1 = "f93655dc73d7a0b4a368e3c0bce296ae035ad76e"

uuid = "cd3eb016-35fb-5094-929b-558a96fad6f3"

version = "1.10.16"

[[deps.HarfBuzz_jll]]

deps = ["Artifacts", "Cairo_jll", "Fontconfig_jll", "FreeType2_jll", "Glib_jll", "Graphite2_jll", "JLLWrappers", "Libdl", "Libffi_jll"]

git-tree-sha1 = "55c53be97790242c29031e5cd45e8ac296dadda3"

git-tree-sha1 = "dba9ddf07f77f60450fe5d2e2beb9854d9a49bd0"

uuid = "8197267c-284f-5f27-9208-e0e47529a953"

version = "0.7.10"

weakdeps = ["Random", "RecipesBase", "Statistics"]

    [deps.IntervalSets.extensions]

    IntervalSetsRandomExt = "Random"

    IntervalSetsRecipesBaseExt = "RecipesBase"

    IntervalSetsStatisticsExt = "Statistics"

    [deps.IntervalSets.weakdeps]

    Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

    RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"

    Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[[deps.IrrationalConstants]]

git-tree-sha1 = "e2222959fbc6c19554dc15174c81bf7bf3aa691c"

uuid = "92d709cd-6900-40b7-9082-c6be49f344b6"

uuid = "56ddb016-857b-54e1-b83d-db4d58db5568"

version = "1.11.0"

[[deps.LoggingExtras]]

deps = ["Dates", "Logging"]

git-tree-sha1 = "f02b56007b064fbfddb4c9cd60161b6dd0f40df3"

uuid = "e6f89c97-d47a-5376-807f-9c37f3926c36"

version = "1.1.0"

[[deps.MIMEs]]

git-tree-sha1 = "1833212fd6f580c20d4291da9c1b4e8a655b128e"

uuid = "6c6e2e6c-3030-632d-7369-2d6c69616d65"

uuid = "0a4f8689-d25c-4efe-a92b-7142dfc1aa53"

version = "0.6.2"

[[deps.MbedTLS]]

deps = ["Dates", "MbedTLS_jll", "MozillaCACerts_jll", "NetworkOptions", "Random", "Sockets"]

git-tree-sha1 = "c067a280ddc25f196b5e7df3877c6b226d390aaf"

uuid = "739be429-bea8-5141-9913-cc70e7f3736d"

version = "1.1.9"

[[deps.MbedTLS_jll]]

deps = ["Artifacts", "Libdl"]

uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1"

uuid = "14a3606d-f60d-562e-9121-12d972cd8159"

version = "2023.12.12"

[[deps.MsgPack]]

deps = ["Serialization"]

git-tree-sha1 = "f5db02ae992c260e4826fe78c942954b48e1d9c2"

uuid = "99f44e22-a591-53d1-9472-aa23ef4bd671"

version = "1.2.1"

[[deps.NaNMath]]

deps = ["OpenLibm_jll"]

git-tree-sha1 = "cc0a5deefdb12ab3a096f00a6d42133af4560d71"

[[deps.OpenLibm_jll]]

deps = ["Artifacts", "Libdl"]

uuid = "05823500-19ac-5b8b-9628-191a04bc5112"

version = "0.8.1+4"

version = "0.8.5+0"

[[deps.OpenSSL]]

deps = ["BitFlags", "Dates", "MozillaCACerts_jll", "OpenSSL_jll", "Sockets"]

git-tree-sha1 = "38cb508d080d21dc1128f7fb04f20387ed4c0af4"

uuid = "4d8831e6-92b7-49fb-bdf8-b643e874388c"

version = "1.4.3"

[[deps.OpenSSL_jll]]

deps = ["Artifacts", "JLLWrappers", "Libdl"]

    [deps.Ratios.extensions]

    RatiosFixedPointNumbersExt = "FixedPointNumbers"

[[deps.RecipesBase]]

deps = ["PrecompileTools"]

git-tree-sha1 = "5c3d09cc4f31f5fc6af001c250bf1278733100ff"

uuid = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"

version = "1.3.4"

[[deps.Reexport]]

git-tree-sha1 = "45e428421666073eab6f2da5c9d310d99bb12f9b"

uuid = "189a3867-3050-52da-a836-e630ba90ab69"

uuid = "73760f76-fbc4-59ce-8f25-708e95d2df96"

version = "0.4.0"

[[deps.SimpleBufferStream]]

git-tree-sha1 = "f305871d2f381d21527c770d4788c06c097c9bc1"

uuid = "777ac1f9-54b0-4bf8-805c-2214025038e7"

version = "1.2.0"

[[deps.SimpleTraits]]

deps = ["InteractiveUtils", "MacroTools"]

git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231"

uuid = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

version = "1.11.0"

[[deps.ThreadPools]]

deps = ["Printf", "RecipesBase", "Statistics"]

git-tree-sha1 = "50cb5f85d5646bc1422aa0238aa5bfca99ca9ae7"

uuid = "b189fb0b-2eb5-4ed4-bc0c-d34c51242431"

version = "2.1.1"

[[deps.TiffImages]]

deps = ["ColorTypes", "DataStructures", "DocStringExtensions", "FileIO", "FixedPointNumbers", "IndirectArrays", "Inflate", "Mmap", "OffsetArrays", "PkgVersion", "ProgressMeter", "SIMD", "UUIDs"]

git-tree-sha1 = "f21231b166166bebc73b99cea236071eb047525b"

    ConstructionBase = "187b0558-2788-49d3-abe0-74a17ed4e7c9"

    InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112"

[[deps.WGLMakie]]

deps = ["Bonito", "Colors", "FileIO", "FreeTypeAbstraction", "GeometryBasics", "Hyperscript", "LinearAlgebra", "Makie", "Observables", "PNGFiles", "PrecompileTools", "RelocatableFolders", "ShaderAbstractions", "StaticArrays"]

git-tree-sha1 = "dca58a991d7ea60ec5adf8c57afee5bf94c2539d"

uuid = "276b4fcb-3e11-5398-bf8b-a0c2d153d008"

version = "0.10.5"

[[deps.WebP]]

deps = ["CEnum", "ColorTypes", "FileIO", "FixedPointNumbers", "ImageCore", "libwebp_jll"]

git-tree-sha1 = "aa1ca3c47f119fbdae8770c29820e5e6119b83f2"

uuid = "e3aaa7dc-3e4b-44e0-be63-ffb868ccd7c1"

version = "0.1.3"

[[deps.WidgetsBase]]

deps = ["Observables"]

git-tree-sha1 = "30a1d631eb06e8c868c559599f915a62d55c2601"

uuid = "eead4739-05f7-45a1-878c-cee36b57321c"

version = "0.1.4"

[[deps.WoodburyMatrices]]

deps = ["LinearAlgebra", "SparseArrays"]

git-tree-sha1 = "c1a7aa6219628fcd757dede0ca95e245c5cd9511"

CommonMark = "a80b9123-70ca-4bc0-993e-6e3bcb318db6"

Format = "1fa38f19-a742-5d3f-a2b9-30dd87b9d5f8"

PlutoUI = "7f904dfe-b85e-4ff6-b463-dae2292396a8"

WGLMakie = "276b4fcb-3e11-5398-bf8b-a0c2d153d008"

# -*- coding: utf-8 -*-

"""

Created on Sun Sep 20 18:13:40 2015

@author: Felix Darvas

compute the Gregory-Laredo algorithm on arival times

This function computes the likelihood of a set of arrival times originating 

from a periodic system rather than constant rate (poisson) background noise

based on

Gregory, P. C. and Thomas. J. Loredo, 1992, 

"A New Method For The Detection Of A Periodic Signal Of Unknown Shape And Period" 

in 

The Astrophysical Journal, Astrophysical J., 398, p.146

inputs:

Tlist    -  list of arrival times, numpy int array

m_max    -  max number of bins typically 12-15, we use 12 as default

w_range  -  frequency range to scan numpy float array of frequency values 

           default is  w_lo=20*pi/T at delta w = pi/T to w_hi=pi*N/T

           where N=#arrival times, T=observation time

ni       - number of integration steps, default ni=10

parallel - use paralel execution - default is off

outut:

O_period - Odds ratio for a periodic process vs. constant rate process

p_period - probability of a periodic process 0<=p_period<=1

m_opt    - optimal bin size 1<= m_opt <=m_max

S        - The probability spectrum 

w        - The frequency range for S

w_peak   - the expected frequency of the process

w_conf   - 95% confidence interval of w_peak

"""

import numpy as np 

import multiprocessing as mp 

import functools

def compute_bin(Tlist,m,w,p):

    n=np.zeros(m,'int')

    j=np.floor(m*np.mod(w*Tlist+p,2*np.pi)/(2*np.pi))

    j.astype(int)

    for u in range(0,m):

        n[u]=np.size(np.extract(j==u,j))

    return n

def compute_factor(N,m,v): 

    # compute m^N *(N+m-1) over N /(2pi v)

    # whic is used to scale the multiplicity function (see eqn.5.25 of the paper )

    # we return the log of the integral scale here to avoid numerical overflow

    f1=N*np.log(m) # log of m^N

    f2=np.sum(np.log(np.arange(1,N+m))) # log of (N+m-1)!

    f3=np.sum(np.log(np.arange(1,m+1))) # log of m!

    f=f1+f3-f2-np.log(2*np.pi*v)

    return f

def precompute_binmult(N):

	# precompute all potential bin factorials

    fbin=np.zeros(int(N)+1)

    for i in range(2,int(N+1)): # n=0 -> define log(n)=0, n=1, log(n)=0, so no need to compute n=0,1

        fbin[i]=fbin[i-1]+np.log(i)

    return fbin

def compute_W_scaled(Tlist,m,w,phase,factor,fbin):

    # compute the scaled multiplicity 1/Wm(w,phase) (see eqn. 5.25)

    # note that for large arrival time numbers the multiplicity becomes

    # excessively large. Since W_m(phase,w) is never needed explicitly,

    # we use the scaled version. 

    # input factor is the log of the actual factor

    n=compute_bin(Tlist,m,w,phase) # find bin histogram for a given bin number m, frequency w and phase 

    f=0    

    for i in range(0,m):

        #f=f+np.sum(np.log(np.arange(2,n[i]+1)))

        f=f+fbin[n[i]]

    y=np.exp(f+factor)

    return y

def compute_Om1(w,Tlist,m,factor,fbin,ni):

    #compute  specific odds-ratios (eqn. 5.25)

    p=np.arange(0,ni)/float(ni)*2*np.pi/m # intgration range, only integrate over a single bin, as values of the integral repeat over bins       

    y=np.zeros(np.size(p),'float') # array to hold values of W_scaled over the integration range

    for i in range(0,np.size(y)):

        y[i]=compute_W_scaled(Tlist,m,w,p[i],factor,fbin)

    return np.trapz(y,p)*m # return intregrated W_Scaled

def compute_Om1wPar(Tlist,m_max,w,fa,fbin,ni): # compute odds-ratios for bins and frequencies

	# parallel version 

    Om1w=np.zeros((m_max,np.size(w)),'float') # odds ratio matrix

    pool=mp.Pool() # use all workers

    for m in range(0,m_max):

        Om1w[m,:]=pool.map(functools.partial(compute_Om1,Tlist=Tlist,m=(m+1),factor=fa[m],fbin=fbin,ni=ni),w)

    return Om1w

def compute_Om1w(Tlist,m_max,w,fa,fbin,ni): # compute odds-ratios for bins and frequencies

    Om1w=np.zeros((m_max,np.size(w)),'float') # odds ratio matrix

    for m in range(0,m_max):

        for wi in range(0,np.size(w)):

            Om1w[m,wi]=compute_Om1(w[wi],Tlist,m+1,fa[m],fbin,ni)

    return Om1w

def compute_GL(Tlist,m_max=12,w_range=None,ni=10,parallel=False):

    # initialize output values

    O_period=None

    p_period=None

    m_opt=None

    S=None

    w=None

    w_peak=None

    w_mean=None

    w_conf=None

    N=float(np.size(Tlist)) # need float value to avoid int/int

    if N>0:

        # compute GL algorithm

        fbin=precompute_binmult(N)

        v=m_max-1

        T=float(np.max(Tlist)) # duration of the observation

        if w_range is None: # use default frequencies

            w_hi=np.pi*N/T # max default frequency

            w_lo=np.minimum(20,N/10)*np.pi/T # min default frequency

            dw=np.pi/T/10 # step size

            w=np.arange(w_lo,w_hi,dw)

            if np.size(w)<2:

                print ("error ")

                raise ValueError('bad arrival time list')

        else: # use user supplied frequency vector

            w=w_range 

            w_hi=np.max(w_range)

            w_lo=np.min(w_range)

        if w_lo==0:

            # cannot have w_lo =0

            print ("minimum frequency cannot be 0!\n")

            return

        fa=np.zeros(m_max)

        for m in range(0,m_max): # precompute factors for each m

            fa[m]=compute_factor(N,m+1,v)

        if parallel:

            Om1w=compute_Om1wPar(Tlist,m_max,w,fa,fbin,ni)

        else:

            Om1w=compute_Om1w(Tlist,m_max,w,fa,fbin,ni)

        pw=1/w/np.log(w_hi/w_lo)

        O1m=np.zeros(m_max)

        for i in range(0,m_max): # intgreate odd ratios across frequencies

            O1m[i]=np.trapz(pw*Om1w[i],w)

        m_opt=np.argmax(O1m) # find optimum bin number, i.e largest odds-ratio    

        S=Om1w[m_opt]/w # compute Spectral probability

        m_opt=m_opt+1 # start bin index with 1

        C=np.trapz(S,w) # compute normalization

        S=S/C # normalized probability

        O_period=np.sum(O1m) # integrated odds ratio

        p_period=O_period/(1+O_period) # likelihood of periodic event

        cdf=np.array(S)

        for i in range(0,np.size(S)):

            cdf[i]=np.trapz(S[0:i],w[0:i])

        wr=np.extract(np.logical_and(cdf>.025, cdf<.975),w)

        w_peak=w[np.argmax(S)]

        w_mean=np.trapz(S*w,w)

        if np.size(wr)>0:

            w_conf=[np.min(wr),np.max(wr)]

        else:

            w_conf=[w_peak,w_peak]

        return O_period,p_period,m_opt,S,w,w_peak,w_mean,w_conf

    else:

        # throw an error

        print ("No valid arrival time array provided!\n")

        return O_period,p_period,m_opt,S,w,w_peak,w_mean,w_conf