<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Goal
Step1:  Load some data
I'm going to work with the data from the combined data sets. The analysis for this data set is in analysis\Cf072115_to_Cf072215b.
The one limitation here is that this data has already cut out the fission chamber neighbors. 
det_df without fission chamber neighbors
Step2:  Specify energy range
Step3:  singles_hist_e_n.npz
Step4:  Load bhp_nn_e for all pairs
I'm going to skip a few steps in order to save memory. This data was produced in analysis_build_bhp_nn_by_pair_1_ns.ipynb and is stored in datap\bhp_nn_by_pair_1ns.npz. Load it now, as explained in the notebook.
Step5:  Set up det_df columns and singles_df
Step6:  Calculate and fill doubles sums
Step7:  Calculate singles sums
Step8:  Calculate W values
Step9:  Condense to angle bin
Step10:  Plot it
Step11:  Save to disk
In order to compare datasets, it would be nice to save these results to disk and reload in another notebook for comparison. These results are pretty easy, format-wise, so I'll just use the built-in pandas methods.
Step12:  Reload
<ASSISTANT_TASK:>
Python Code:
import os
import sys
import matplotlib.pyplot as plt
import numpy as np
import imageio
import pandas as pd
import seaborn as sns
sns.set(style='ticks')
sys.path.append('../scripts/')
import bicorr as bicorr
import bicorr_e as bicorr_e
import bicorr_plot as bicorr_plot
import bicorr_sums as bicorr_sums
import bicorr_math as bicorr_math
%load_ext autoreload
%autoreload 2
Explanation: Goal: Correct for singles rate with $W$ calculation
In order to correct for differences in detection efficiencies and solid angles, we will divide all of the doubles rates by the singles rates of the two detectors as follows:
$ W_{i,j} = \frac{D_{i,j}}{S_i*S_j}$
This requires calculating $S_i$ and $S_j$ from the cced files. I need to rewrite my analysis from the beginning, or write another function that parses the cced file.
In this file, I will import the singles and bicorr data and calculate all  $D_{i,j}$, $S_i$, $S_j$, and $W_{i,j}$. 
This notebook does the analysis in energy space.
End of explanation
det_df = bicorr.load_det_df('../meas_info/det_df_pairs_angles.csv')
det_df.head()
chList, fcList, detList, num_dets, num_det_pairs = bicorr.build_ch_lists()
dict_pair_to_index, dict_index_to_pair, dict_pair_to_angle = bicorr.build_dict_det_pair(det_df)
num_fissions = 2194651200.00
Explanation: Load some data
I'm going to work with the data from the combined data sets. The analysis for this data set is in analysis\Cf072115_to_Cf072215b.
The one limitation here is that this data has already cut out the fission chamber neighbors. 
det_df without fission chamber neighbors
End of explanation
e_min = 0.62
e_max = 12
Explanation: Specify energy range
End of explanation
singles_hist_e_n, e_bin_edges, dict_det_to_index, dict_index_to_det = bicorr_e.load_singles_hist_both(filepath = '../analysis/Cf072115_to_Cf072215b/datap/',plot_flag=True, show_flag=True)
bicorr_plot.plot_singles_hist_e_n(singles_hist_e_n, e_bin_edges, show_flag=False, clear_flag=False)
for e in [e_min, e_max]:
    plt.axvline(e,c='r') 
plt.show()
singles_hist_e_n.shape
Explanation: singles_hist_e_n.npz
End of explanation
bhm_e, e_bin_edges, note = bicorr_e.load_bhm_e('../analysis/Cf072115_to_Cf072215b/datap')
bhm_e.shape
bhp_e = np.zeros((len(det_df),len(e_bin_edges)-1,len(e_bin_edges)-1))
bhp_e.shape
for index in det_df.index.values: # index is same as in `bhm`
    bhp_e[index,:,:] = bicorr_e.build_bhp_e(bhm_e,e_bin_edges,pair_is=[index])[0]
