Python debug_.debug_函数代码示例

 def __init__(self, svd_data, sv_list, timebase, min_dphase = -np.pi, phase_pairs = None):
     # NOTE I'd prefer not to duplicate info here which is in svd_data - should be able to refer to that, once sqlalchemy is hooked in
     #self.topo_channels = svd_data.topo_channels
     self.channels = svd_data.channels
     #self.svs = sv_list
     if len(sv_list) == 1:
         self.a12 = 0
     else:
         # this ratio was inverted accidentally at first
         self.a12 = svd_data.svs[sv_list[1]]/svd_data.svs[sv_list[0]]
     self._binary_svs = list2bin(sv_list)
     # peak frequency for fluctuation structure
     self.timebase = timebase
     self.freq, self.freq_elmt = peak_freq(svd_data.chronos[sv_list[0]], self.timebase)
     self.phase_pairs = phase_pairs
     self.t0 = timebase[0]
     # singular value filtered signals
     self.signal = np.dot(np.transpose(svd_data.topos[sv_list,:]),
                        np.dot(np.diag(svd_data.svs.take(sv_list)), svd_data.chronos[sv_list,:]))
     # phase differences between nearest neighbour channels
     self.dphase, self.fourier_values = self._get_dphase(min_dphase=min_dphase, get_fourier_value = 1)
     self.p = np.sum(svd_data.svs.take(sv_list)**2)/svd_data.E
     self.H = svd_data.H
     self.E = svd_data.E
     debug_(pyfusion.DEBUG, 4, key='FlucStruc')
     super(FlucStruc, self).__init__()

def regenerate_dim(x):
    """ assume x in ns since epoch from the current time """
    msg = None  # msg allows us to see which shot/diag was at fault
    diffs = np.diff(x)
    # bincount needs a positive input and needs an array with N elts where N is the largest number input
    small = (diffs > 0) & (diffs < 1000000)
    sorted_diffs = np.sort(diffs[np.where(small)[0]])
    counts = np.bincount(sorted_diffs)
    bigcounts, bigvals = myhist(diffs[np.where(~small)[0]])

    if pyfusion.VERBOSE>0:
        print('[[diff, count],....]')
        print('small:', [[argc, counts[argc]] for argc in np.argsort(counts)[::-1][0:5]])
        print('big or negative:', [[bigvals[argc], bigcounts[argc]] for argc in np.argsort(bigcounts)[::-1][0:10]])

    dtns = 1 + np.argmax(counts[1:])  # skip the first position - it is 0
    # wgt0 = np.where(sorted_diffs > 0)[0]  # we are in ns, so no worry about rounding
    histo = plt.hist if pyfusion.DBG() > 1 else np.histogram
    cnts, vals = histo(x, bins=200)[0:2]
    # ignore the two end bins - hopefully there will be very few there
    wmin = np.where(cnts[1:-1] < np.max(cnts[1:-1]))[0]
    if len(wmin)>0:
        print('**********\n*********** Gap in data > {p:.2f}%'.format(p=100*len(wmin)/float(len(cnts))))
    x01111 = np.ones(len(x))  # x01111 will be all 1s except for the first elt.
    x01111[0] = 0
    errcnt = np.sum(bigcounts) + np.sum(np.sort(counts)[::-1][1:])
    if errcnt>0 or (pyfusion.VERBOSE > 0): 
        msg = str('** repaired length of {l:,}, dtns={dtns:,}, {e} erroneous utcs'
              .format(l=len(x01111), dtns=dtns, e=errcnt))

    fixedx = np.cumsum(x01111)*dtns
    wbad = np.where((x - fixedx)>1e8)[0]
    fixedx[wbad] = np.nan
    debug_(pyfusion.DEBUG, 3, key="repair", msg="repair of W7-X scrambled Langmuir timebase") 
    return(fixedx, msg)

