<jupyter_start><jupyter_text>Gather all tails<jupyter_code>all_tails = {'B': {0: {}, 1: {}},
             'R': {0: {}, 1: {}},
             'P': {0: {}, 1: {}, 2: {}, 3: {}}}
for fname in os.listdir(DARK_DIR):
    if '_25_' not in fname:
        continue
    frame = int(fname.split('_')[2])
    channel = fname.split('.')[0][-1]
    n_amps = 4 if channel == 'P' else 2
    print(fname, channel, )
    for amp in range(n_amps):
        data = get_data(frame, channel, amp)
        objs = get_objs(data)
        tails, peak_vals = get_tails(data, objs)
        s_tails, _ = get_tails(data, objs, serial=True)
        all_tails[channel][amp][frame] = {'tails_parallel': tails,
                                          'tails_serial': s_tails,
                                          'loc_y': objs.ypeak,
                                          'loc_x': objs.xpeak,
                                          'peak_vals': peak_vals}
for x in all_tails['R'][1].values():
    plt.plot(x['tails_serial'].T)
plt.xlabel('Pixels from hit')
plt.ylabel('Value')
def nmad(arr):
    med = np.median(arr)
    mad = np.median(np.abs(arr-med))
    return 1.4826 * mad
def fit_cte(channel, amp, serial=False, verbose=False):
    peak_vals = np.concatenate([v['peak_vals'] for k, v in all_tails[channel][amp].items()])
    if serial:
        locs = np.concatenate([v['loc_x'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_serial'][:, :4] for k, v in all_tails[channel][amp].items()]),
                             axis=1)
    else:
        locs = np.concatenate([v['loc_y'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_parallel'][:, :4] for k, v in all_tails[channel][amp].items()]),
                             axis=1)

    if serial:
        bins = np.arange(0, 1050, 64)
    else:
        bins = np.arange(0, 5000, 256)
 
    median = np.array([np.median((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))])
                       for i in range(len(bins)-1)])
    nmads = np.array([nmad((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))])
                      for i in range(len(bins)-1)])
    cts = np.array([len((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))])
                    for i in range(len(bins)-1)])
    bins = (bins[1:]+bins[:-1])/2.
    nan_cut = cts==0
    median = median[~nan_cut]
    nmads = nmads[~nan_cut]/cts[~nan_cut]
    bins = bins[~nan_cut]
    
    def chisq(m, b):
        line = m * bins + b
        diff = median - line
        return np.sum(diff**2/nmads**2)
    
    m = Minuit(chisq, pedantic=False, print_level=verbose,
               m=1e-6, b=0,
               limit_m=(-1, 1), limit_b=(-1, 1))
    m.migrad()
    return m.values

