latools/latools.py
def calibrate(self, analytes=None, drift_correct=True,
                  srms_used=['NIST610', 'NIST612', 'NIST614'],
                  zero_intercept=True, n_min=10, reload_srm_database=False):
        """
        Calibrates the data to measured SRM values.

        Assumes that y intercept is zero.

        Parameters
        ----------  
        analytes : str or iterable
            Which analytes you'd like to calibrate. Defaults to all.
        drift_correct : bool
            Whether to pool all SRM measurements into a single calibration,
            or vary the calibration through the run, interpolating
            coefficients between measured SRMs.
        srms_used : str or iterable
            Which SRMs have been measured. Must match names given in
            SRM data file *exactly*.
        n_min : int
            The minimum number of data points an SRM measurement
            must have to be included.

        Returns
        -------
        None
        """
        if analytes is None:
            analytes = self.analytes[self.analytes != self.internal_standard]
        elif isinstance(analytes, str):
            analytes = [analytes]

        if not hasattr(self, 'srmtabs'):
            self.srm_id_auto(srms_used=srms_used, n_min=n_min, reload_srm_database=reload_srm_database)

        # make container for calibration params
        if not hasattr(self, 'calib_params'):
            gTime = self.stdtab.gTime.unique()
            self.calib_params = pd.DataFrame(columns=pd.MultiIndex.from_product([analytes, ['m']]),
                                            index=gTime)

        calib_analytes = self.srmtabs.index.get_level_values(0).unique()

        if zero_intercept:
            fn  = lambda x, m: x * m
        else:
            fn = lambda x, m, c: x * m + c

        for a in calib_analytes:
            if zero_intercept:
                if (a, 'c') in self.calib_params:
                    self.calib_params.drop((a, 'c'), 1, inplace=True)
            if drift_correct:
                for g in self.stdtab.gTime.unique():
                    ind = idx[a, :, :, g]
                    if self.srmtabs.loc[ind].size == 0:
                        continue
                    # try:
                    meas = self.srmtabs.loc[ind, 'meas_mean']
                    srm = self.srmtabs.loc[ind, 'srm_mean']
                    # TODO: replace curve_fit with Sambridge's 2D likelihood function for better uncertainty incorporation.
                    merr = self.srmtabs.loc[ind, 'meas_err']
                    serr = self.srmtabs.loc[ind, 'srm_err']
                    sigma = np.sqrt(merr**2 + serr**2)

                    if len(meas) > 1:
                        # multiple SRMs - do a regression
                        p, cov = curve_fit(fn, meas, srm, sigma=sigma)
                        pe = unc.correlated_values(p, cov)                
                        self.calib_params.loc[g, (a, 'm')] = pe[0]
                        if not zero_intercept:
                            self.calib_params.loc[g, (a, 'c')] = pe[1]
                    else:
                        # deal with case where there's only one datum
                        self.calib_params.loc[g, (a, 'm')] = (un.uarray(srm, serr) / 
                                                              un.uarray(meas, merr))[0]
                        if not zero_intercept:
                            self.calib_params.loc[g, (a, 'c')] = 0

                    # This should be obsolete, because no-longer sourcing locator from calib_params index.
                    # except KeyError:
                    #     # If the calibration is being recalculated, calib_params
                    #     # will have t=0 and t=max(uTime) values that are outside
                    #     # the srmtabs index.
                    #     # If this happens, drop them, and re-fill them at the end.
                    #     self.calib_params.drop(g, inplace=True)
            else:
                ind = idx[a, :, :, :]
                meas = self.srmtabs.loc[ind, 'meas_mean']
                srm = self.srmtabs.loc[ind, 'srm_mean']
                merr = self.srmtabs.loc[ind, 'meas_err']
                serr = self.srmtabs.loc[ind, 'srm_err']
                sigma = np.sqrt(merr**2 + serr**2)
                
                if len(meas) > 1:
                    p, cov = curve_fit(fn, meas, srm, sigma=sigma)
                    pe = unc.correlated_values(p, cov)                
                    self.calib_params.loc[:, (a, 'm')] = pe[0]
                    if not zero_intercept:
                        self.calib_params.loc[:, (a, 'c')] = pe[1]
                else:
                    self.calib_params.loc[:, (a, 'm')] = (un.uarray(srm, serr) / 
                                                          un.uarray(meas, merr))[0]
                    if not zero_intercept:
                        self.calib_params.loc[:, (a, 'c')] = 0

        # if fill:
        # fill in uTime=0 and uTime = max cases for interpolation
        if self.calib_params.index.min() == 0:
            self.calib_params.drop(0, inplace=True)
            self.calib_params.drop(self.calib_params.index.max(), inplace=True)
        self.calib_params.loc[0, :] = self.calib_params.loc[self.calib_params.index.min(), :]
        maxuT = np.max([d.uTime.max() for d in self.data.values()])  # calculate max uTime
        self.calib_params.loc[maxuT, :] = self.calib_params.loc[self.calib_params.index.max(), :]
        # sort indices for slice access
        self.calib_params.sort_index(1, inplace=True)
        self.calib_params.sort_index(0, inplace=True)

        # calculcate interpolators for applying calibrations
        self.calib_ps = Bunch()
        for a in analytes:
            # TODO: revisit un_interp1d to see whether it plays well with correlated values. 
            # Possible re-write to deal with covariance matrices?
            self.calib_ps[a] = {'m': un_interp1d(self.calib_params.index.values,
                                                self.calib_params.loc[:, (a, 'm')].values)}
            if not zero_intercept:
                self.calib_ps[a]['c'] = un_interp1d(self.calib_params.index.values,
                                                    self.calib_params.loc[:, (a, 'c')].values)

        with self.pbar.set(total=len(self.data), desc='Applying Calibrations') as prog:
            for d in self.data.values():
                d.calibrate(self.calib_ps, analytes)
                prog.update()

        # record SRMs used for plotting
        markers = 'osDsv<>PX'  # for future implementation of SRM-specific markers.
        if not hasattr(self, 'srms_used'):
            self.srms_used = set(srms_used)
        else:
            self.srms_used.update(srms_used)
        self.srm_mdict = {k: markers[i] for i, k in enumerate(self.srms_used)}

        self.stages_complete.update(['calibrated'])
        self.focus_stage = 'calibrated'

        return