def process_chirps(dd, shots=None, Ns=None, Ms=None, maxd = 4, dfdtunit=1000, minlen=2,
                   plot=0, hold=0, clear_dfdt=False, debug=0, verbose=0):
    """ M and N default to all values present
    >>> 
    """
    # add islead, lead, dfdt entries
    if plot>0 and hold==0: pl.clf()  # plot(hold=0) runs slower here
    num = len(dd['shot'])
    if not dd.has_key('indx'): raise LookupError('dd needs and indx!')
    indx = dd['indx']
    if clear_dfdt or not dd.has_key('dfdt'):
        print('******* initialising dfdt data ************')
        dd.update({'islead': np.zeros(num, dtype=int)})
        dd.update({'lead': np.zeros(num, dtype=int)})
        dd['lead'][:]= -1
        dd.update({'dfdt': np.zeros(num, dtype=np.float32)})
        dd['dfdt'][:]=np.nan
    # extract the shot, M and N

    warn("need to choose a 'nan' equivalent in process_chirps")
    for k in ['shot', 'M', 'N']: 

        w = np.where(dd[k]!= -4)[0]  # right now it is -4 - but effect on colorscale?
        exec("all_{0}=np.array(dd['{0}'])[w]".format(k))

    if Ms == None: 
        Ms = np.unique(all_M) 
        warn('defaulting Ms to {0}'.format(Ms))
    if Ns == None: 
        Ns = np.unique(all_N) 
        warn('defaulting Ns to {0}'.format(Ns))
    debug_(debug,2,key='start_chirps')
    for shot_ in shots:    
        if verbose>0: print("\n{0}:".format(shot_)),
        for N_ in Ns:
            if verbose>0: print(" {0}:".format(N_)),
            for M_ in Ms:
                inds = np.where((shot_ == all_shot) & 
                             (M_ == all_M) & (N_ == all_N))[0]
    # extract, but order in time
                tord=np.argsort(dd['t_mid'][inds])
                if verbose > 4:
                    print("shot {0}, len tord={1}".format(shot_,len(tord)))
                for k in dd.keys(): 
                    if not hasattr(dd[k],'keys'):
                        exec("{0}=np.array(dd['{0}'])[inds[tord]]".format(k))
                (islead, lead, dfdt) = \
                    process_chirps_m_shot(t_mid, freq, indx, maxd=maxd, 
                                          dfdtunit=dfdtunit, minlen=minlen,
                                          plot=plot, debug=debug)
    # write them back into the dict dd
                if len(islead) != 0:            
                    if verbose > 2: print("indx={0}".format(indx))
                    dd['islead'][indx] = islead
                    dd['lead'][indx] = lead
                    dd['dfdt'][indx] = dfdt
                debug_(debug, 3, key='write_back')
    print("minlen = {m}, maxd = {d}, dfdtunit={df}".
          format(m=minlen, d=maxd, df=dfdtunit))
    if plot>0: pl.show()