def fit_cte_all(channel, amp, serial=False, verbose=False):
    peak_vals = np.concatenate([v['peak_vals'] for k, v in all_tails[channel][amp].items()])
    if serial:
        locs = np.concatenate([v['loc_x'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_serial'][:, :4] for k, v in all_tails[channel][amp].items()]), axis=1)
    else:
        locs = np.concatenate([v['loc_y'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_parallel'][:, :4] for k, v in all_tails[channel][amp].items()]), axis=1)

    if serial:
        bins = np.arange(0, 1050, 64)
    else:
        bins = np.arange(0, 5000, 256)
    
    def chisq(m, b):
        line = m * locs + b
        diff = tail_counts/peak_vals - line
        return np.sum(diff**2)
    
    m = Minuit(chisq, pedantic=False, print_level=verbose,
               m=1e-6, b=0,
               limit_m=(-1, 1), limit_b=(-1, 1))
    m.migrad()
    return m.values
def cte_plot(channel, amp, serial=False, all_vals=False):
    peak_vals = np.concatenate([v['peak_vals'] for k, v in all_tails[channel][amp].items()])
    if serial:
        locs = np.concatenate([v['loc_x'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_serial'] for k, v in all_tails[channel][amp].items()]), axis=1)
    else:
        locs = np.concatenate([v['loc_y'] for k, v in all_tails[channel][amp].items()])
        tail_counts = np.sum(np.concatenate([v['tails_parallel'][:, :4] for k, v in all_tails[channel][amp].items()]), axis=1)
    if all_vals:
        plt.scatter(locs, tail_counts/peak_vals, s=1, alpha=0.2, label='{} channel, amp {}'.format(channel, amp))
    if serial:
        bins = np.arange(0, 1050, 64)
    else:
        bins = np.arange(0, 5000, 256)
    median = [np.median((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))]) for i in range(len(bins)-1)]
    nmads = np.array([nmad((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))]) for i in range(len(bins)-1)])
    cts = np.array([len((tail_counts/peak_vals)[np.where((locs>bins[i]) & (locs<bins[i+1]))]) for i in range(len(bins)-1)])
    nmads /= np.sqrt(cts)
    plt.errorbar((bins[1:]+bins[:-1])/2, median,
                 xerr=(bins[1:]-bins[:-1])/2, yerr=nmads,
                 color='C3', linewidth=0, elinewidth=1, marker='.', label='{} channel, amp {}'.format(channel, amp))
    fit = fit_cte(channel, amp, serial)
    plt.plot(bins, fit['m']*bins + fit['b'], color='k')
    plt.legend(loc='upper right')
#     plt.plot(bins, fit_all['m']*bins + fit_all['b'], color='k', alpha=0.3)
    if serial:
        plt.xlim(0, 1028)
    else:
        plt.xlim(0, 4096)
#     plt.ylim(-0.3, 0.3)
plt.figure(figsize=(12, 6))

plt.subplot(2, 4, 1)
cte_plot('B', 0)
plt.ylabel('Fraction of counts in tail')

plt.subplot(2, 4, 2)
cte_plot('B', 1)

plt.subplot(2, 4, 3)
cte_plot('R', 0)

plt.subplot(2, 4, 4)
cte_plot('R', 1)

plt.subplot(2, 4, 5)
cte_plot('P', 0)
plt.ylabel('Fraction of counts in tail')
plt.xlabel('Y-position')

plt.subplot(2, 4, 6)
cte_plot('P', 1)
plt.xlabel('Y-position')

plt.subplot(2, 4, 7)
cte_plot('P', 2)
plt.xlabel('Y-position')

plt.subplot(2, 4, 8)
cte_plot('P', 3)
plt.xlabel('Y-position')

# plt.savefig('figures/single_night_example_parallel.png', bbox_inches='tight')
plt.savefig('figures/single_night_example_parallel_medians.png', bbox_inches='tight')
plt.figure(figsize=(12, 6))

plt.subplot(2, 4, 1)
cte_plot('B', 0, True)
plt.ylabel('Fraction of counts in tail')

plt.subplot(2, 4, 2)
cte_plot('B', 1, True)

plt.subplot(2, 4, 3)
cte_plot('R', 0, True)

plt.subplot(2, 4, 4)
cte_plot('R', 1, True)

plt.subplot(2, 4, 5)
cte_plot('P', 0, True)
plt.ylabel('Fraction of counts in tail')
plt.xlabel('X-position')

plt.subplot(2, 4, 6)
cte_plot('P', 1, True)
plt.xlabel('X-position')

plt.subplot(2, 4, 7)
cte_plot('P', 2, True)
plt.xlabel('X-position')

plt.subplot(2, 4, 8)
cte_plot('P', 3, True)
plt.xlabel('X-position')
# plt.savefig('figures/single_night_example_serial.png', bbox_inches='tight')
plt.savefig('figures/single_night_example_serial_medians.png', bbox_inches='tight')
fit_cte('R', 1, False)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1929328621908127, 'fraction_numerical': 0.029540636042402826, 'mean_word_length': 3.0854503464203233, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '651194', 'license': 'MIT', 'n_tokens_mistral': 3044, 'n_tokens_neox': 2784, 'n_words': 540, 'path': 'cosmic_rays/do_cte_one_night.ipynb', 'repo_name': 'snfactory/cte', 'size': '7075'}