~jb753/clusterfunc

329a1e819f1eff6740208af1bc0e8844bc649d34 — James Brind 1 year, 1 month ago efd00d7
Free N by bisection
1 files changed, 37 insertions(+), 9 deletions(-)

M clusterfunc.py
M clusterfunc.py => clusterfunc.py +37 -9
@@ 63,8 63,6 @@ def normalised_N_fixed(Dst, Dmax, ERmax, N, check=True):
    # Number of points for uniform spacing at Dst
    Nunif = np.ceil(1./Dst)

    print(N1a, N1b, N, L1a)

    if N > Nunif:
        raise ClusteringException(
        f'Too many points N={N} for a uniform spacing of Dst={Dst}.'


@@ 156,6 154,34 @@ def normalised_N_fixed(Dst, Dmax, ERmax, N, check=True):

    return x

def normalised_N_free(Dst, Dmax, ERmax, mult=8, check=True):

    # Find high and low guesses for n, N = mult*n + 1
    nlow = np.round(1./Dmax/mult).astype(int)
    nhigh = np.round(1./Dst/mult).astype(int)

    # Loop until converged
    while nhigh-nlow > 1:

        # Bisect high and low guesses
        nnew = (nhigh+nlow)//2
        Nnew = nnew * mult + 1

        # Attempt to cluster
        try:
            x = normalised_N_fixed(Dst, Dmax, ERmax, Nnew, check=False)
            # new n is valid, use as high limiit
            nhigh = nnew
        except ClusteringException:
            # new n is not valid, use as low limiit
            nlow = nnew

    if check:
        check_normalised(x, Dst, Dmax, ERmax)

    return x



def cumsum0(x, axis=None):
    return np.insert(np.cumsum(x, axis=axis), 0, 0.0, axis=axis)


@@ 169,19 195,21 @@ import matplotlib.pyplot as plt
Dst = 0.001
Dmax = 0.08
ER = 1.1
N = 56
N = 57

xx = normalised_N_fixed(Dst, Dmax, ER, N)
x1 = normalised_N_fixed(Dst, Dmax, ER, N)
x2 = normalised_N_free(Dst, Dmax, ER)

fig, ax = plt.subplots()
ax.plot(xx,'k-x')
ax.plot(x1,'-x')
ax.plot(x2,'-x')
ax.set_title('x')
plt.show()

fig, ax = plt.subplots()
ax.set_title('dx')
ax.plot(np.diff(xx),'k-x')
plt.show()
# fig, ax = plt.subplots()
# ax.set_title('dx')
# ax.plot(np.diff(xx),'k-x')
# plt.show()

# fig, ax = plt.subplots()
# ax.plot(dx,'k-x')