def plot_shot(DA, sh, ax=None, diags = None, debug=0, fontsize=None):
    """ more flexible - no need to check for errors
    """
    if fontsize != None:     
        pl.rcParams['legend.fontsize']=fontsize

    if diags == None:
        diags = 'i_p,w_p,flat_level,NBI,ech,p_frac'

    inds=np.where(sh==DA.da['shot'])[0]
    pl.rcParams['legend.fontsize']='small'
    if ax == None: ax = pl.gca()
    #(t_mid,w_p,dw_pdt,dw_pdt2,b_0,ech,NBI,p_frac,flat_level)=9*[None]
    t = DA.da['t_mid'][inds]
    b_0 = DA.da['b_0'][inds]
    if (len(np.shape(diags)) == 0): diags = diags.split(',')
    for diag in diags:
        if DA.da.has_key(diag):
            dat = DA.da[diag][inds] ; lab = diag; linestyle='-'
            if diag in 'p_frac': 
                dat=dat*100
                lab+="*100"
            elif diag in 'ech': 
                dat=dat*10
                lab+="*10"
            elif 'flat_level' in diag: 
                dat = 30+200*dat
                linestyle='--'
                
            if diag == 'p_frac': linestyle = ':'    
            ax.plot(t,dat, linestyle=linestyle, label=lab)
    pl.legend()
    debug_(debug,1,key='plot_shot')

    pl.title("{s} {b}T".format(s=sh,b=b_0[0]))

    def do_fetch(self):
        channel_length = int(self.length)
        outdata=np.zeros(1024*2*256+1)
        ##  !! really should put a wrapper around gethjdata to do common stuff          
        #  outfile is only needed if the direct passing of binary won't work
        #  with tempfile.NamedTemporaryFile(prefix="pyfusion_") as outfile:
        ierror, getrets = gethjdata.gethjdata(self.shot,channel_length,self.path,
                                              verbose=VERBOSE, opt=1, ierror=2,
                                              outdata=outdata, outname='')
        if ierror != 0:
            raise LookupError('hj Okada style data not found for {s}:{c}'.format(s=self.shot, c=self.path))

        ch = Channel(self.path,
                     Coords('dummy', (0,0,0)))

        # the intent statement causes the out var to be returned in the result lsit
        # looks like the time,data is interleaved in a 1x256 array
        # it is fed in as real*64, but returns as real*32! (as per fortran decl)
        debug_(pyfusion.DEBUG, 4, key='Heliotron_fetch', msg='after call to getdata')
        # timebase in secs (ms in raw data) - could add a preferred unit?
        # this is partly allowed for in savez_compressed, newload, and
        # for plotting, in the config file.
        # important that the 1e-3 be inside the Timebase()
        output_data = TimeseriesData(timebase=Timebase(1e-3 * getrets[1::2]),
                                 signal=Signal(getrets[2::2]), channels=ch)
        output_data.meta.update({'shot':self.shot})
        if pyfusion.VERBOSE>0: print('HJ config name',self.config_name)
        output_data.config_name = self.config_name         
        stprms = get_static_params(shot=self.shot,signal=self.path)
        if len(list(stprms)) == 0:  # maybe this should be ignored - how do we use it?
            raise LookupError(' failure to get params for {shot}:{path}'
                              .format(shot=self.shot, path=self.path))
        output_data.params = stprms
        return output_data

    def do_fetch(self):
        # Allow for movement of Mirnov signals from A14 to PXI crate
        if pyfusion.VERBOSE > 1:
            print("LHDfetch - timeseries")
        chnl = int(self.channel_number)
        dggn = self.diag_name
        # the clever "-" thing should only be used in members of multi signal diagnostics.
        # so I moved it to base.py
        # dggn = (self.diag_name.split('-'))[-1]  # remove -
        debug_(pyfusion.DEBUG, level=4, key="LHD fetch debug")

        if dggn == "FMD":
            if self.shot < 72380:
                dggn = "SX8O"
                if chnl != 0:
                    chnl = chnl + 33
                    if self.shot < 26208:
                        chnl = chnl + 1

        filename_dict = {"diag_name": dggn, "channel_number": chnl, "shot": self.shot}
        self.basename = path.join(self.filepath, data_fileformat % filename_dict)

        files_exist = path.exists(self.basename + ".dat") and path.exists(self.basename + ".prm")
        if not files_exist:
            if pyfusion.VERBOSE > 3:
                print("RETR: retrieving %d chn %d to %s" % (self.shot, int(chnl), self.filepath))
            tmp = retrieve_to_file(diagg_name=dggn, shot=self.shot, subshot=1, channel=int(chnl), outdir=self.filepath)
            if not path.exists(self.basename + ".dat") and path.exists(self.basename + ".prm"):
                raise Exception, "something is buggered."

        return fetch_data_from_file(self)

