Browse Source

Reject compton measurements on the border of the energy points

Убрал из рассмотрения комптоновские измерения, которые пришлись на границу между двумя
энергетическими точками.
compton 3 years ago
parent
commit
aec4ea3225
2 changed files with 37 additions and 14 deletions
  1. 31 9
      compton_combiner.py
  2. 6 5
      compton_filter.py

+ 31 - 9
compton_combiner.py

@@ -167,7 +167,7 @@ class Combiner():
         df['run_duration'] = df['run_stop'] - df['run_start']
         df['run_in_measurement'] = df['common_duration']/df['run_duration']
         df = df.sort_values(by='run_in_measurement', ascending=False).drop_duplicates(subset='run').sort_values(by='run')
-        df = df.drop(['run_duration', 'common_duration', 'run_start', 'run_stop'], axis=1)
+        df = df.drop(['run_duration', 'common_duration'], axis=1) #, 'run_start', 'run_stop'
         return df
         
     def __del__(self):
@@ -258,8 +258,8 @@ def calculate_point(comb_df: pd.DataFrame, runs_df: pd.DataFrame, compton_df: pd
     
     Returns
     -------
-    dict
-        average parameters on this DataFrame
+    dict, pd.DataFrame
+        average parameters on this DataFrame, clean dataFrame
     """
         
     if (len(comb_df) == 1) and pd.isnull(comb_df.iloc[0].at['compton_start']):
@@ -305,7 +305,7 @@ def calculate_point(comb_df: pd.DataFrame, runs_df: pd.DataFrame, compton_df: pd
         }, res_df
         
     
-    df = comb_df.copy()
+    df = comb_df.loc[~comb_df.compton_start.isna()].copy()
     df.spread_std = np.where(df.spread_std < 1e-4, 1e-4, df.spread_std)
     
     df = df[df.e_std > 0]
@@ -313,7 +313,7 @@ def calculate_point(comb_df: pd.DataFrame, runs_df: pd.DataFrame, compton_df: pd
     # std_energy = np.sqrt(1/np.sum((df.luminosity/df.luminosity.mean())/df.e_std**2))
     
     good_criterion = np.abs((df.e_mean - mean_energy)/np.sqrt(df.e_mean.std()**2 + df.e_std**2)) < 5
-    good_criterion = good_criterion
+    #print(df[~good_criterion].elabel.value_counts())
     df = df[good_criterion]
     
     m = Minuit(Likelihood(df.e_mean, df.e_std, df.luminosity), mean=df.e_mean.mean(), sigma=df.e_mean.std())
@@ -339,6 +339,25 @@ def calculate_point(comb_df: pd.DataFrame, runs_df: pd.DataFrame, compton_df: pd
         'comment': '',
     }
     return res_dict, df
+    
+def process_intersected_compton_meas(combined_df: pd.DataFrame) -> pd.DataFrame:
+    """Replaces compton measurements writed on the border of two energy points on NaNs
+    """
+    
+    energy_point_borders = combined_df[['point_idx', 'elabel', 'run_start', 'run_stop']].groupby(['point_idx'], dropna=True).agg(
+        elabel_start_time=('run_start', 'min'), elabel_stop_time=('run_stop', 'max'),
+    )
+    df_comb = combined_df.set_index('point_idx').join(energy_point_borders, how='left')
+    
+    df_comb['comptonmeas_in_elabel'] = (df_comb[['elabel_stop_time', 'compton_stop']].min(axis=1) - df_comb[['elabel_start_time', 'compton_start']].max(axis=1))/(df_comb['compton_stop'] - df_comb['compton_start'])
+    #print(df_comb['comptonmeas_in_elabel'].dropna().sort_values())
+    df_comb = df_comb.query('comptonmeas_in_elabel < 0.7')
+    border_comptons = df_comb.compton_start.values
+    #print(combined_df.compton_start.isin(border_comptons).sum())
+    #print(combined_df.loc[combined_df.compton_start.isin(border_comptons)].elabel.value_counts())
+    
+    combined_df.loc[combined_df.compton_start.isin(border_comptons), ['compton_start', 'compton_stop', 'e_mean', 'e_std', 'spread_mean', 'spread_std']] = np.nan
+    return combined_df
 
 def process_combined(combined_df: pd.DataFrame, runs_df: pd.DataFrame, compton_df: pd.DataFrame, pics_folder: Optional[str] = None, rdb: Optional[RunsDBHandler] = None) -> pd.DataFrame:
     
@@ -358,6 +377,8 @@ def process_combined(combined_df: pd.DataFrame, runs_df: pd.DataFrame, compton_d
     combined_df['luminosity'] = combined_df['luminosity'].fillna(0)
     
     combined_df['point_idx'] = np.cumsum(~np.isclose(combined_df.elabel, combined_df.elabel.shift(1), atol=1e-4))
+    combined_df = process_intersected_compton_meas(combined_df)
+    
     combined_df = combined_df.groupby(['point_idx', 'compton_start'], dropna=False).agg(
         elabel=('elabel', 'min'), elabel_test=('elabel', 'max'),
         run_first=('run', 'min'), run_last=('run', 'max'),
@@ -368,7 +389,7 @@ def process_combined(combined_df: pd.DataFrame, runs_df: pd.DataFrame, compton_d
         spread_mean=('spread_mean', 'min'), spread_mean_test=('spread_mean', 'max'),
         spread_std=('spread_std', 'min'), spread_std_test=('spread_std', 'max'),
     ).reset_index().set_index('point_idx')
-    # return combined_df
+    #return combined_df
     
     result_df = pd.DataFrame(columns=['energy_point', 'first_run', 'last_run', 'mean_energy', 'mean_energy_stat_err', 'mean_energy_sys_err', 'mean_spread', 'mean_spread_stat_err', 'used_lum', 'comment'])
     
@@ -385,7 +406,7 @@ def process_combined(combined_df: pd.DataFrame, runs_df: pd.DataFrame, compton_d
             half_timedelta = (plt_table.compton_stop - plt_table.compton_start)/2 
             time = plt_table.compton_start + half_timedelta
             dlt0, total_time = timedelta(days=1), plt_table.compton_stop.max() - plt_table.compton_stop.min()
-            timelim = [plt_table.compton_stop.min() - 0.05*total_time, plt_table.compton_stop.max() + 0.05*total_time]
+            timelim = [plt_table.compton_start.min() - 0.05*total_time, plt_table.compton_stop.max() + 0.05*total_time]
             
             fig, ax = plt.subplots(1, 1, dpi=120, tight_layout=True)
             ax.errorbar(time, plt_table.e_mean, xerr=half_timedelta, yerr=plt_table.e_std, fmt='.')
@@ -396,7 +417,8 @@ def process_combined(combined_df: pd.DataFrame, runs_df: pd.DataFrame, compton_d
             ax.tick_params(axis='x', labelrotation=45)
             ax.xaxis.set_major_locator(locator)
             ax.xaxis.set_major_formatter(formatter)
-            ax.set(title=f'{res_dict["energy_point"]}, E = {res_dict["mean_energy"]:.3f} ± {res_dict["mean_energy_stat_err"]:.3f} ± {res_dict["mean_energy_sys_err"]:.3f} MeV', xlabel='Time, NSK', ylabel='Energy, [MeV]', xlim=timelim)
+            ax.set(title=f'{res_dict["energy_point"]}, E = {res_dict["mean_energy"]:.3f} ± {res_dict["mean_energy_stat_err"]:.3f} ± {res_dict["mean_energy_sys_err"]:.3f} MeV', 
+                   xlabel='Time, NSK', ylabel='Energy, [MeV]', xlim=timelim)
             plt.savefig(f'{pics_folder}/{res_dict["first_run"]}_{res_dict["energy_point"]}.png', transparent=True)
             plt.close()
     
@@ -484,4 +506,4 @@ def main():
     return 
 
 if __name__ == "__main__":
-    main()
+    main()

+ 6 - 5
compton_filter.py

@@ -276,7 +276,7 @@ class CalibrdbHandler(PostgreSQLHandler):
             WHERE sid = %(sid)s AND comment = %(comment)s
             """
             for x in new_rows:
-                season_point = (comment if comment is not None else '') + '_' + str(x[3])
+                season_point = (comment if comment is not None else '') + '_' + str(x[4]) + '_' + str(x[3])
                 dict_row = {
                     'sid': sid,
                     'createdby': 'lxeuser',
@@ -289,9 +289,9 @@ class CalibrdbHandler(PostgreSQLHandler):
                 self.cur.execute(update_query, dict_row)
         
         insert_query = """INSERT INTO clbrdata (sid, createdby, time, begintime, endtime, comment, data) VALUES %s"""
-        comment_creator = lambda x: f'{comment if comment is not None else ""}_{str(x[3])}'
+        comment_creator = lambda x: f'{comment if comment is not None else ""}_{str(x[4])}_{str(x[3])}'
         insert_rows = list(map(lambda x: (sid, 'lxeuser', x[0], x[1], x[2], comment_creator(x), x[3:]), new_rows))        
-        execute_values(self.cur, insert_query, insert_rows, fetch=False)        
+        execute_values(self.cur, insert_query, insert_rows, fetch=False)  
         
         drop_query = f"""
             DELETE FROM clbrdata a
@@ -300,8 +300,9 @@ class CalibrdbHandler(PostgreSQLHandler):
                 a.sid = {sid}
                 AND a.cid > b.cid
                 AND a.sid = b.sid
-                AND a.comment = b.comment
+                AND a.begintime = b.begintime
         """
+        #AND a.comment = b.comment
         self.cur.execute(drop_query)
         
         logging.info(f"Inserted {len(insert_rows)} rows into table: {system}/{algo}/{name}/{version}")
@@ -398,4 +399,4 @@ def main():
     
 # python scripts/compton_filter.py --season cmd3_2021_2 --config database.ini --update
 if __name__ == "__main__":
-    main()
+    main()