def calibrate(self, analytes=None, drift_correct=True,
                  srms_used=['NIST610', 'NIST612', 'NIST614'],
                  zero_intercept=True, n_min=10, reload_srm_database=False):
        """
        Calibrates the data to measured SRM values.

        Assumes that y intercept is zero.

        Parameters
        ----------  
        analytes : str or iterable
            Which analytes you'd like to calibrate. Defaults to all.
        drift_correct : bool
            Whether to pool all SRM measurements into a single calibration,
            or vary the calibration through the run, interpolating
            coefficients between measured SRMs.
        srms_used : str or iterable
            Which SRMs have been measured. Must match names given in
            SRM data file *exactly*.
        n_min : int
            The minimum number of data points an SRM measurement
            must have to be included.

        Returns
        -------
        None
        """
        if analytes is None:
            analytes = self.analytes[self.analytes != self.internal_standard]
        elif isinstance(analytes, str):
            analytes = [analytes]

        if not hasattr(self, 'srmtabs'):
            self.srm_id_auto(srms_used=srms_used, n_min=n_min, reload_srm_database=reload_srm_database)

        # make container for calibration params
        if not hasattr(self, 'calib_params'):
            gTime = self.stdtab.gTime.unique()
            self.calib_params = pd.DataFrame(columns=pd.MultiIndex.from_product([analytes, ['m']]),
                                            index=gTime)

        calib_analytes = self.srmtabs.index.get_level_values(0).unique()

        if zero_intercept:
            fn  = lambda x, m: x * m
        else:
            fn = lambda x, m, c: x * m + c

        for a in calib_analytes:
            if zero_intercept:
                if (a, 'c') in self.calib_params:
                    self.calib_params.drop((a, 'c'), 1, inplace=True)
            if drift_correct:
                for g in self.stdtab.gTime.unique():
                    ind = idx[a, :, :, g]
                    if self.srmtabs.loc[ind].size == 0:
                        continue
                    # try:
                    meas = self.srmtabs.loc[ind, 'meas_mean']
                    srm = self.srmtabs.loc[ind, 'srm_mean']
                    # TODO: replace curve_fit with Sambridge's 2D likelihood function for better uncertainty incorporation.
                    merr = self.srmtabs.loc[ind, 'meas_err']
                    serr = self.srmtabs.loc[ind, 'srm_err']
                    sigma = np.sqrt(merr**2 + serr**2)

                    if len(meas) > 1:
                        # multiple SRMs - do a regression
                        p, cov = curve_fit(fn, meas, srm, sigma=sigma)
                        pe = unc.correlated_values(p, cov)                
                        self.calib_params.loc[g, (a, 'm')] = pe[0]
                        if not zero_intercept:
                            self.calib_params.loc[g, (a, 'c')] = pe[1]
                    else:
                        # deal with case where there's only one datum
                        self.calib_params.loc[g, (a, 'm')] = (un.uarray(srm, serr) / 
                                                              un.uarray(meas, merr))[0]
                        if not zero_intercept:
                            self.calib_params.loc[g, (a, 'c')] = 0

                    # This should be obsolete, because no-longer sourcing locator from calib_params index.
                    # except KeyError:
                    #     # If the calibration is being recalculated, calib_params
                    #     # will have t=0 and t=max(uTime) values that are outside
                    #     # the srmtabs index.
                    #     # If this happens, drop them, and re-fill them at the end.
                    #     self.calib_params.drop(g, inplace=True)
            else:
                ind = idx[a, :, :, :]
                meas = self.srmtabs.loc[ind, 'meas_mean']
                srm = self.srmtabs.loc[ind, 'srm_mean']
                merr = self.srmtabs.loc[ind, 'meas_err']
                serr = self.srmtabs.loc[ind, 'srm_err']
                sigma = np.sqrt(merr**2 + serr**2)
                
                if len(meas) > 1:
                    p, cov = curve_fit(fn, meas, srm, sigma=sigma)
                    pe = unc.correlated_values(p, cov)                
                    self.calib_params.loc[:, (a, 'm')] = pe[0]
                    if not zero_intercept:
                        self.calib_params.loc[:, (a, 'c')] = pe[1]
                else:
                    self.calib_params.loc[:, (a, 'm')] = (un.uarray(srm, serr) / 
                                                          un.uarray(meas, merr))[0]
                    if not zero_intercept:
                        self.calib_params.loc[:, (a, 'c')] = 0

        # if fill:
        # fill in uTime=0 and uTime = max cases for interpolation
        if self.calib_params.index.min() == 0:
            self.calib_params.drop(0, inplace=True)
            self.calib_params.drop(self.calib_params.index.max(), inplace=True)
        self.calib_params.loc[0, :] = self.calib_params.loc[self.calib_params.index.min(), :]
        maxuT = np.max([d.uTime.max() for d in self.data.values()])  # calculate max uTime
        self.calib_params.loc[maxuT, :] = self.calib_params.loc[self.calib_params.index.max(), :]
        # sort indices for slice access
        self.calib_params.sort_index(1, inplace=True)
        self.calib_params.sort_index(0, inplace=True)

        # calculcate interpolators for applying calibrations
        self.calib_ps = Bunch()
        for a in analytes:
            # TODO: revisit un_interp1d to see whether it plays well with correlated values. 
            # Possible re-write to deal with covariance matrices?
            self.calib_ps[a] = {'m': un_interp1d(self.calib_params.index.values,
                                                self.calib_params.loc[:, (a, 'm')].values)}
            if not zero_intercept:
                self.calib_ps[a]['c'] = un_interp1d(self.calib_params.index.values,
                                                    self.calib_params.loc[:, (a, 'c')].values)

        with self.pbar.set(total=len(self.data), desc='Applying Calibrations') as prog:
            for d in self.data.values():
                d.calibrate(self.calib_ps, analytes)
                prog.update()

        # record SRMs used for plotting
        markers = 'osDsv<>PX'  # for future implementation of SRM-specific markers.
        if not hasattr(self, 'srms_used'):
            self.srms_used = set(srms_used)
        else:
            self.srms_used.update(srms_used)
        self.srm_mdict = {k: markers[i] for i, k in enumerate(self.srms_used)}

        self.stages_complete.update(['calibrated'])
        self.focus_stage = 'calibrated'

        return