def find_clipped(sigs, fact):
    """ look for digitizer or amplifier saturation in all raw 
    signals, and return offending indices.  Digitizer saturation
    is easy, amplifier saturation is softer - need to be careful 

    sigs can be a signal, array or a list of those
    Note that restore_sin() can be used on the sweep, rather than omitting
    the points.
    """
    if isinstance(sigs,(list, tuple)):
        goodbits = np.ones(len(sigs[0]),dtype=np.bool)
        for sig in sigs:
            goodbits = goodbits & find_clipped(sig, fact)
            debug_(debug, 5, key='find_clipped')

        return(goodbits)
    
    cnts, vals = np.histogram(sigs, 50)   # about 100ms for 250k points
    
    if (np.sum(cnts[0]) > fact * np.sum(cnts[1])):
        wunder = (sigs < vals[1])
    else:
        wunder = np.zeros(len(sigs),dtype=np.bool)

    if (np.sum(cnts[-1]) > fact * np.sum(cnts[-2])):
        wover =  (sigs > vals[-2])
    else:
        wover = np.zeros(len(sigs),dtype=np.bool)

    return(~wunder & ~wover)

def get_coords_for_channel(channel_name=None, **kwargs):
    config_dict = kwargs.copy()
    if channel_name:
        config_dict.update(get_config_as_dict("Diagnostic", channel_name))
    coord_name = "dummy"
    coord_values = (0.0, 0.0, 0.0)
    transforms = []
    for k in config_dict.keys():
        if k.startswith("coords_"):
            coord_name = k[7:]
            coord_values = tuple(map(float, config_dict[k].split(",")))
    coords_instance = Coords(coord_name, coord_values)
    if "coord_transform" in config_dict:
        transform_list = pyfusion.config.pf_options("CoordTransform", config_dict["coord_transform"])
        for transform_name in transform_list:
            # this seems to return all the globals too
            transform_class_str = pyfusion.config.pf_get(
                "CoordTransform", config_dict["coord_transform"], transform_name
            )
            # so tyr to exclude the globals
            if pyfusion.config.has_option("global", transform_name):
                continue
            transform_class = import_from_str(transform_class_str)
            # if not hasattr(transform_class, 'output_coord'):
            #    raise Exception('??')
            coords_instance.load_transform(transform_class)
    debug_(pyfusion.DEBUG, 1, key=["coord", "device_name"])
    return coords_instance

    def _get_dphase(self, min_dphase = -np.pi):
        """
        remap to [min_dphase, min_dphase+2pi]
        """
        phases = np.array([self._get_single_channel_phase(i) for i in range(self.signal.shape[0])])
        debug_(pyfusion.DEBUG, 2, key='_get_dphase')
        if self.phase_pairs != None:
            tmp = []
            for a,b in self.phase_pairs:
                ind_a = self.channels.get_channel_index(a)
                ind_b = self.channels.get_channel_index(b)
                tmp.append(phases[ind_b]-phases[ind_a])
            d_phases=remap_periodic(np.array(tmp), min_val=min_dphase)
                
        else:
            d_phases = remap_periodic(phases[1:]-phases[:-1], min_val = min_dphase)
        #d_phase_dataset = OrderedDataSet(ordered_by="channel_1.name")
        d_phase_dataset = BaseOrderedDataSet('d_phase_%s' %datetime.now())
        ## append then sort should be faster than ordereddataset.add() [ fewer sorts()]
        if self.phase_pairs != None:
            for i,phase_pair in enumerate(self.phase_pairs):
                ind_a = self.channels.get_channel_index(phase_pair[0])
                ind_b = self.channels.get_channel_index(phase_pair[1])
                d_phase_dataset.append(FloatDelta(self.channels[ind_a],self.channels[ind_b],d_phases[i]))
        else:
            for i, d_ph in enumerate(d_phases):
                d_phase_dataset.append(FloatDelta(self.channels[i], self.channels[i+1], d_ph))

        #d_phase_dataset.sort()
        return d_phase_dataset

