This file describes the software modules used for SOIR archiving. (written by E. NEEFS) 1. NON-LINEARITY CORRECTION SOFTWARE ------------------------------------ (bring internal SOIR level 0.1 data to internal SOIR level 0.2 data) ###################################################################### ## Project # SOIR # ## File # produce_lev0_2/non_lin_correction.py # ###################################################################### ## Author # Roland Clairquin # ## E-Mail # rolandcl@oma.be # ## Company # BIRA-IASB # ###################################################################### def adccode_to_charge( x ): if x < 6000: y=-109.4112717552833+x*(0.3281672408563101+\ x*(-0.0003846513541535442+x*(2.869226627796301E-07+\ x*(-1.381722060516796E-10+x*(4.459643046851159E-14+\ x*(-9.752279474228916E-18+x*(1.426792904826683E-21+\ x*(-1.337703563748429E-25+x*(7.266297806363216E-30+\ x*-1.738835026549852E-34))))))))) else: y = 6.0634764 + 0.02184421*x return y bkgcodes_vs_integtime = [663,663,679,693, 706,721,738,755,772,790,808,827,846,866,886,908,930,952,975,1000, 1024,1050,1077,1104,1134,1164,1194,1225,1257,1289,1323,1357,1391, 1427,1463,1500,1536,1574,1611,1650,1688,1727,1766,1806,1846,1886, 1926,1966,2008,2048,2089,2131,2173,2215,2257,2299,2340,2383,2426, 2469,2511,2555,2599,2641,2684,2729,2772,2815,2860,2903,2947,2992, 3035,3080,3125,3168,3213,3257,3302,3346,3391,3437,3481,3527,3572, 3616,3661,3706,3752,3797,3842,3887,3933,3977,4022,4068,4113,4159, 4205,4250,4296,4342,4387,4432,4479,4524,4570,4616,4661,4707,4753, 4799,4844,4891,4936,4982,5028,5075,5121,5166,5212,5259,5305,5350, 5396,5442,5488,5534,5581,5627,5672,5719,5765,5811,5858,5903,5950, 5995,6042,6088,6134,6182,6227,6274,6319,6366,6412,6458,6504,6551, 6597] ... for line in ifd: ## If we are in the data zone ############################################################ if data_zone_flag: data = line.split(',') if len(data) < 9: data_zone_flag = False ofd.write( line ) else: ## Write the first column (unchanged) ofd.write( data[0] + ',' ) for i in range(1,9): order_ix = (i-1) / bins_per_order adc_bkg = non_lin_correction.bkgcodes_vs_integtime [ integ_times_ms[order_ix] ] n_accum = (dcbf+1) * ( nracs[order_ix]-1 ) / 2 raw_value = float(data[i] ## Corrects a possible overflow of the 20-bit accumulator ## in the FPGA if raw_value < -50000.0: raw_value += 1048576.0 ## Correction of the detector non-linearity adc_sig = float(data[i]) / n_accum ads_sig_plus_bkg = adc_sig + adc_bkg charge = non_lin_correction.adccode_to_charge( ads_sig_plus_bkg ) - integ_times_ms[order_ix] ## Write the corrected pixel ofd.write( "%0.6f," % (charge) ) ofd.write("\n") ... 2. WAVENUMBER CORRECTION SOFTWARE --------------------------------- ## Corrections from Arnaud Mahieux (E-mail 13/12/2007) def wn_corr_function( binning_fact, n_of_bins, str_date ): corrections = [ [] for i in range(n_of_bins) ] if binning_fact == 4: if str_date < "20070804": # pixVSpixcorr_fullscan_140_binning_04.dat # 9.3251286611e+0 1.0589168767e+0 -1.1283944205e-4 4.2612617637e-7 # 9.8349173816e+0 1.0472960365e+0 -2.8012841041e-5 2.4683001803e-7 # 9.4664850174e+0 1.0580933699e+0 -9.9878156613e-5 3.8249107023e-7 # 9.6655591559e+0 1.0565595264e+0 -8.4763534379e-5 3.4485635641e-7 # 9.9936245729e+0 1.0511133316e+0 -6.1668396712e-5 3.2548128607e-7 # 1.0448645788e+1 1.0386135544e+0 4.5452379099e-5 6.6899950885e-8 if n_of_bins == 8: for x in range(320): val = 9.3251286611 + 1.0589168767 * x - 1.1283944205e-004 * x**2 + 4.2612617637e-007 * x**3 corrections[0].append( val ) val = 9.3251286611 + 1.0589168767 * x - 1.1283944205e-004 * x**2 + 4.2612617637e-007 * x**3 corrections[1].append( val ) val = 9.8349173816 + 1.0472960365 * x - 2.8012841041e-005 * x**2 + 2.4683001803e-007 * x**3 corrections[2].append( val ) val = 9.4664850174 + 1.0580933699 * x - 9.9878156613e-005 * x**2 + 3.8249107023e-007 * x**3 corrections[3].append( val ) val = 9.6655591559 + 1.0565595264 * x - 8.4763534379e-005 * x**2 + 3.4485635641e-007 * x**3 corrections[4].append( val ) val = 9.9936245729 + 1.0511133316 * x - 6.1668396712e-005 * x**2 + 3.2548128607e-007 * x**3 corrections[5].append( val ) val = 10.448645788 + 1.0386135544 * x + 4.5452379099e-005 * x**2 + 6.6899950885e-008 * x**3 corrections[6].append( val ) val = 10.448645788 + 1.0386135544 * x + 4.5452379099e-005 * x**2 + 6.6899950885e-008 * x**3 corrections[7].append( val ) elif n_of_bins == 4: for x in range(320): val = 9.8349173816 + 1.0472960365 * x - 2.8012841041e-005 * x**2 + 2.4683001803e-007 * x**3 corrections[0].append( val ) val = 9.4664850174 + 1.0580933699 * x - 9.9878156613e-005 * x**2 + 3.8249107023e-007 * x**3 corrections[1].append( val ) val = 9.6655591559 + 1.0565595264 * x - 8.4763534379e-005 * x**2 + 3.4485635641e-007 * x**3 corrections[2].append( val ) val = 9.9936245729 + 1.0511133316 * x - 6.1668396712e-005 * x**2 + 3.2548128607e-007 * x**3 corrections[3].append( val ) else: # pixVSpixcorr_fullscan_470_binning_04.dat # 1.2495153766e+001 1.0467746258e+000 -2.6829518463e-005 2.5356765471e-007 # 1.2426527569e+001 1.0464762560e+000 -1.1630159206e-005 1.9959751196e-007 # 1.2526660570e+001 1.0430118883e+000 1.7858397790e-005 1.2457837106e-007 # 1.2043549187e+001 1.0537009941e+000 -5.9477423843e-005 2.9492967477e-007 # 1.2005238771e+001 1.0564214556e+000 -8.0628266962e-005 3.2619518389e-007 # 1.2057862214e+001 1.0553827878e+000 -7.6692983298e-005 3.2429320975e-007 if n_of_bins == 8: for x in range(320): val = 12.495153766 + 1.0467746258 * x - 2.6829518463e-005 * x**2 + 2.5356765471e-007 * x**3 corrections[0].append( val ) val = 12.495153766 + 1.0467746258 * x - 2.6829518463e-005 * x**2 + 2.5356765471e-007 * x**3 corrections[1].append( val ) val = 12.426527569 + 1.0464762560 * x - 1.1630159206e-005 * x**2 + 1.9959751196e-007 * x**3 corrections[2].append( val ) val = 12.526660570 + 1.0430118883 * x + 1.7858397790e-005 * x**2 + 1.2457837106e-007 * x**3 corrections[3].append( val ) val = 12.043549187 + 1.0537009941 * x - 5.9477423843e-005 * x**2 + 2.9492967477e-007 * x**3 corrections[4].append( val ) val = 12.005238771 + 1.0564214556 * x - 8.0628266962e-005 * x**2 + 3.2619518389e-007 * x**3 corrections[5].append( val ) val = 12.057862214 + 1.0553827878 * x - 7.6692983298e-005 * x**2 + 3.2429320975e-007 * x**3 corrections[6].append( val ) val = 12.057862214 + 1.0553827878 * x - 7.6692983298e-005 * x**2 + 3.2429320975e-007 * x**3 corrections[7].append( val ) elif n_of_bins == 4: for x in range(320): val = 12.426527569 + 1.0464762560 * x - 1.1630159206e-005 * x**2 + 1.9959751196e-007 * x**3 corrections[0].append( val ) val = 12.526660570 + 1.0430118883 * x + 1.7858397790e-005 * x**2 + 1.2457837106e-007 * x**3 corrections[1].append( val ) val = 12.043549187 + 1.0537009941 * x - 5.9477423843e-005 * x**2 + 2.9492967477e-007 * x**3 corrections[2].append( val ) val = 12.005238771 + 1.0564214556 * x - 8.0628266962e-005 * x**2 + 3.2619518389e-007 * x**3 corrections[3].append( val ) elif binning_fact == 12: if str_date < "20060908": # pixVSpixcorr_fullscan_106_binning_12.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) elif str_date < "20061116": # pixVSpixcorr_fullscan_140_binning_12.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) elif str_date < "20070110": # pixVSpixcorr_fullscan_209_binning_12.dat # 1.0891223136e+1 1.0468377964e+0 -2.1635349008e-5 2.3014038635e-7 # 1.0884423892e+1 1.0495689756e+0 -3.2619055528e-5 2.4004439949e-7 for x in range(320): val = 10.891223136 + 1.0468377964 * x - 2.1635349008e-005 * x**2 + 2.3014038635e-007 * x**3 corrections[0].append( val ) val = 10.884423892 + 1.0495689756 * x - 3.2619055528e-005 * x**2 + 2.4004439949e-007 * x**3 corrections[1].append( val ) elif str_date < "20070804": # pixVSpixcorr_fullscan_264_binning_12.dat # 1.0891223136e+1 1.0468377964e+0 -2.1635349008e-5 2.3014038635e-7 # 1.0884423892e+1 1.0495689756e+0 -3.2619055528e-5 2.4004439949e-7 for x in range(320): val = 10.891223136 + 1.0468377964 * x - 2.1635349008e-005 * x**2 + 2.3014038635e-007 * x**3 corrections[0].append( val ) val = 10.884423892 + 1.0495689756 * x - 3.2619055528e-005 * x**2 + 2.4004439949e-007 * x**3 corrections[1].append( val ) else: # pixVSpixcorr_fullscan_470_binning_12.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) elif binning_fact == 16: if str_date < "20060908": # pixVSpixcorr_fullscan_106_binning_16.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) if str_date < "20061116": # pixVSpixcorr_fullscan_140_binning_16.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) if str_date < "20070110": # pixVSpixcorr_fullscan_209_binning_16.dat # 1.0891223136e+1 1.0468377964e+0 -2.1635349008e-5 2.3014038635e-7 # 1.0884423892e+1 1.0495689756e+0 -3.2619055528e-5 2.4004439949e-7 for x in range(320): val = 10.891223136 + 1.0468377964 * x - 2.1635349008e-005 * x**2 + 2.3014038635e-007 * x**3 corrections[0].append( val ) val = 10.884423892 + 1.0495689756 * x - 3.2619055528e-005 * x**2 + 2.4004439949e-007 * x**3 corrections[1].append( val ) if str_date < "20070804": # pixVSpixcorr_fullscan_264_binning_16.dat # 1.0891223136e+1 1.0468377964e+0 -2.1635349008e-5 2.3014038635e-7 # 1.0884423892e+1 1.0495689756e+0 -3.2619055528e-5 2.4004439949e-7 for x in range(320): val = 10.891223136 + 1.0468377964 * x - 2.1635349008e-005 * x**2 + 2.3014038635e-007 * x**3 corrections[0].append( val ) val = 10.884423892 + 1.0495689756 * x - 3.2619055528e-005 * x**2 + 2.4004439949e-007 * x**3 corrections[1].append( val ) else: # pixVSpixcorr_fullscan_470_binning_16.dat # 1.1581630262e+1 1.0473801863e+0 -2.4637866889e-5 2.3680104252e-7 # 1.1394038164e+1 1.0545550974e+0 -6.6628118196e-5 3.0043205490e-7 for x in range(320): val = 11.581630262 + 1.0473801863 * x - 2.4637866889e-005 * x**2 + 2.3680104252e-007 * x**3 corrections[0].append( val ) val = 11.394038164 + 1.0545550974 * x - 6.6628118196e-005 * x**2 + 3.0043205490e-007 * x**3 corrections[1].append( val ) else: raise RuntimeError, "No corrections data" return corrections ... from wn_corr_function import wn_corr_function ## Spectra Viewer Corrections Contants Le_MINUS_OFFSET =2274.75 Le_PLUS_OFFSET =2256.41 Le_MINUS_SLOPE =22.52135922 Le_PLUS_SLOPE =22.34019417 ORDER_BASE =101 class WavenumberCorr: def __init__( self, tm, str_date ): binning_fact, nof_bins = tm.get_binning() self.corrections = wn_corr_function ( binning_fact,nof_bins,str_date ) def get_wavenumbers( self, order ): Le_minus = (order - ORDER_BASE) * Le_MINUS_SLOPE + Le_MINUS_OFFSET; Le_plus = (order - ORDER_BASE) * Le_PLUS_SLOPE + Le_PLUS_OFFSET; slope = (Le_minus - Le_plus) / 319; wavenumbers = [] for corr_bin in self.corrections: wavenumbers.append(map( lambda x: x * slope + Le_plus,corr_bin )) return wavenumbers ... 3. FULL SUN REFERENCING SOFTWARE -------------------------------- ## Read the TMs in the reference zone ###################################################################### regr_data_dico = {} regr_flag = False tm_ix = 0 for file_name in csv_files: if file_name[2:16] == beg_regr: regr_flag = True if regr_flag: tm = TmCsv( source_dir + '/' + file_name ) aofs_list = tm.get_aofs_list() spectrum_list = tm.get_spectrum_list() for aofs, spectrum in zip( aofs_list, spectrum_list ): if not regr_data_dico.has_key( aofs ): n_of_bins = spectrum.get_nof_bins() regr_data_dico[aofs] = [[([],[]) for i in range(320)] for j in range(n_of_bins)] curve_tab = regr_data_dico[aofs] for bin_ix in range(n_of_bins): for pix_ix in range(320): curve_tab[bin_ix][pix_ix][0].append( float(tm_ix) ) curve_tab[bin_ix][pix_ix][1].append( spectrum[pix_ix][bin_ix]) ... ## Execute the linear regression for each pixel ###################################################################### failed = False statistics = [] regr_dico = copy.deepcopy( regr_data_dico ) for aofs in regr_dico.keys(): for bin_ix in range(n_of_bins): for pix_ix in range(320): try: xs, ys = regr_dico[aofs][bin_ix][pix_ix] coeff_a, coeff_b = linear_regression( xs, ys ) regr_dico[aofs][bin_ix][pix_ix] = (coeff_a, coeff_b) if pix_ix == 150: cc = abs( corr_coeff( xs, ys ) ) dy = delta_y(xs, ys) err_rel = dy / moy(ys) statistics.append( (order_names_dico[aofs], cc, dy, err_rel)) except ZeroDivisionError: failed = True if failed: self.write_report( submeas_name, "Une erreur est survenue lors de la regression lineaire" ) continue ... ###################################################################### ## Produce the .csv corrected data ###################################################################### ## For each TM ###################################################################### occult_flag = False tm_ix = 0 for file_name in csv_files: ## Si premiere TM de la zone occultation if file_name[2:16] == beg_occult: occult_flag = True ## If we are in the occultation zone ################################################################## if occult_flag: tm = TmCsv( source_dir + '/' + file_name ) aofs_list = tm.get_aofs_list() order_name_list = map( lambda x: order_names_dico[x], aofs_list ) wavenumbers_list = map( lambda x: wn_correction.get_wavenumbers( orders_dico[x] ), aofs_list ) regr_list = map( lambda x: regr_dico[x], aofs_list ) filename_list = map( lambda x: "%s/%s_order%s.csv" % ( dest_dir, file_name[:-4], x), order_name_list ) tm.write_corrected( filename_list, regr_list, tm_ix, order_name_list, wavenumbers_list, bad_pixels, log_infos, histo ) ..