def calibrate ( self , analytes = None , drift_correct = True , srms_used = [ 'NIST610' , 'NIST612' , 'NIST614' ] , zero_intercept = True , n_min = 10 , reload_srm_database = False ) : if analytes is None : analytes = self . analytes [ self . analytes != self . internal_standard ] elif isinstance ( analytes , str ) : analytes = [ analytes ] if not hasattr ( self , 'srmtabs' ) : self . srm_id_auto ( srms_used = srms_used , n_min = n_min , reload_srm_database = reload_srm_database ) # make container for calibration params if not hasattr ( self , 'calib_params' ) : gTime = self . stdtab . gTime . unique ( ) self . calib_params = pd . DataFrame ( columns = pd . MultiIndex . from_product ( [ analytes , [ 'm' ] ] ) , index = gTime ) calib_analytes = self . srmtabs . index . get_level_values ( 0 ) . unique ( ) if zero_intercept : fn = lambda x , m : x * m else : fn = lambda x , m , c : x * m + c for a in calib_analytes : if zero_intercept : if ( a , 'c' ) in self . calib_params : self . calib_params . drop ( ( a , 'c' ) , 1 , inplace = True ) if drift_correct : for g in self . stdtab . gTime . unique ( ) : ind = idx [ a , : , : , g ] if self . srmtabs . loc [ ind ] . size == 0 : continue # try: meas = self . srmtabs . loc [ ind , 'meas_mean' ] srm = self . srmtabs . loc [ ind , 'srm_mean' ] # TODO: replace curve_fit with Sambridge's 2D likelihood function for better uncertainty incorporation. merr = self . srmtabs . loc [ ind , 'meas_err' ] serr = self . srmtabs . loc [ ind , 'srm_err' ] sigma = np . sqrt ( merr ** 2 + serr ** 2 ) if len ( meas ) > 1 : # multiple SRMs - do a regression p , cov = curve_fit ( fn , meas , srm , sigma = sigma ) pe = unc . correlated_values ( p , cov ) self . calib_params . loc [ g , ( a , 'm' ) ] = pe [ 0 ] if not zero_intercept : self . calib_params . loc [ g , ( a , 'c' ) ] = pe [ 1 ] else : # deal with case where there's only one datum self . calib_params . loc [ g , ( a , 'm' ) ] = ( un . uarray ( srm , serr ) / un . uarray ( meas , merr ) ) [ 0 ] if not zero_intercept : self . calib_params . loc [ g , ( a , 'c' ) ] = 0 # This should be obsolete, because no-longer sourcing locator from calib_params index. # except KeyError: #     # If the calibration is being recalculated, calib_params #     # will have t=0 and t=max(uTime) values that are outside #     # the srmtabs index. #     # If this happens, drop them, and re-fill them at the end. #     self.calib_params.drop(g, inplace=True) else : ind = idx [ a , : , : , : ] meas = self . srmtabs . loc [ ind , 'meas_mean' ] srm = self . srmtabs . loc [ ind , 'srm_mean' ] merr = self . srmtabs . loc [ ind , 'meas_err' ] serr = self . srmtabs . loc [ ind , 'srm_err' ] sigma = np . sqrt ( merr ** 2 + serr ** 2 ) if len ( meas ) > 1 : p , cov = curve_fit ( fn , meas , srm , sigma = sigma ) pe = unc . correlated_values ( p , cov ) self . calib_params . loc [ : , ( a , 'm' ) ] = pe [ 0 ] if not zero_intercept : self . calib_params . loc [ : , ( a , 'c' ) ] = pe [ 1 ] else : self . calib_params . loc [ : , ( a , 'm' ) ] = ( un . uarray ( srm , serr ) / un . uarray ( meas , merr ) ) [ 0 ] if not zero_intercept : self . calib_params . loc [ : , ( a , 'c' ) ] = 0 # if fill: # fill in uTime=0 and uTime = max cases for interpolation if self . calib_params . index . min ( ) == 0 : self . calib_params . drop ( 0 , inplace = True ) self . calib_params . drop ( self . calib_params . index . max ( ) , inplace = True ) self . calib_params . loc [ 0 , : ] = self . calib_params . loc [ self . calib_params . index . min ( ) , : ] maxuT = np . max ( [ d . uTime . max ( ) for d in self . data . values ( ) ] ) # calculate max uTime self . calib_params . loc [ maxuT , : ] = self . calib_params . loc [ self . calib_params . index . max ( ) , : ] # sort indices for slice access self . calib_params . sort_index ( 1 , inplace = True ) self . calib_params . sort_index ( 0 , inplace = True ) # calculcate interpolators for applying calibrations self . calib_ps = Bunch ( ) for a in analytes : # TODO: revisit un_interp1d to see whether it plays well with correlated values.  # Possible re-write to deal with covariance matrices? self . calib_ps [ a ] = { 'm' : un_interp1d ( self . calib_params . index . values , self . calib_params . loc [ : , ( a , 'm' ) ] . values ) } if not zero_intercept : self . calib_ps [ a ] [ 'c' ] = un_interp1d ( self . calib_params . index . values , self . calib_params . loc [ : , ( a , 'c' ) ] . values ) with self . pbar . set ( total = len ( self . data ) , desc = 'Applying Calibrations' ) as prog : for d in self . data . values ( ) : d . calibrate ( self . calib_ps , analytes ) prog . update ( ) # record SRMs used for plotting markers = 'osDsv<>PX' # for future implementation of SRM-specific markers. if not hasattr ( self , 'srms_used' ) : self . srms_used = set ( srms_used ) else : self . srms_used . update ( srms_used ) self . srm_mdict = { k : markers [ i ] for i , k in enumerate ( self . srms_used ) } self . stages_complete . update ( [ 'calibrated' ] ) self . focus_stage = 'calibrated' return
Calibrates the data to measured SRM values.

        Assumes that y intercept is zero.

        Parameters
        ----------  
        analytes : str or iterable
            Which analytes you'd like to calibrate. Defaults to all.
        drift_correct : bool
            Whether to pool all SRM measurements into a single calibration,
            or vary the calibration through the run, interpolating
            coefficients between measured SRMs.
        srms_used : str or iterable
            Which SRMs have been measured. Must match names given in
            SRM data file *exactly*.
        n_min : int
            The minimum number of data points an SRM measurement
            must have to be included.

        Returns
        -------
        None
Calibrates the data to measured SRM values .