def svd(input_data):
    from timeseries import SVDData
    svddata = SVDData(input_data.timebase, input_data.channels, linalg.svd(input_data.signal, 0))
    svddata.history = input_data.history
    svddata.scales = input_data.scales # need to pass it on to caller
    if pyfusion.DEBUG>4: print("input_data.scales",input_data.scales)
    debug_(pyfusion.DEBUG, key='svd',msg='about to return from svd')
    return svddata

 def reg_item(plot_method):
     for cl_name in class_names:
         if cl_name not in plot_reg:
             plot_reg[cl_name] = [plot_method]
         else:
             plot_reg[cl_name].append(plot_method)
     debug_(pyfusion.DEBUG, 2, key='register plot')
     return plot_method

    def to_sqlalchemy(self,db = 'sqlite:///:memory:',n_recs=1000, chunk=1000):
        """ Write to an sqlachemy database 
            chunk = 500: 2000 34 element recs/sec to (big) sqllite file, 
                    1600/sec to mysql
            'mysql://[email protected]/junk'  (mysql need snans cvtd to nul
        """        
        import sqlalchemy as SA
        def cvt(val):
            if np.isnan(val): return(None)
            elif np.issubdtype(val, np.float):
                return(float(val))
            elif np.issubdtype(val, np.int):
                return(int(val))
            else:
                if self.debug>0: print('unrecognised type {t} in cvt'
                                       .format(t=type(val)))
                return(None)

        # define the table
        self.engine = SA.create_engine(db, echo=self.debug>2)
        self.metadata = SA.MetaData()
        self.fs_table = SA.Table('fs_table', self.metadata)
        (dbkeys,dbtypes)=([],[])
        for k in self.da.keys():
            arr = self.da[k]
            typ = None
            print(k)
            if hasattr(arr,'dtype'):
                if np.issubdtype(arr.dtype, np.int): 
                    typ = SA.Integer
                elif np.issubdtype(arr.dtype, np.float): 
                    typ = SA.Float
                else:
                    print('unknown dtype {d}'.format(d=arr.dtype))
                    
                if typ != None: # if it gets here, it is recognised
                    dbkeys.append(k)
                    dbtypes.append(typ)
                    self.metadata.tables['fs_table'].append_column(SA.Column(k, typ))
                    debug_(self.debug, 1)

            if self.debug>0: print(self.metadata.tables)

        if len(dbkeys)==0: return('nothing to create')
        self.metadata.create_all(self.engine)
        conn=self.engine.connect()
        if self.len > n_recs: print('Warning - only storing n_rec = {n} records'
                                    .format(n=n_recs))
        for c in range(0,min(n_recs,len(self.da[dbkeys[0]])),chunk):
            print(c, min(c+chunk, self.len))
            lst = []
            for i in range(c,min(c+chunk, min(self.len,n_recs))):
                dct = {}
                for (k,key) in enumerate(dbkeys): 
                    dct.update({key: cvt(self.da[key][i])})
                lst.append(dct)    
            if self.debug>0: print(lst)
            conn.execute(self.fs_table.insert(),lst)

