# coding: utf-8 import progressbar import calendar import numpy as np from astropy.io import fits import pandas as pd from sqlalchemy import create_engine night_int = 20150721 def fjd(datetime_inst): """ convert datetime instance to FJD """ return calendar.timegm(datetime_inst.utctimetuple()) / (24.*3600.) def create_DB_connection(): from ConfigParser import SafeConfigParser config = SafeConfigParser() config.optionxform = str # this make the parsing case sensitive config.read('config.ini') factdb = create_engine( "mysql+mysqldb://{user}:{pw}@{host}/{db}".format( user=config.get('database', 'user'), pw=config.get('database', 'password'), host=config.get('database', 'host'), db=config.get('database', 'database'), ) ) return factdb def get_all_data_runs_from_run_db(factdb): keys = [ 'fRunID', 'fNight', 'fRunStart', 'fRunStop', 'fZenithDistanceMean', 'fAzimuthMean', ] sql_query = """SELECT {comma_sep_keys} FROM RunInfo WHERE fRunTypeKey=1 ORDER BY fNight; """ sql_query = sql_query.format( comma_sep_keys=', '.join(keys), ) data_runs = pd.read_sql_query( sql_query, factdb, parse_dates=['fRunStart', 'fRunStop'], ) return data_runs def get_trigger_rates(night_int, base_path='/fact/aux/'): night_string = str(night_int) fits_file = fits.open( base_path+'{y}/{m}/{d}/{n}.FTM_CONTROL_TRIGGER_RATES.fits'.format( n=night_string, y=night_string[0:4], m=night_string[4:6], d=night_string[6:8], ) ) trigger_rates = fits_file[1].data return trigger_rates def add_ratio_and_more_to_dataframe(data_runs): # create new columns and pre-assign some hopefully good NULL-like-values data_runs['number_of_measurements'] = 0 data_runs['median_of_board_rates'] = np.nan data_runs['std_dev_of_board_rates'] = np.nan data_runs['fBoardTriggerRateRatioAboveThreshold'] = np.nan trigger_rates = None last_night = None progress = progressbar.ProgressBar(widgets=[progressbar.Bar('=', '[', ']'), ' ', progressbar.Percentage(), ' ', progressbar.ETA()]) for dataframe_index in progress(data_runs.index): this_run = data_runs.ix[dataframe_index] if last_night != this_run['fNight']: try: trigger_rates = get_trigger_rates(this_run['fNight']) except (IOError, ValueError): continue last_night = this_run['fNight'] try: mask = ( (trigger_rates['Time'] > fjd(this_run['fRunStart'])) * (trigger_rates['Time'] < fjd(this_run['fRunStop'])) ) except ValueError: continue if mask.sum() == 0: continue try: this_board_rates = trigger_rates['BoardRate'][mask] except KeyError: continue N = this_board_rates.shape[0] data_runs.ix[dataframe_index, 'number_of_measurements'] = N left, med, right = np.percentile( this_board_rates, [50-20, 50, 50+20]) something_like_width = (right - left)/2. # 40% of the area of the normal distrubution # fall between -0.52*sigma and +0.52*sigma. # The estimation of sigma in a distorted CDF # is less affected near the median. (but of course less accurate as well.) std_dev = something_like_width / 0.52 median_board_rates = np.median(this_board_rates, axis=1) over = (median_board_rates > med+3*std_dev).sum() data_runs.ix[dataframe_index, 'median_of_board_rates'] = med data_runs.ix[dataframe_index, 'std_dev_of_board_rates'] = std_dev data_runs.ix[dataframe_index, 'fBoardTriggerRateRatioAboveThreshold'] = float(over)/N def main(night_int): factdb = create_DB_connection() data_runs = get_all_data_runs_from_run_db(factdb) add_ratio_and_more_to_dataframe(data_runs) return data_runs if __name__ == '__main__': data_runs = main(night_int) data_runs.to_csv('burst_ratio.csv')