bicorr_plot.bhp_e_plot(np.sum(bhp_e,axis=0),e_bin_edges, show_flag=True)
Explanation: Load bhp_nn_e for all pairs
I'm going to skip a few steps in order to save memory. This data was produced in analysis_build_bhp_nn_by_pair_1_ns.ipynb and is stored in datap\bhp_nn_by_pair_1ns.npz. Load it now, as explained in the notebook.
End of explanation
det_df.head()
det_df = bicorr_sums.init_det_df_sums(det_df)
det_df.head()
singles_e_df = bicorr_sums.init_singles_e_df(dict_index_to_det)
singles_e_df.head()
Explanation: Set up det_df columns and singles_df
End of explanation
bhp_e.shape
det_df, energies_real = bicorr_sums.fill_det_df_doubles_e_sums(det_df, bhp_e, e_bin_edges, e_min, e_max, True)
det_df.head()
bicorr_plot.counts_vs_angle_all(det_df, save_flag=False)
Explanation: Calculate and fill doubles sums
End of explanation
singles_e_df.head()
bicorr_plot.Sd_vs_ch_all(singles_e_df, save_flag=False)
det_df = bicorr_sums.fill_det_df_singles_sums(det_df, singles_e_df)
det_df.head()
plt.figure(figsize=(4,4))
ax = plt.gca()
sc = ax.scatter(det_df['d1'],det_df['d2'],s=13,marker='s',
                edgecolor = 'none', c=det_df['Cd'],cmap='viridis')
ax.set_xlabel('d1 channel')
ax.set_ylabel('d2 channel')
ax.set_title('Doubles counts')
cbar = plt.colorbar(sc,fraction=0.043,pad=0.1)
cbar.set_label('Doubles counts')
plt.show()
plt.figure(figsize=(4,4))
ax = plt.gca()
sc = ax.scatter(det_df['d1'],det_df['d2'],s=13,marker='s',
                edgecolor = 'none', c=det_df['Sd1'],cmap='viridis')
ax.set_xlabel('d1 channel')
ax.set_ylabel('d2 channel')
ax.set_title('D1 singles counts')
cbar = plt.colorbar(sc,fraction=0.043,pad=0.1)
cbar.set_label('D1 singles counts')
plt.show()
plt.figure(figsize=(4,4))
ax = plt.gca()
sc = ax.scatter(det_df['d1'],det_df['d2'],s=13,marker='s',
                edgecolor = 'none', c=det_df['Sd2'],cmap='viridis')
ax.set_xlabel('d1 channel')
ax.set_ylabel('d2 channel')
ax.set_title('D2 singles counts')
cbar = plt.colorbar(sc,fraction=0.043,pad=0.1)
cbar.set_label('D2 Singles counts')
plt.show()
Explanation: Calculate singles sums
End of explanation
det_df = bicorr_sums.calc_det_df_W(det_df)
det_df.head()
plt.figure(figsize=(4,4))
ax = plt.gca()
sc = ax.scatter(det_df['d1'],det_df['d2'],s=13,marker='s',
                edgecolor = 'none', c=det_df['W'],cmap='viridis')
ax.set_xlabel('d1 channel')
ax.set_ylabel('d2 channel')
ax.set_title('W')
cbar = plt.colorbar(sc,fraction=0.043,pad=0.1)
cbar.set_label('W')
plt.show()
chIgnore = [1,17,33]
det_df_ignore = det_df[~det_df['d1'].isin(chIgnore) & ~det_df['d2'].isin(chIgnore)]
bicorr_plot.W_vs_angle_all(det_df_ignore, save_flag=False)
bicorr_plot.W_vs_angle_all?
Explanation: Calculate W values
End of explanation
angle_bin_edges = np.arange(8,190,10)
print(angle_bin_edges)
by_angle_df = bicorr_sums.condense_det_df_by_angle(det_df_ignore, angle_bin_edges)
by_angle_df.head()
Explanation: Condense to angle bin
End of explanation
bicorr_plot.W_vs_angle(det_df_ignore, by_angle_df, save_flag=False)
Explanation: Plot it
End of explanation
singles_e_df.to_csv('singles_e_df_filled.csv')
det_df.to_csv(r'det_df_e_filled.csv')
by_angle_df.to_csv(r'by_angle_e_df.csv')
Explanation: Save to disk
In order to compare datasets, it would be nice to save these results to disk and reload in another notebook for comparison. These results are pretty easy, format-wise, so I'll just use the built-in pandas methods.
End of explanation
det_df_filled = pd.read_csv(r'det_df_e_filled.csv',index_col=0)
det_df_filled.head()
chIgnore = [1,17,33]
det_df_ignore = det_df_filled[~det_df_filled['d1'].isin(chIgnore) & ~det_df_filled['d2'].isin(chIgnore)]
det_df_ignore.head()
singles_e_df_filled = pd.read_csv(r'singles_e_df_filled.csv',index_col=0)
singles_e_df_filled.head()
by_angle_e_df = pd.read_csv(r'by_angle_e_df.csv',index_col=0)
by_angle_e_df.head()
bicorr_plot.W_vs_angle(det_df_ignore, by_angle_e_df, save_flag=False)
Explanation: Reload
End of explanation
<END_TASK>