def plot_shot(da, sh=None, ax=None, diags = None, marker=None, extra_diags=None, debug=0, quiet=True, fontsize=None, hold=1, **kwargs):
    """ more flexible - no need to check for errors
    Also initial version tailored to LHD, but a fudge used to detect HJ
    extra_diags are simply added to the default list.  This way a default list
    can be used simultaneously with a one or more extra diags.
    """

    if fontsize is not None:     
        pl.rcParams['legend.fontsize']=fontsize

    if diags is None:
        if 'MICRO01' in da.keys():
            diags = 'MICRO01,DIA135'
        else:
            diags = 'i_p,w_p,flat_level,NBI,ech,p_frac'

    if extra_diags is not None:
        diags += "," + extra_diags

    if sh is None: sh = da['shot'][0]

    inds=np.where(sh==da['shot'])[0]
    pl.rcParams['legend.fontsize']='small'
    if ax is None: 
        if hold==0: pl.plot(hold=0)
        ax = pl.gca()
    #(t_mid,w_p,dw_pdt,dw_pdt2,b_0,ech,NBI,p_frac,flat_level)=9*[None]
    t = da['t_mid'][inds]
    b_0 = safe_get(da, 'b_0',inds)

    if (len(np.shape(diags)) == 0): diags = diags.split(',')
    for (i,diag) in enumerate(diags):
        lab = diag   # set default label, marker
        if marker is None: marker = '-' 

        (err, dat) = eval_diag(da=da, inds=inds, diag=diag, debug=debug)
        if quiet and err is not None:
            continue

        if diag in 'p_frac': 
            dat=dat*100
            lab+="*100"
        elif diag in 'ech': 
            dat=dat*10
            lab+="*10"
        elif 'flat_level' in diag: 
            dat = 30+200*dat
            marker='--'

        if diag == 'p_frac': 
            marker = ':'    

        ax.plot(t,dat, marker, label=lab, **kwargs)
        ## no hold keyword in ax.plot #, hold=(hold & (i==0)))
    pl.legend()
    debug_(debug,1,key='plot_shot')

    pl.title("{s} {b:.3f}T".format(s=sh,b=b_0[0]))

    def fetch(self):
        """Fetch each  channel and combine into  a multichannel instance
        of :py:class:`~pyfusion.data.timeseries.TimeseriesData`.

        :rtype: :py:class:`~pyfusion.data.timeseries.TimeseriesData`
        """
 
        ## initially, assume only single channel signals
        # this base debug breakpoint will apply to all flavours of acquisition
        debug_(pyfusion.DEBUG, level=2, key='base_multi_fetch')
        ordered_channel_names = self.ordered_channel_names()
        data_list = []
        channels = ChannelList()
        timebase = None
        meta_dict={}
        if hasattr(self, 't_min') and hasattr(self, 't_max'):
            t_range = [float(self.t_min), float(self.t_max)]
        else:
            t_range = []
        for chan in ordered_channel_names:
            sgn = 1
            if chan[0]=='-': sgn = -sgn
            bare_chan = (chan.split('-'))[-1]
            fetcher_class = import_setting('Diagnostic', bare_chan, 'data_fetcher')
            tmp_data = fetcher_class(self.acq, self.shot,
                                     config_name=bare_chan).fetch()

            if len(t_range) == 2:
                tmp_data = tmp_data.reduce_time(t_range)
            channels.append(tmp_data.channels)
            # two tricky things here - tmp.data.channels only gets one channel hhere
            #config_name for a channel is attached to the multi part -
            #we need to move it to the particular channel 
            # was  channels[-1].config_name = chan
            channels[-1].config_name = tmp_data.config_name
            meta_dict.update(tmp_data.meta)
            #print(tmp_data.signal[-1], sgn)
            tmp_data.signal = sgn * tmp_data.signal
            #print(tmp_data.signal[-1], sgn)
            if timebase == None:
                timebase = tmp_data.timebase
                data_list.append(tmp_data.signal)
            else:
                if hasattr(self, 'skip_timebase_check') and self.skip_timebase_check == 'true':
                    data_list.append(tmp_data.signal)
                else:
                    try:
                        assert_array_almost_equal(timebase, tmp_data.timebase)
                        data_list.append(tmp_data.signal)
                    except:
                        raise
        signal=Signal(data_list)
        output_data = TimeseriesData(signal=signal, timebase=timebase,
                                     channels=channels)
        #output_data.meta.update({'shot':self.shot})
        output_data.meta.update(meta_dict)
        return output_data

 def test_local_access(self):
     # This tests access to a small test tree included in this directory
     shot = -1
     debug_(pyfusion.DEBUG, level=3, key='local_access',msg='enter test_local_access')
     tree_path = os.path.join(TEST_DATA_PATH, 'test_tree')
     from pyfusion.acquisition.utils import get_acq_from_config
     test_acq_class = get_acq_from_config('test_local_tree')
     test_acq = test_acq_class('test_local_tree', test_tree_path=tree_path)
     #test_acq = pyfusion.getAcquisition('test_tree', 
     test_data = test_acq.getdata(shot, 'test_signal')