data Package

data Package

DA_datamining Module

class pyfusion.data.DA_datamining.DA(fileordict, debug=0, verbose=0, load=0, limit=None, mainkey=None)

Class to handle and save data in a special dictionary of arrays referred to hereafter as a “DA”.

• Can deal with databases larger than memory, by using load=0
• Faster to use if load=1, but if you subselect by using extract you get the speed for large data sets (once extract is done).
• Extract can be used over and over to get different data set selections.

Parameters:

• fileordict – An .npz file containing a DA object or a dictionary of arrays sharing a common first dimension, including the result of a loadtxt(dtype=...) command. The filename is processed for env vars ~/ etc, but sometimes this seems to substitute the path of the DA module? (bug)
• load – 1 will immediately load into memory, 0 will defer load allowing some operations (but slowly) without ocnsuming memory.
• mainkey – The main key, not necessarily a unique identifier - e.g it can be shot.
• limit – Decimates the data when loaded into memory (via load=1). It is the most effective space saver, but you need to reload if more (or a different subselection of data) is needed. The alternative is to downselect by using ‘extract=’ (but this applies only to the variables extracted into namespace (e.g. locals())

Returns:

A DA object as described above

Raises:

KeyError, ValueError, LookupError –

Experimental new feature allows use of the DA object itself as a

dictionary (e.g. DA59[‘shot’]).

For more info type help(DA)

Note: This is my prototype of google style python sphinx docstrings - based on

http://www.sphinx-doc.org/en/stable/ext/example_google.html Had to include ‘sphinx.ext.napoleon’ in documentation.conf.py to get the parameters on separate lines.

append(dd)

append the data arrays in dd to the data arrays in self - i.e.

extend the existing arrays. Typical use is in serial processing of a range od shots.

See also append_to_DA_file to add an extra variable

copyda(force=False)

make a deepcopy of self.da typically dd=DAxx.copy() instead of dd=DAxx.da - which will make dd and DAxx.da the same thing (not usually desirable)

extract(dictionary=False, varnames=None, inds=None, limit=None, strict=0, masked=1, debug=0)

extract the listed variables into the dictionary (local by default) selecting those at indices <inds> (all be default variables must be strings, either an array, or separated by commas

if the dictionary is False, return them in a tuple instead Note: returning a list requires you to make the order consistent

if varnames is None - extract all.

e.g. if da is a dictionary or arrays da = DA(‘mydata.npz’) da.extract(‘shot,beta’) plot(shot,beta)

(shot,beta,n_e) = da.extract([‘shot’,’beta’,’n_e’], inds=np.where(da[‘beta’]>3)[0]) # makes a tuple of 3 arrays of data for high beta. Note syntax of where()! It is evaluted in your variable space.

to extract one var, need trailing ”,” (tuple notation) e.g.

(allbeta,) = D54.extract(‘beta’,locals())

which can be abbreviated to

allbeta, = D54.extract(‘beta’,locals())

hist(key, bins=50, nanval=-0.01, percentile=99, label=None)

plot a histogram of Te or resid etc, replacing Nans or infs with nanval, and considering only up to the <percentile>th percentile DA(‘LP20160310_9_L57__amoeba21_1.2_2k.npz’).hist(‘resid’)

Examples

>>> da.hist('resid',percentile=97,label='{k}: {fn}')
>>> da.hist('resid',percentile=97,label='{k}: {actual_fit_params}')
>>> da.hist('resid',percentile=97,label='{k}: {i_diag} {actual_fit_params}')

info(verbose=None)

load(sel=None)

make_attributes()

make each element of the dictionary an attribute of the DA object This is very convenient for operations on more than one dataset e.g. plot(da.im2 - daold.im2) Is this python 3 compatible? It seems to work fine for continuum 3.5.1

plot(key, xkey='t', sharey=1, select=None, sharex='col', masked=1, ebar=1, marker='', elinewidth=0.3, **kwargs)

Plot the member ‘key’ of the Dictionary of Arrays

masked [1] 1 show only ‘unmasked’ points ebar [1] the number of points per errorbar - None will suppress

save(filename, verbose=None, sel=None, use_dictionary=False, tempdir=None, zipopt=-1)

Save as an npz file, using an incremental method, which only uses as much /tmp space as required by each var at a time. Select which to save with sel: if sel is None, save all except for use_dictionary below. If use_dictionary is a valid dictionary, save the values of ANY AND ONLY the LOCAL variables whose names are in the keys for this set. So if you have extracted a subset, and you specify use_dictionary=locals(), only that subset is saved (both in array length, and variables chosen). Beware locals that are not your variables - e.g. mtrand.beta To avoid running out of space on tmp, or to speed up zip - Now included as an argument (Note that the normal os.putenv() doesn’t seem to write to THIS environment use the fudge below - careful - no guarantees) os.environ.__setitem__(‘TMPDIR’,os.getenv(‘HOME’)) actually - this seems OK os.environ[‘IGETFILE’]=’/data/datamining/myView/bin/linux/igetfile’

reload tempfile tempfile.gettempdir() also (‘ZIPOPT’,’”-1”’) (Now incorporated into args, not tested) ** superseded by zlib.Z_DEFAULT_COMPRESSION 0–9 (or -1 for default)

to_sqlalchemy(db='sqlite:///:memory:', mytable='fs_table', n_recs=1000, newfmts={}, chunk=1000)

update(new_dict, check=True)

Add a new variable to the dictionary. Better than simply updating dd, as it allows length check and updates the list of keys.

write_arff(filename, keys=[])

keys is a list of keys to include, and empty list includes all

class pyfusion.data.DA_datamining.Masked_DA(valid_keys=[], baseDA=None, mask=None)

A virtual sub dictionary of a DA, contained in the ‘masked’ attribute

and returning applicable (valid_keys) elements, masked by DA.da[‘mask’] to have Nans in the positions where mask = False or 0

An important side effect is to add the mask array to the main dictionary Probably should NOT be a subclass - we don’t want to do unnecessary copying.

Parameters:

• valid_keys – keys to which mask should be applied.
• mask – An array (usualy 2D) of the same shape as the data, is usually set at a later stage, when the quality or error criteria are evaluated.
• baseDA – not sure why this is required, because this object is usually attached to an existing DA - needs some thought.

Example

>>> from pyfusion.data.DA_datamining import Masked_DA, DA
>>> mydDA=DA('20160310_9_L57',load=1)  # needs to be loaded for this operation
>>> da.masked=Masked_DA(valid_keys=['Te','I0','Vp'], baseDA=myDA)
>>> myDA.da['mask']=-myDA['resid']/abs(myDA['I0'])<.35
>>> clf();plot(myDA.masked['Te']);ylim(0,100)

keys()

Return the keys for elements with validity masking only. To see all elements, including raw data, use keys() on the parent DA. object.

pyfusion.data.DA_datamining.append_to_DA_file(filename, new_dict, force=False)

Adds a new variable to the file - more like ‘update’ than append

Opens filename with mode=a, after checking if the indx variables align force=1 ignores checks for consistent length c.f. the var shot.

Works with a DA file, in contrast to DA.append() which extends a DA

Parameters:

• filename – file to append to
• new_dict – dictionary with new data to append
• force – try to continue if there is a mismatch error

Returns:

no return - side effect is to add a new variable to a DA file

Raises:

ValueError – if new arrays don’t match old.

Example

>>> append_to_DA_file('DAX.npz',dict(N=dd['N'])     # simple
>>> append_to_DA_file('foo.npz',dict((k, mydict[k]) for k in ['N','M']))

Not a member of the class DA, because the class has memory copies of the file, so it would be confusing.

pyfusion.data.DA_datamining.da(filename='300_small.npz', dd=True)

return a da dictionary (used to be called dd - not the DA object) mainly for automated tests of example files.

pyfusion.data.DA_datamining.info_to_bytes(inf)

Eventually this should deal with bytes/unicode compatibility

pyfusion.data.DA_datamining.mylen(ob)

return the length of an array or dictionary for diagnostic info

pyfusion.data.DA_datamining.process_file_name(filename)

Allow shell shortcuts such as ~/ and env var expansion - Note: not tested in windows

pyfusion.data.DA_datamining.report_mem(prev_values=None, msg=None)

LPextra Module

pyfusion.data.LPextra.estimate_error(self, result, debug=1)

pyfusion.data.LPextra.lpfilter(t, i, v, lpf, debug=1, plot=0)

Remove v is used only to determine the number of cycles

amoeba Module

# Amoeba uses the simplex method of Nelder and Mead to maximize a # function of 1 or more variables. # # Copyright (C) 2005 Thomas R. Metcalf # # This program is free software; you can redistribute it and/or modify # it under the terms of the GNU General Public License as published by # the Free Software Foundation; either version 2 of the License, or # (at your option) any later version. # # This program is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # GNU General Public License for more details. # # You should have received a copy of the GNU General Public License # along with this program; if not, write to the Free Software # Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.

pyfusion.data.amoeba.amoeba(var, scale, func, ftolerance=0.0001, xtolerance=0.0001, itmax=500, data=None)

Use the simplex method to maximize a function of 1 or more variables.

Input:

var = the initial guess, a list with one element for each variable scale = the search scale for each variable, a list with one

element for each variable.

func = the function to maximize.

Optional Input:

ftolerance = convergence criterion on the function values (default = 1.e-4) xtolerance = convergence criterion on the variable values (default = 1.e-4) itmax = maximum number of iterations allowed (default = 500). data = data to be passed to func (default = None).

Output:

(varbest,funcvalue,iterations) varbest = a list of the variables at the maximum. funcvalue = the function value at the maximum. iterations = the number of iterations used.

• Setting itmax to zero disables the itmax check and the routine will run until convergence, even if it takes forever.
• Setting ftolerance or xtolerance to 0.0 turns that convergence criterion off. But do not set both ftolerance and xtolerance to zero or the routine will exit immediately without finding the maximum.
• To check for convergence, check if (iterations < itmax).

The function should be defined like func(var,data) where data is optional data to pass to the function.

Example:

>>> import amoeba
>>> def afunc(var,data=None): return 1.0-var[0]*var[0]-var[1]*var[1]
>>> print amoeba.amoeba([0.25,0.25],[0.5,0.5],afunc)
([0.0, 0.0], 1.0, 17)

Version 1.0 2005-March-28 T. Metcalf

1.1 2005-March-29 T. Metcalf - Use scale in simsize calculation.

• Use func convergence and x convergence rather than func convergence or x convergence.

1.2 2005-April-03 T. Metcalf - When contracting, contract the whole

simplex.

bdb added doctest

base Module

Base classes for data.

class pyfusion.data.base.BaseCoordTransform

Bases: object

Base class does nothing useful at the moment

input_coords = 'base_input'

output_coords = 'base_output'

transform(coords)

class pyfusion.data.base.BaseData(*args, **kwargs)

Bases: object

Base class for handling processed data.

In general, specialised subclasses of BaseData will be used to handle processed data rather than BaseData itself.

Usage: ..........

save()

class pyfusion.data.base.BaseDataSet(label='')

Bases: object

add(item)

copy()

pop()

remove(item)

save()

update(item)

class pyfusion.data.base.BaseOrderedDataSet(label='')

Bases: object

append(item)

save()

class pyfusion.data.base.Channel(name, coords, source='', parent_device='not specified')

Bases: object

save()

applicable only to ORM db

class pyfusion.data.base.ChannelList(*args)

Bases: list

get_channel_index(channel_name, bounds=None)

repopulate()

save()

class pyfusion.data.base.Coords(default_coords_name, default_coords_tuple, **kwargs)

Bases: object

Stores coordinates with an interface for coordinate transforms.

add_coords(**kwargs)

load_from_config(**kwargs)

load_transform(transform_class)

save()

class pyfusion.data.base.DataSet(label='')

Bases: pyfusion.data.base.BaseDataSet

copy(input_data, *args, **kwargs)

normalise(input_data, *args, **kwargs)

reduce_time(input_data, *args, **kwargs)

remove_noncontiguous(input_data, *args, **kwargs)

segment(input_data, *args, **kwargs)

subtract_mean(input_data, *args, **kwargs)

class pyfusion.data.base.DynamicDataSet(label='')

Bases: pyfusion.data.base.BaseDataSet

class pyfusion.data.base.FloatDelta(channel_1, channel_2, delta, **kwargs)

Bases: pyfusion.data.base.BaseData

class pyfusion.data.base.MetaMethods

Bases: type

Metaclass which provides filter and plot methods for data classes.

class pyfusion.data.base.OrderedDataSet(label='')

Bases: pyfusion.data.base.BaseOrderedDataSet

class pyfusion.data.base.OrderedDataSetItem(item, index)

Bases: object

class pyfusion.data.base.PfMetaData

Bases: dict

pyfusion.data.base.get_coords_for_channel(channel_name=None, **kwargs)

pyfusion.data.base.get_history_args_string(*args, **kwargs)

pyfusion.data.base.history_reg_method(method)

Wrapper for filter and plot methods which updates the data history.

pyfusion.data.base.orm_load_basedata(man)

pyfusion.data.base.orm_load_basedataset(man)

pyfusion.data.base.orm_load_baseordereddataset(man)

pyfusion.data.base.orm_load_channel(man)

pyfusion.data.base.orm_load_channel_map(man)

pyfusion.data.base.orm_load_channellist(man)

pyfusion.data.base.orm_load_dataset(man)

pyfusion.data.base.orm_load_dynamic_dataset(man)

pyfusion.data.base.orm_load_floatdelta(man)

pyfusion.data.base.setup_coords(man)

check_crosstalk Module

convenience Module

pyfusion.data.convenience.between(var, lower, upper=None, closed=True)

return whether var is between lower and upper includes end points if closed=True alternative call is between(var, range) e.g. between(x, [1, 2])

pyfusion.data.convenience.broaden(inds, data=None, dw=1)

broaden a set of indices or data in width by dw each side

pyfusion.data.convenience.btw(var, lower, upper=None, closed=True)

return whether var is between lower and upper includes end points if closed=True alternative call is between(var, range) e.g. between(x, [1, 2])

pyfusion.data.convenience.bw(var, lower, upper=None, closed=True)

return whether var is between lower and upper includes end points if closed=True alternative call is between(var, range) e.g. between(x, [1, 2])

pyfusion.data.convenience.decimate(data, limit=None, fraction=None)

reduce the number of items to a limit or by a fraction returns the same data every call decimation is regular, not random

pyfusion.data.convenience.his(xa, tabs=False, sort=-1, dfmt=None, total=True, maxsort=None)

print the counts and fraction of xa binned as integers sort=1,-1 sorts by most frequent (first, last), maxsort sets the maximum number kept after the sort

pyfusion.data.convenience.inds_from_list(var, lst)

return a list of indices of var which match the items in lst. This is to replace shot in [90000,90001] etc which can’t be used in where() (gets msg ”...... Use a.any() or a.all()”)

>>> inds_from_list([1,4,5,1], [1,2,3])
array([0, 3])

pyfusion.data.convenience.inlist(var, lst)

return a list of True/False of var which match the items in lst. This is to replace shot in [90000,90001] etc which would get msg ”...... Use a.any() or a.all()”

>>> inlist([1,4,5,1], [1,2,3])
array([ True, False, False,  True], dtype=bool)

pyfusion.data.convenience.whr(*args)

short for where, doesn’t need np. or [0] at end, and echoes number

evalexpr_script Module

filters Module

Some un-pythonic code here (checking instance type inside function). Need to figure out a better way to do this.

python3 issues:

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/filters.py:65: DeprecationWarning: classic int division

nice = [2**p * n/16 for p in range(minp2,maxp2) for n in [16, 18, 20, 24, 27]]

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/utils.py:120: DeprecationWarning: classic int division

ipks = find_peaks(np.abs(FT)[0:ns/2], minratio = minratio, debug=1)

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/filters.py:443: DeprecationWarning: classic int division

twid = 2*(1+max(n_pb_low - n_sb_low,n_sb_hi - n_pb_hi)/2)

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/base.py:57: UserWarning: defaulting taper to 1 as band edges are sharp /home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/filters.py:508: DeprecationWarning: classic int division

if np.mod(NA,2)==0: mask[:NA/2:-1] = mask[1:(NA/2)] # even

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/filters.py:490: DeprecationWarning: classic int division

low_mid = n_pb_low - twid/2

/home/bdb112/pyfusion/mon121210/pyfusion/pyfusion/data/filters.py:491: DeprecationWarning: classic int division

high_mid = n_pb_hi + twid/2

pyfusion.data.filters.copy(input_data)

safe (deep) copy of camplete data

pyfusion.data.filters.correlate(input_data, index_1, index_2, **kwargs)

pyfusion.data.filters.cps(a, b)

pyfusion.data.filters.downsample(input_data, skip=10, chan=None, copy=False)

Good example of filter that changes the size of the data.

class pyfusion.data.filters.dummysig(tb, sig)

timebase is given as an array in seconds 30 (10 lots of three) signals are generated: as scaled copies 1,2,3x

pyfusion.data.filters.filter_fourier_bandpass(input_data, passband, stopband, taper=None, debug=None)

Note: Is MUCH (2.2x faster) more efficient to use real ffts, (implemented April) Use a Fourier space taper/tophat or pseudo gaussian filter to perform narrowband filtering (much narrower than butterworth). Problem is that bursts may generate ringing. This should be better with taper=2, but it is not clear

See the __main__ code below for nice test facilities twid is the width of the transition from stop to pass (not impl.?) >>> tb = timebase(np.linspace(0,20,512)) >>> w = 2*np.pi* 1 # 1 Hertz >>> dat = dummysig(tb,np.sin(w*tb)*(tb<np.max(tb)/3)) >>> fop = filter_fourier_bandpass(dat,[0.9,1.1],[0.8,1.2],debug=1).signal[0]

Testing can be done on the dummy data set generated after running filters.py e.g. (with pyfusion,DEBUG=2 make_mask(512, [0.8,.93], [0.9,.98],dat,2) # medium sharp shoulder fopmed = filter_fourier_bandpass(dat,[9.5,10.5],[9,11],debug=1).signal[0] # very sharp shoulders fopsharp = filter_fourier_bandpass(dat,[9.1,10.9],[9,11],debug=1)

pyfusion.data.filters.flucstruc(input_data, min_dphase=-3.141592653589793, group=<function fs_group_geometric>, method='rms', separate=True, label=None, segment=0, segment_overlap=1.0)

If segment is 0, then we dont segment the data (assume already done)

pyfusion.data.filters.fs_group_geometric(input_data, max_energy=1.0)

no filtering implemented yet we don’t register this as a filter, because it doesn’t return a Data or DataSet subclass # TODO: write docs for how to use max_energy - not obvious if using flucstruc() filter...

pyfusion.data.filters.fs_group_threshold(input_data, threshold=0.7)

no filtering implemented yet we don’t register this as a filter, because it doesn’t return a Data or DataSet subclass

pyfusion.data.filters.get_optimum_time_range(input_data, new_time_range, try_more=0.0002)

This grabs a few more (or a few less, if enough not available in input_data) points than requested so that the FFT is more efficient. try_more is the number of samples to increase/decrease as a fraction of the input data.

Note: For FFTW, it is more efficient to zero pad to a nice number above, even if it is a long way away. This is always true for Fourier filtering, in which case you never see the zeros. For general applications, the zeros might be confusing if you forget they have been put there.

pyfusion.data.filters.integrate(input_data, baseline=[], delta_t=0.01, chan=None, copy=False)

Return the time integral of a signal, with baseline optionally removed before integration

baseline = None No removal

= [] sloping baseline for delta_t at either end of the data = [t0, t1, t2, t3] as above, but t values are explicit = [t0, t1] as above, but simple constant removal.

Perhaps the first sample will be a Nan - no,the samples are all displaced by one If we used the trapeziodal method, the first and last samples are ‘incorrect’ and Nans may be appropriate

pyfusion.data.filters.make_mask(NA, norm_passband, norm_stopband, input_data, taper)

works well now, except that the stopband is adjusted to be symmetric about the passband (take the average of the differences The problem with crashes (zero mask) was solved by shifting the mask before and after integrating, also a test for aliasing (on the mask before integration).

pyfusion.data.filters.next_nice_number(N)

return the next highest power of 2 including nice fractions (e.g. 2**n *5/4) takes about 10us - should rewrite more carefully to calculate starting from smallest power of 2 less than N, but hard to do better >>> print(next_nice_number(256), next_nice_number(257)) (256, 288)

Have to be careful this doen’t take more time than it saves!

pyfusion.data.filters.normalise(input_data, method=None, separate=False)

method=None -> default, method=0 -> DON’T normalise

pyfusion.data.filters.reduce_time(input_data, new_time_range, fftopt=0)

reduce the time range of the input data in place(copy=False) or the returned Dataset (copy=True - default at present). if fftopt, then extend time if possible, or if not reduce it so that ffts run reasonably fast. Should consider moving this to actual filters? But this way users can obtain optimum fft even without filters. The fftopt is only visited when it is a dataset, and this isn’t happening

pyfusion.data.filters.register(*class_names)

pyfusion.data.filters.remove_baseline(input_data, baseline=None, delta_t=0.01, chan=None, copy=False)

Remove a tilted baseline from a signal if the baseline is 4 elements, correct at two points (mid point of those intervals) baseline in the same units as the timebase

pyfusion.data.filters.remove_noncontiguous(input_dataset)

pyfusion.data.filters.segment(input_data, n_samples, overlap=1.0)

Break into segments length n_samples.

Overlap of 2.0 starts a new segment halfway into previous, overlap=1 is no overlap. overlap should divide into n_samples. Probably should consider a nicer definition such as in pyfusion 0

n_samples < 1 implies a time interval, which is adjusted to suit fft otherwise n_samples is taken literally and not adjusted. fractional n_samples>1 allows the step size to be fine-tuned.

pyfusion.data.filters.sp_filter_butterworth_bandpass(input_data, passband, stopband, max_passband_loss, min_stopband_attenuation, btype='bandpass')

** Warning - fails for a single signal in the enumerate step.

This actually does ALL butterworth filters - just select bptype and use scalars instead of [x,y] for the passband.

e.g df=data.sp_filter_butterworth_bandpass(2e3,4e3,2,20,btype=’lowpass’)

pyfusion.data.filters.subtract_mean(input_data)

pyfusion.data.filters.svd(input_data)

fits Module

pyfusion.data.fits.fitbiexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits.fitgauss(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits.fitmirexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits.gaus(x)

mis-spell to avoid clash

pyfusion.data.fits.mybiexp(x, params)

pyfusion.data.fits.mygauss(x, params)

pyfusion.data.fits.mymirexp(x, params)

pyfusion.data.fits.residuals(params, fn, x, y, yerrs=None)

fits0 Module

pyfusion.data.fits0.fitbiexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits0.fitgauss(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits0.fitmirexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits0.gaus(x)

mis-spell to avoid clash

pyfusion.data.fits0.mybiexp(x, params)

pyfusion.data.fits0.mygauss(x, params)

pyfusion.data.fits0.mymirexp(x, params)

pyfusion.data.fits0.residuals(params, x, y, fn)

fits1 Module

pyfusion.data.fits1.fitbiexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits1.fitgauss(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits1.fitmirexp(x_data, y_data, yerrs=None, ax=None, guess=None, pltkwargs={})

pyfusion.data.fits1.gaus(x)

mis-spell to avoid clash

pyfusion.data.fits1.mybiexp(x, params)

pyfusion.data.fits1.mygauss(x, params)

pyfusion.data.fits1.mymirexp(x, params)

pyfusion.data.fits1.residuals(params, fn, x, y, yerrs=None)

histogramHD Module

class pyfusion.data.histogramHD.CoordHD(dims, debug=0)

get(*indices)

get(1,2,3) args not a tuple

set(indices, val)

set((1,2,3),49.)

class pyfusion.data.histogramHD.CoordHDs(dims, debug=0)

get(indices)

string version doesn’t require the tuple notation (*)

set(indices, val)

set(inds,49.)

pyfusion.data.histogramHD.find_eps(x, value=None)

find the smallest number that will always exceed the representational accuracy of x in the range around value.

pyfusion.data.histogramHD.histogramHD(d, bins=None, method='safe')

make a histogram of data too high in dimension to use histogramdd, by successive calls to histogramdd(). The simplest implementation is to assume all bins are equal and in the range -pi..pi subdivided into nbins bins.

First version will use the coo_utils for simplicity - not sure how efficient lookup is though - and there is no obvious way to make an array larger than memory

Second - use a dictionary to store the results - 120ns lookup in a 100k dict Assume we can do a rank R array with histogramdd Should really store the bin boundaries for later reference

plots Module

Note, plots (this file) doesn’t have unittests

Problems with checkbuttons working since 2011-2012? Temporary fix is to block in svdplot when hold==1. This enables the function, then you kill the window to move on (see plot_svd.py) Also attempted to use subplots here to tidy up putting additonal graphs on top. See plots_1.py and svd_plots1.py - need to sovle subplot pars problem bdb made format more conforming. 2013

class pyfusion.data.plots.Energy(energy_list, initial_list)

add(elmt)

sub(elmt)

pyfusion.data.plots.findZero(i, x, y1, y2)

pyfusion.data.plots.fsplot_phase(input_data, closed=True, ax=None, hold=0, offset=0, block=False)

plot the phase of a flucstruc, optionally inserting the first point at the end (if closed=True). Applies to closed arrays (e.g complete 2pi). Until Feb 2013, this version did not yet attempt to take into account angles, or check that adjacent channels are adjacent (i.e. ch2-ch1, ch2-c2 etc). Channel names are taken from the fs and plotted abbreviated

1/1/2011: TODO This appears to work only for database = None config 1/17/2011: bdb: May be fixed - I had used channel instead of channel.name

pyfusion.data.plots.join_ends(inarray, add_2pi=False, add_360deg=False, add_lenarray=False, add_one=False)

used in old code, needs clean up....

pyfusion.data.plots.mydiff(t, y)

numpy diff in a wrapper to give correct time and units

pyfusion.data.plots.myiden(y)

pyfusion.data.plots.myiden2(t, y)

pyfusion.data.plots.overlap(interval1, interval2)

return true if there is any overlap between the two intervals probably belongs in convenience or utils

pyfusion.data.plots.plot_fs_groups(input_data, fs_grouper)

Show what would be generated by the supplied fs_grouper function.

An example of a grouper funcion is pyfusion.data.filters.fs_group_threshold

The grouper function takes an SVDData instance as its first input and returns an iterable where each element is a set of singular value indices which define a flucstruc.

pyfusion.data.plots.plot_signals(input_data, filename=None, downsamplefactor=1, n_columns=1, hspace=None, sharey=False, sharex=True, ylim=None, xlim=None, marker='None', decimate=0, markersize=2, linestyle=True, labelfmt='{short_name} {units}', filldown=True, suptitle='shot {shot}', raw_names=False, labeleg='False', color='b', fun=<function myiden>, fun2=None, **kwargs)

Plot a figure full of signals using n_columns[1],

sharey [=1] “gangs” y axes - sim for sharex - sharex=None stops this sharey: 2 gangs all but first (top) axis x axes are ganged by default: see Note:

fun, fun2: optionally plot a function of the signal. fun= refers to a function of one variable. fun2 is a function (t,x), such as one that returns a different timebase (diff should do this) if fun2 is given (a function of t and sig), then fun is ignored

labelfmt[“{short_name} {units}”] controls the channel labels.

The default ignores the shot and uses an abbreviated form of the channel name. If the short form is very short, it becomes a y label. A full version is “{name} {units}” and if > 8 chars, will become the x label. Even longer is “Shot={shot}, k_h={kh}, {name}”

labeleg: If ‘True’ put label in legend - else use this str as a legend lab linestyle: default of True means use ‘-‘ if no marker, and nothing if a

marker is given e.g. marker of ‘.’ and markersize<1 will produce “shaded” waveforms, good to see harmonic structure even without zooming in (need to adjust markersize or plot size for best results.

raw_names uses “digitiser” names, otherwise use names from the config file Note = sharex that to allow implicit overlay by using the same subplot specs the sharex must be the same between main and overlay - hence the use of explicit sharex = None

suptitle by default refers to the shot number

pyfusion.data.plots.plot_spectrogram(input_data, windowfn=None, units='kHz', channel_number=0, filename=None, coloraxis=None, noverlap=0, NFFT=None, title=None, **kwargs)

title will be auto generated, if supplied, include ‘+’ to include the auto-generated part

pyfusion.data.plots.posNegFill(x, y1, y2)

pyfusion.data.plots.register(*class_names)

pyfusion.data.plots.set_axis_if_OK(ax, xlims, ylims)

set axes if they containg at least a bit of the plot originally meant to be used when forcing specific axes, to allow for attempt to force axes suitable for LHD on HJ data for example

pyfusion.data.plots.svdplot(input_data, fmax=None, hold=0)

plots_1 Module

Note, plots (this file) doesn’t have unittests

class pyfusion.data.plots_1.Energy(energy_list, initial_list)

add(elmt)

sub(elmt)

pyfusion.data.plots_1.findZero(i, x, y1, y2)

pyfusion.data.plots_1.fsplot_phase(input_data, closed=True, ax=None, hold=0)

plot the phase of a flucstruc, optionally calculating the extra difference e.g. MP6-MP1 NOT just replicating the last point at the beginning (if closed=True). This version does not yet attempt to take into account angles, or check that adjacent channels are adjacent (i.e. ch2-ch1, ch2-c2 etc). Channel names are taken from the fs and plotted abbreviated

1/1/2011: TODO This appears to work only for database=None config 1/17/2011: bdb: May be fixed - I had used channel instead of channel.name

pyfusion.data.plots_1.join_ends(inarray, add_2pi=False, add_360deg=False, add_lenarray=False, add_one=False)

used in old code, needs clean up....

pyfusion.data.plots_1.plot_fs_groups(input_data, fs_grouper)

Show what would be generated by the supplied fs_grouper function.

An example of a grouper funcion is pyfusion.data.filters.fs_group_threshold

The grouper function takes an SVDData instance as its first input and returns an iterable where each element is a set of singular value indices which define a flucstruc.

pyfusion.data.plots_1.plot_signals(input_data, filename=None, downsamplefactor=1, n_columns=1, hspace=None, sharey=False, sharex=True, ylim=None, xlim=None, marker='None', markersize=0.3, linestyle=True, labelfmt='%(short_name)s', filldown=True)

Plot a figure full of signals using n_columns[1],

sharey [=1] “gangs” y axes x axes are always ganged labelfmt[“%(short_name)s”] controls the channel labels.

The default ignores the shot and uses an abbreviated form of the channel name. If the short form is very short, it becomes a y label. A full version is “%(name)s” and if > 8 chars, will become the x label. Even longer is “Shot=%(shot), k_h=%(kh)s, %(name)s”

linestyle: default of True means use ‘-‘ if no marker, and nothing if a

marker is given e.g. marker of ‘.’ and markersize<1 will produce “shaded” waveforms, good to see harmonic structure even without zooming in (need to adjust markersize or plot size for best results.

pyfusion.data.plots_1.plot_spectrogram(input_data, windowfn=None, units='kHz', channel_number=0, filename=None, coloraxis=None, noverlap=0, NFFT=None, **kwargs)

pyfusion.data.plots_1.posNegFill(x, y1, y2)

pyfusion.data.plots_1.register(*class_names)

pyfusion.data.plots_1.svdplot(input_data, fmax=None, axs=None, hold=0)

axs is a set of four axes

process_swept_Langmuir Module

Complete processing of swept Langmuir data, using two classes. New version allows choice of algorithm and adds more params This supersedes the example in examples/process_swept_Langmuir.py

pyfusion v0.7.0 - estimate error, include a lpfilter, resid normed to I0 This version assumes that all required sweep data are in a multi-channel diagnostic (dataset)

• the other version (should not be in the repo) grabs sweepV data as required

See main README.rst for more recent changes

• see also the docs of the main method process_swept_Langmuir()
• The convenience script examples/run_process_Langmuir runs a list of shots for both limiter segments

with suitable inputs

Example of step by step operation and tuning:

>>> run pyfusion/data/process_swept_Langmuir
>>> LP.process_swept_Langmuir()
>>> LP.write_DA('20160310_9_L57')
>>> # then to tune mask:
>>> from pyfusion.data.DA_datamining import Masked_DA, DA
>>> myDA=DA('20160310_9_L57',load=1)  # just to be sure it is a good DA
>>> myDA.masked=Masked_DA(['Te','I0','Vf'], baseDA=myDA)
>>> myDA.da['mask']=(myDA['resid']/abs(myDA['I0'])<.35) & (myDA['nits']<100)
>>> clf();plot(myDA.masked['Te']);ylim(0,100)

pyfusion.data.process_swept_Langmuir.AC(x)

class pyfusion.data.process_swept_Langmuir.LPfitter(i, v, fit_params=None, plot=1, verbose=1, parent=None, debug=1)

Mainly used by process_swept_Langmuir() but also stands alone. Assumes the i and v are in arrays already. See :py:meth:LPfitter.__init__ for more info

This is the code presently used to default fit params

>>> self.fit_params = dict(maxits=300, alg='leastsq', ftol=1e-4, xtol=1e-4, Lnorm=1.2, track_ratio=1.2, lpf=None)

Safer to type this to see the actual defaults if the code has changed:

>>> from pyfusion.data.process_swept_Langmuir import LPfitter
>>> LPfitter(None, None).fit_params

LPchar(v, Te, Vf, I0)

error_fun(vars, data)

return RMS error

fit(init=None, plot=None, fit_params=None, default_init={'I0': None, 'Te': 50, 'Vf': 15})

Perform a curve fit operation with the given parameter and initial value. ‘None’ invokes the fit_params determined at creation of the class, and for ‘init’, the default_init

The final parameters actually used are saved in actual_params

plotchar(v, Te, Vf, I0, alpha=1, col='g', linewidth=2)

plot a Langmuir characteristic

residuals(params, x, y)

class pyfusion.data.process_swept_Langmuir.Langmuir_data(shot, i_diag, v_diag, dev_name='W7X', debug=0, plot=1, verbose=0, params=None)

get the fits for a multi (or single) cpt diagnostic using segments If a dictionary params is supplied, use these to process the data, otherwise just read in the data.

Eternal question: Should I pass parameters or set them in the object? Answer for now - if you want them in the record of analysis, pass them through process_swept_Langmuir

obvious candidates for object: debug, plot, select obvious candidates for args: Lnorm, initial_TeVfI0, clipfact, rest_swp

If params contains a filename entry, save the result as that name

Final solution: save the ACTUAL params used in the object, and then save THOSE in when writing, not the onbes entered.

See process_swept_Langmuir() for the processing arguments

Create a Langmuir Data oject for later processing

Parameters:

• shot – shot number
• i_diag – name of probe corrent diagnostic - at present should be a multi channel
• v_diag – sweep voltage diagnostic - should include all voltages referred to by all the corrent channels config as sweepV
• dev_name –
• debug –
• plot – can be set later in process_langmuir
• params – if set, can cause the process step to be executed after loading

fit_swept_Langmuir_seg_chan(m_sig, i_sig, v_sig, i_segds, channame, fit_params=None, clipfact=5, initial_TeVfI0=None, plot=None)

Lowish level routine - takes care of i clipping and but not leakage or restoring sweep.

Returns fitted params including error

fit_swept_Langmuir_seg_multi(m_seg, i_seg, v_seg, clipfact=5, initial_TeVfI0=None, fit_params=None, plot=None)

get_iprobe(leakage=None, t_comp=None)

The main purpose is to subtract leakage currents Will use the full data, as the t_range is meant to be plasma interval returns a copy of the measured courremt, overwritten with the corrected iprobe

prepare_sweeps(rest_swp='auto', sweep_freq=500, Vpp=180, clip_level_minus=-88)

extracts sweep voltage data for all probes, call restore_sin if necessary Initially, the voltage data is assumed to be in a dictionary

process_constant_Langmuir(t_range=None, t_comp=[0, 0.1], leakage=None, threshold=0.01, threshchan=12, filename=None, amu=1, plot=None, return_data=False, suffix='')

This function should analyse the shots with constant Voltage, usually measuring I_sat. It will try to use as much of process_swept_Langmuir as possible and the pyfusion environment.

process_swept_Langmuir(t_range=None, t_comp=[0, 0.1], fit_params={'lpf': None, 'esterr': 1, 'maxits': 200, 'alg': 'leastsq'}, initial_TeVfI0={'I0': None, 'Te': 50, 'Vf': 15}, dtseg=0.004, overlap=1, rest_swp='auto', clipfact=5, clip_iprobe=None, leakage=None, threshold=0.01, threshchan=12, filename=None, amu=1, plot=None, return_data=False, sweep_freq=500, suffix='')

Process the I, V probe data in the Langmuir_data object

Parameters:

• fit_params – a dictionary of parameters used in the fitting stage (passed to the LP_fitter class) Most of these have reasonable defaults (for W7X limiter probes), and the actual values used are attached to the Langmuir_data object and saved in the the DA_ file.
• filename – The output (DA) file name, defaults to a pattern containing the shot number etc. If a user entered value contains a ‘*’, then that pattern is inserted in its place e.g. filename=’test_*’ generates a DA file test_LP20160308_23_L53.npz.
• t_range – time range in secs over which data is processed - None uses a minimum current value as follows.
• threshold – processing starts once this value is exceeded in the following channel (assuming t_range=None).
• threshchan – the current channel number (base 0) used to detect plasma start.
• overlap – degree of overlap of the analysed segments: 1 is no overlap, 2 means half the segment is common to neighboring data.
• dtseg – The length of a segment over shich analysis is performed, either in seconds (float) or samples (int)
• rest_swp – Restore clipped sweep V if True. If None, for March 10 only.
• leakage – A complex conductance to represent the crosstalk between voltage and current channels. If None it is automatically calculated for the interval t_comp. To set to zero, use [0,0]

Returns:

A DA object as described above

Raises:

ValueError –

==> results[time,probe,quantity] plot = 1 : V-I and data if there are not too many. plot >= 2 : V-I curves plot >= 3 : ditto + time plot plot >= 4 : ditto + all I-V iterations

Can send parameters through the init or process - either way they are all

recorded in actual_params

clip_iprobe = [-0.015, .02] # used to check if a resistive term is affecting Te

Start by processing in the ideal order: - logical, but processes more data than necessary fix up the voltage sweep compensate I detect plasma reduce time segment and process

Faster order - reduces time range earlier: (maybe implement later) detect plasma time range with quick and dirty method restrict time range, but keep pre-shot data accessible for evaluation of success of removal

write_DA(filename)

pyfusion.data.process_swept_Langmuir.find_clipped(sigs, clipfact)

look for digitizer or amplifier saturation in all raw signals, and return offending indices. Detecting 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.

pyfusion.data.process_swept_Langmuir.fixup(da=None, locs=None, newmask=None, suppress_ne=None, newpath='/tmp')

some miscellaneous post processing mainly to set the time offset of ECH need to run plot_LP2d.py first. sets t_zero, and tweaks the Nan coverage

Parameters:

• da – a dictionary of arrays (DA) file.
• locs – a dictionary containing the variable dtprobe, usually locals()

pyfusion.data.process_swept_Langmuir.tog(colls=None, hide_all=0)

toggle lines on and off, unless hide_all=1

pyfusion.data.process_swept_Langmuir.tryone(var, data, mrk='r', w=4)

evaluate fit and show curve - meant for manual use tryone([Te,Vf,I0],dict(v=sweepV,i=iprobe))

process_swept_Langmuir_0 Module

process_swept_Langmuir_1 Module

This version grabs sweepV data as required - the other version assumes that all required sweep data are in a multi-channel diagnostic (dataset)

pyfusion.data.process_swept_Langmuir_1.AC(x)

pyfusion.data.process_swept_Langmuir_1.LPfit(v, i, init={'Vp': 15, 'I0': None, 'Te': 50}, plot=None, maxits=None)

class pyfusion.data.process_swept_Langmuir_1.LPfitter(i, v, pnorm=1, maxits=100)

LPchar(v, Te, Vp, I0)

error_fun(vars, data)

return RMS error

plotchar(v, Te, Vp, I0, alpha=1, col='g', linewidth=2)

class pyfusion.data.process_swept_Langmuir_1.Langmuir_data(shot, i_diag, dev_name='W7X', debug=1, plot=1)

get the fits for a multi (or single) cpt diagnostic using segments

fit_one_swept_Langmuir(m_seg, i_seg, v_seg, clipfact=5)

Lowish level routine - takes care of i clipping and but not leakage or restoring sweep.

Returns fitted params including error

get_iprobe(leakage=None, t_comp=[0, 0.1])

main purpose is to subtract leakage currents Will use the full data, as the t_range is meant to be plasma interval returns a copy of the measured courremt, overwritten with the corrected iprobe

prepare_sweeps(rest_swp=None)

extracts sweep voltage data for all probes, call restore_sin if necessary Initially, the voltage data is assumed to be in a dictionary

process_swept_Langmuir(t_range=None, t_comp=[0, 0.1], dtseg=0.001, overlap=1, rest_swp=None, clipfact=5, leakage=None)

pyfusion.data.process_swept_Langmuir_1.find_clipped(sigs, clipfact)

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.

pyfusion.data.process_swept_Langmuir_1.tog(colls=None, hide_all=0)

toggle lines on and off, unless hide_all=1

pyfusion.data.process_swept_Langmuir_1.tryone(var, data, mrk='r', w=4)

evaluate fit and show curve - meant for manual use tryone([Te,Vp,I0],dict(v=sweepV,i=iprobe))

pyfusion_sigproc Module

signal processing peculiar to plasma devices - more general stuff in signal_processing.py Simple criterion is if it imports pyfusion it should be here, if general, in signal_processing. Important to separate from routines that initialise database, as they can’t be recompiled/reloaded easily during debugging.

pyfusion.data.pyfusion_sigproc.find_shot_times(dev, shot, activity_indicator=None, debug=0)

Note: This is inside a try/except - errors will just skip over!! fixme From the channel specified in the expression “activity_indicator”, determine the beginning and end of pulse. A suitable expression is hard to find. For example, density usually persists too long after the shot, and sxrays appear a little late in the shot. The magnetics may be useful if magnet power supply noise could be removed. (had trouble with lhd 50628 until adj threshold ?start and end were at 0.5-0.6 secs ) >>> import pyfusion >>> sh=pyfusion.core.get_shot(15043,activity_indicator=”MP4”) >>> print(‘start=%.3g, end=%.3g’ % (sh.pulse_start, sh.pulse_end) ) start=177, end=218 >>> sh=pyfusion.core.get_shot(33372,activity_indicator=”MP4”) >>> print(‘start=%.3g, end=%.3g’ % (sh.pulse_start, sh.pulse_end) ) start=168, end=290

restore_sin Module

subroutine to restore the clipped part of the sinusoid - see examples/restore_sin for the old script version

pyfusion.data.restore_sin.analytic_signal(x)

A short-cut assuming that the incoming signal is reasonable e.g. fairly pure sinusoid. So far this has no strategy to minimize fourier transform time.

pyfusion.data.restore_sin.restore_sin(data, t_range=None, chan=None, method=2, sweep_freq=500, Vpp=180, clip_level_minus=-88, verbose=1)

Restore the clipped part of the sinusoid - see examples/restore_sin for the old script version

sweep_freq: sinusoidal sweep freq in Hz Vpp: pp value of signal before it was clipped clip_level_minus: value to ensure even soft clipping is excluded method: so far only 2

save_compress Module

ported from the old pyfusion:

A “smart compression” replacement for savez, assuming data is quantised. The quantum is found, and the data replaced by a product of and integer sequence and quantum with offset. delta encoding is optional and often saves space. The efficiency is better than bz2 of ascii data for individual channels, and a little worse if many channels are lumped together with a common timebase in the bz2 ascii format, because save_compress stores individual timebases. $ wc -c /f/python/local_data/027/27999_P* 2300176 total

At the moment (2010), save_compress is not explicitly implemented - it is effected by calling discretise_signal() with a filename argument.

July 2009 - long-standing error in delta_encode_signal fixed (had not been usable before)

may 2016 - version 104 works around an error caused by W7X corrupted timebases (all 0s and nans)

pyfusion.data.save_compress.discretise_array(arrin, eps=0, bits=0, maxcount=0, verbose=None, delta_encode=False)

Return an integer array and scales etc in a dictionary - the dictionary form allows for added functionaility. If bits=0, find the natural accuracy. eps defaults to 3e-6, and is the error relative to the largest element, as is maxerror.

pyfusion.data.save_compress.discretise_signal(timebase=None, signal=None, parent_element=array(0), eps=0, verbose=0, params={}, delta_encode_time=True, delta_encode_signal=False, filename=None)

a function to return a dictionary from signal and timebase, with relative accuracy eps, optionally saving if filename is defined. Achieves a factor of >10x on MP1 signal 33373 using delta_encode_time=True Delta encode on signal is not effective for MP and MICROFAST (.005% worse) Probably should eventually separate the file write from making the dictionary. Intended to be aliased with loadz, args slightly different. There is no dependence on pyfusion. Version 101 adds time_unit_in_seconds version 102 adds utc, raw, 103 after correction of probes 11-20 Note: changed to parent_element=array(0) by default - not sure what this is!

pyfusion.data.save_compress.newload(filename, verbose=0)

Intended to replace load() in numpy This is being used with nan data. The version in data/base.py is closer to python 3 compatible, but can’t deal with the nans yet.

pyfusion.data.save_compress.save_compress(timebase=None, signal=None, filename=None, *args, **kwargs)

save a signal and timebase into a compress .npz file. See arglist of discretise_signal. Example:

>>> sig=[1,2,1,2] ; tb=[1,2,3,4] # need this only for later comparison
>>> save_compress(timebase=tb, signal=sig, filename='junk')
>>> readback=newload('junk.npz',verbose=0)
>>> if (readback['signal'] != sig).any(): print 'error in save/restore'

pyfusion.data.save_compress.test_compress(file=None, verbose=0, eps=0, debug=False, maxcount=0)

Used in developing the save compress routines. Not tested since then

>>> test_compress()

Looks like it only saves the time series, not the rest.

pyfusion.data.save_compress.try_discretise_array(arr, eps=0, bits=0, deltar=None, verbose=0, delta_encode=False)

Return an integer array and scales etc in a dictionary - the dictionary form allows for added functionality. If bits=0, find the natural accuracy. eps defaults to 1e-6

shot_range Module

pyfusion.data.shot_range.next_shot(shot)

pyfusion.data.shot_range.shot_gte(shot1, shot2)

pyfusion.data.shot_range.shot_range(shot_from, shot_to)

signal_processing Module

Boyd’s python for stand alone, general signal processing, try to be “efficient”. Stuff that imports pyfusion should live in pyfusion_sigproc.py. Put stuff here so that recompilation doesn’t require restarting pyfusion.

pyfusion.data.signal_processing.analytic_phase(x, t=None, subint=None)

gets the phase from an amazing variety of signals http://en.wikipedia.org/wiki/Analytic_signal subinterval idea is not debugged and is probably unnecessary may shorten data?

pyfusion.data.signal_processing.cross_correl(x1, x2, nsmooth=21, n_times=3)

<x1.x2>/sqrt(<x1.x1> * <x2,x2>) averaged over nsmooth points n_times Ideally extract raw data from pyfusion signals, but that makes it pyfusion specfic - bad.

pyfusion.data.signal_processing.powerof2(n, near=False)

pyfusion.data.signal_processing.smooth(data, n_smooth=3, timebase=None, causal=False, indices=False, keep=False)

An efficient top hat smoother based on the IDL routine of that name. The use of cumsum-shift(cumsum) means that execution time is 2xN flops compared to 2 x n_smooth x N for a convolution. If supplied with a timebase, the shortened timebase is returned as the first of a tuple.

causal – If true, the smoothed signal never preceded the input,

otherwise, the smoothed signal is “centred” on the input (for n_smooth odd) and close (1/2 timestep off) for even

indices – if true, return the timebase indices instead of the times data = (timebase, data) is a shorthand way to pass timebase n_smooth - apply recursively if an array e.g. n_smooth=[33,20,14,11]

removes 3rd, 5th, 7th, 9th harmonics fraction values are interpreted as timeintervals.

>>> smooth([1,2,3,4],3)
array([ 2.,  3.])
>>> smooth([1.,2.,3.,4.,5.],3)
array([ 2.,  3.,  4.])
>>> smooth([1,2,3,4,5],timebase=array([1,2,3,4,5]),n_smooth=3, causal=False)
(array([2, 3, 4]), array([ 2.,  3.,  4.]))
>>> smooth([0,0,0,3,0,0,0],timebase=[1,2,3,4,5,6,7],n_smooth=3, causal=True)
([3, 4, 5, 6, 7], array([ 0.,  1.,  1.,  1.,  0.]))
>>> smooth([0,0,0,3,0,0,0],timebase=[1,2,3,4,5,6,7],n_smooth=3, causal=True, indices=True)
([2, 3, 4, 5, 6], array([ 0.,  1.,  1.,  1.,  0.]))
>>> smooth([0,   0,   0,   0,   5,   0,    0,  0,   0,   0,   0], 5, keep=1)
array([ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  0.,  0., -1., -1.])

Last example: keep=1:

Better to throw the partially cooked ends away, but if you want to

keep them use keep=True. THis is useful for quick filtering applications so that original and filtered signals are easily compared without worrying about timebase

pyfusion.data.signal_processing.smooth_n(data, n_smooth=3, timebase=None, causal=False, iter=3, keep=False, indices=False)

Apply smooth “iter” times. [ smooth() doc follows: ]

An efficient top hat smoother based on the IDL routine of that name.

The use of cumsum-shift(cumsum) means that execution time is 2xN flops compared to 2 x n_smooth x N for a convolution. If supplied with a timebase, the shortened timebase is returned as the first of a tuple.

causal – If true, the smoothed signal never preceded the input,

otherwise, the smoothed signal is “centred” on the input (for n_smooth odd) and close (1/2 timestep off) for even

indices – if true, return the timebase indices instead of the times data = (timebase, data) is a shorthand way to pass timebase n_smooth - apply recursively if an array e.g. n_smooth=[33,20,14,11]

removes 3rd, 5th, 7th, 9th harmonics fraction values are interpreted as timeintervals.

>>> smooth([1,2,3,4],3)
array([ 2.,  3.])
>>> smooth([1.,2.,3.,4.,5.],3)
array([ 2.,  3.,  4.])
>>> smooth([1,2,3,4,5],timebase=array([1,2,3,4,5]),n_smooth=3, causal=False)
(array([2, 3, 4]), array([ 2.,  3.,  4.]))
>>> smooth([0,0,0,3,0,0,0],timebase=[1,2,3,4,5,6,7],n_smooth=3, causal=True)
([3, 4, 5, 6, 7], array([ 0.,  1.,  1.,  1.,  0.]))
>>> smooth([0,0,0,3,0,0,0],timebase=[1,2,3,4,5,6,7],n_smooth=3, causal=True, indices=True)
([2, 3, 4, 5, 6], array([ 0.,  1.,  1.,  1.,  0.]))
>>> smooth([0,   0,   0,   0,   5,   0,    0,  0,   0,   0,   0], 5, keep=1)
array([ 0.,  0.,  1.,  1.,  1.,  1.,  1.,  0.,  0., -1., -1.])

Last example: keep=1:

Better to throw the partially cooked ends away, but if you want to

keep them use keep=True. THis is useful for quick filtering applications so that original and filtered signals are easily compared without worrying about timebase

pyfusion.data.signal_processing.splot(signal, *args, **kwargs)

simple wrapper to plot to accept tuples of form (timebase, data)

pyfusion.data.signal_processing.test_analytic_phase(verbose=3)

>>> test_analytic_phase()

specgram_fftw3 Module

A faster replacement for matplotlib specgram using float32 fftw3 Real-Complex routines. So far about 4x faster for no overlap, 3x faster with overlap Doesn’t include time from return to prompt until display is updated. and times do not have SIMD enabled (fft is not the big cost)

Strategy:

1/ Use fftw3 -> factor of 2-3. (see multiprocessing) 2/ Use Real-Complex 32 bit precision ffts. another factor of 2 (looks like psd already does this) 3/ Take log, abs and do window on a composite array for speed, and avoid sqrt

this is effective as chunks are small (128-1024) so numpy overhead per call of a few us can be important.

Expect about 2.3us (E4300, simd, 1thread, and 1.5 lenny) for a float32 real 512 pt fft

–> so for 2**20 array, 2000*2.3 = 5ms!

However the pyfftw overhead of about 10us has in impact for small bins like 512. The more direct interface of the other fftw3 could help here, but you would need to do copies on input and output - these take 1.5us each for float32, 512 bytes, so would be a nett win of 5us*512 =2.5ms for 512x512 case (15%) and 40ms for 2**22 (8192) = 13%

'FFTW_ESTIMATE', 'FFTW_MEASURE' (default), 'FFTW_PATIENT' and 'FFTW_EXHAUSTIVE'.

Improvements: Use the Advanced interface plan_many to get all short samples done all at once. Hopefully threads are allocated one chunk at a time.

Timings for no overlap, E4300 (power: no diff) clf() means preceded by clf, hold=1, using original (no simd) ubuntu libraries.

specgram This

512x512 512x8192 512x512 512x8192 512x32768

clf() 96 804 69 360 hold=0 125 765 62 351 hold=1 43 695 17.7 307 1486

Sun 21: 360ms for hold=0, 512*8192? instead of 351 cpu is showing 10% when idle, and 86 ms in imshow - no gain by reboot (don’t forget to change hold= in the imshow line)

With FFTW fft instead of numpy replacement, now 150 for ffts, 198 total (hold=0) (noverlap=512, NFFT=1024 - 426 (incl imshow) 336-360 before)

This compares with the raw fft time (noverlap=0) of 8k*2.3 = 18ms Summary of times:

window mult inp= F.execute allout= log imshow

7ms 16ms 15ms 11ms 8ms 77ms 50ms 184 cf 197actual

See /usr/lib/pymodules/python2.7/matplotlib/mlab.py

Can test animate speed by im_obj.set_data(im) # 50ms (std) 80ms 1200x800 im_obj.set_data(im[::5,::20]) # 30ms 60ms 1200x800 time for ii in range(10):draw()

pyfusion.data.specgram_fftw3.all_tests_specgram(slow=False)

compare with pylab specgram for random noise data size=2**12 (and 2*22 for slow=True)

pyfusion.data.specgram_fftw3.specgram(x, NFFT=256, Fs=2, Fc=0, detrend=None, window=<function hanning>, noverlap=128, cmap=None, hold=None, dtype=<type 'numpy.float32'>, threads=1, im_obj=None, interpolation='nearest', pylab_scaling=True, fast=True)

Return and optionally plot the spectrogram of the data in x, using fftw3 to

optimise the speed. Note that using fftw3, most of the time is spent elsewhere, so perhaps a cython implementation would be more efficient.

So far only coded for float32 and float64

Slight increase in speed if the previous image_object is given in the arg list.

Tests below on E4300, power on, 1 thread

• 6 secs til image for specgram of 16M points, NFFT=512, noverlap=128 std win

• 1.8 secs til first image fast=1, im_obj=im_obj (1 sec more for float64)

• e.g. (note this is a wide dynamic range (17 decades of power) so there

  will be slight errors visible in the float32 version.

  >>> time figure(3);s=specgram(arange(2**24)**n,Fs=1e6,NFFT=512,im_obj=s[3],dtype=float64);n=(n+1)%3
  

  matplotlib specgram has good dyn range for arange2*24, but error of ~ -8 for random

pyfusion.data.specgram_fftw3.test_specgram(x=array([ 0.11371453, 0.58468941, 0.58165886, ..., 0.23541809, 0.90983162, 0.537599 ]), NFFT=256, noverlap=128, dtype='float32', rtol=1e-06, atol=3e-07, verbose=1)

compare with pylab specgram for random noise atol=2e-7 shows errors < 1 in 10,000 times - verbose not yet implemented don’t run too many of these (i.e. > 1000)- eats into RAM

tests Module

Tests for the Data and related classes.

class pyfusion.data.tests.CheckChannelList

Bases: pyfusion.test.tests.PfTestBase

test_channel_list()

test_channellist_ORM()

class pyfusion.data.tests.CheckChannels

Bases: pyfusion.test.tests.PfTestBase

Make sure that arguments passed to Channel() appear as attributes.

test_channel_class()

class pyfusion.data.tests.CheckChannelsSQL

Bases: pyfusion.test.tests.PfTestBase

test_channels_SQL()

class pyfusion.data.tests.CheckCoordinates

Bases: pyfusion.test.tests.PfTestBase

Check that we can add and transform coordinates.

test_add_coords()

test_coord_transform()

class pyfusion.data.tests.CheckDataHistory

Bases: pyfusion.test.tests.PfTestBase

testFilteredDataHistory_copy()

make sure that _copy version does NOT alter original

testFilteredDataHistory_nocopy()

testNewData()

class pyfusion.data.tests.CheckDataSet

Bases: pyfusion.test.tests.PfTestBase

test_dataset()

test_dataset_filters_2()

class pyfusion.data.tests.CheckDataSetLabels

Bases: pyfusion.test.tests.PfTestBase

test_baseordereddataset_label()

test_dataset_label()

class pyfusion.data.tests.CheckFilterCopy

Bases: pyfusion.test.tests.PfTestBase

Check that by default, data filters alter a copy of the input data object not the object itself.

dev = False

test_dataset_filter_copy()

test_dataset_filter_nocopy()

test_timeseries_filter_copy()

test_timeseries_filter_nocopy()

class pyfusion.data.tests.CheckFilterMetaClass

Bases: pyfusion.test.tests.PfTestBase

test_new_filter()

class pyfusion.data.tests.CheckFilters

Bases: pyfusion.test.tests.PfTestBase

test_reduce_time_dataset()

test_reduce_time_filter_multi_channel()

test_reduce_time_filter_multi_channel_attached_method()

test_reduce_time_filter_single_channel()

class pyfusion.data.tests.CheckFloatDelta

Bases: pyfusion.test.tests.PfTestBase

delta phase data class.

test_ORM_floatdelta()

check that floatdelta can be saved to database

test_d_phase()

class pyfusion.data.tests.CheckFlucstrucPhases

Bases: pyfusion.test.tests.PfTestBase

Test code to replicate bug found by Shaun

Pyfusion uses different syntax depending if sql is enabled.

test_flucstruc_phases(PfTestCase)

class pyfusion.data.tests.CheckFlucstrucs

Bases: pyfusion.test.tests.PfTestBase

test_ORM_flucstrucs()

check that flucstrucs can be saved to database

test_SVDData_class()

test_flucstruc_phases()

test_flucstruc_signals()

test_svd_data()

class pyfusion.data.tests.CheckGetCoords

Bases: pyfusion.test.tests.PfTestBase

test_get_coords_for_channel_config()

class pyfusion.data.tests.CheckNormalise

Bases: pyfusion.test.tests.PfTestBase

test_multichannel_fakedata()

test_single_channel_fakedata()

class pyfusion.data.tests.CheckNumpyFilters

Bases: pyfusion.test.tests.PfTestBase

test_correlate()

class pyfusion.data.tests.CheckOrderedDataSet

Bases: pyfusion.test.tests.PfTestBase

test the ordered dataset

test_ordered_dataset_ORM()

class pyfusion.data.tests.CheckPlotMethods

Bases: pyfusion.test.tests.PfTestBase

test_svd_plot()

class pyfusion.data.tests.CheckPlotSignals

Bases: pyfusion.test.tests.PfTestBase

broken = True

mds = True

slow = True

test_plot_signals()

class pyfusion.data.tests.CheckRemoveNonContiguousFilter

Bases: pyfusion.test.tests.PfTestBase

test_remove_noncontiguous()

class pyfusion.data.tests.CheckSciPyFilters

Bases: pyfusion.test.tests.PfTestBase

test_sp_filter_butterworth_bandpass()

class pyfusion.data.tests.CheckSegmentFilter

Bases: pyfusion.test.tests.PfTestBase

test_dataset()

test_multi_channel_timeseries()

test_single_channel_timeseries()

class pyfusion.data.tests.CheckShotFlucstrucs

Bases: pyfusion.test.tests.PfTestBase

test_shot_flucstrucs()

Just check that the number of flucstrucs is the same whether we use flucstruc directly on the shot data with the segment kwarg or whether we explicitly call the segment method.

class pyfusion.data.tests.CheckSignal

Bases: pyfusion.test.tests.PfTestBase

Test Signal class.

test_base_class()

test_n_channels()

test_n_samples()

class pyfusion.data.tests.CheckStoredMetaData

Bases: pyfusion.test.tests.PfTestBase

test_stored_metadata_data()

metadata should be stored to data instances, rather than datasets - this might be slower, but more likely to guarantee data is kept track of.

class pyfusion.data.tests.CheckStoredMetaDataForDataSets

Bases: pyfusion.test.tests.PfTestBase

test_stored_metadata_datasets()

Make sure metadata attached to dataset classes is saved to sql.

class pyfusion.data.tests.CheckSubtractMeanFilter

Bases: pyfusion.test.tests.PfTestBase

Test mean subtraction filter for timeseries data.

test_remove_mean_dataset()

test_remove_mean_multichanel()

test_remove_mean_single_channel()

class pyfusion.data.tests.CheckTimebase

Bases: pyfusion.test.tests.PfTestBase

Test Timebase class.

test_timebase()

test_timebase_slice()

class pyfusion.data.tests.CheckTimeseriesData

Bases: pyfusion.test.tests.PfTestBase

Test timeseries data

testBaseClasses()

test_timebase_and_coords()

class pyfusion.data.tests.CheckUtils

Bases: pyfusion.test.tests.PfTestBase

Test the helper functions in pyfusion.data.utils.py

test_peak_freq()

test_remap_periodic()

class pyfusion.data.tests.DummyCoordTransform

Bases: pyfusion.data.base.BaseCoordTransform

Minimal coordinate transform class for testing.

Transforms cylindrical coordinates (a,b,c) to dummy coordinates (2a,3b,4c)

input_coords = 'cylindrical'

output_coords = 'dummy'

transform(coords)

pyfusion.data.tests.get_multimode_test_data(channels=[<pyfusion.data.base.Channel object at 0x2b036ef76710>, <pyfusion.data.base.Channel object at 0x2b036ef76790>, <pyfusion.data.base.Channel object at 0x2b036ef76810>, <pyfusion.data.base.Channel object at 0x2b036ef76890>, <pyfusion.data.base.Channel object at 0x2b036ef76910>, <pyfusion.data.base.Channel object at 0x2b036ef76990>, <pyfusion.data.base.Channel object at 0x2b036ef76a10>, <pyfusion.data.base.Channel object at 0x2b036ef76a90>, <pyfusion.data.base.Channel object at 0x2b036ef76b10>, <pyfusion.data.base.Channel object at 0x2b036ef76b90>], timebase=Timebase([ 0.00000000e+00, 1.00000000e-05, 2.00000000e-05, 3.00000000e-05, 4.00000000e-05, 5.00000000e-05, 6.00000000e-05, 7.00000000e-05, 8.00000000e-05, 9.00000000e-05, 1.00000000e-04, 1.10000000e-04, 1.20000000e-04, 1.30000000e-04, 1.40000000e-04, 1.50000000e-04, 1.60000000e-04, 1.70000000e-04, 1.80000000e-04, 1.90000000e-04, 2.00000000e-04, 2.10000000e-04, 2.20000000e-04, 2.30000000e-04, 2.40000000e-04, 2.50000000e-04, 2.60000000e-04, 2.70000000e-04, 2.80000000e-04, 2.90000000e-04, 3.00000000e-04, 3.10000000e-04, 3.20000000e-04, 3.30000000e-04, 3.40000000e-04, 3.50000000e-04, 3.60000000e-04, 3.70000000e-04, 3.80000000e-04, 3.90000000e-04, 4.00000000e-04, 4.10000000e-04, 4.20000000e-04, 4.30000000e-04, 4.40000000e-04, 4.50000000e-04, 4.60000000e-04, 4.70000000e-04, 4.80000000e-04, 4.90000000e-04, 5.00000000e-04, 5.10000000e-04, 5.20000000e-04, 5.30000000e-04, 5.40000000e-04, 5.50000000e-04, 5.60000000e-04, 5.70000000e-04, 5.80000000e-04, 5.90000000e-04, 6.00000000e-04, 6.10000000e-04, 6.20000000e-04, 6.30000000e-04, 6.40000000e-04, 6.50000000e-04, 6.60000000e-04, 6.70000000e-04, 6.80000000e-04, 6.90000000e-04, 7.00000000e-04, 7.10000000e-04, 7.20000000e-04, 7.30000000e-04, 7.40000000e-04, 7.50000000e-04, 7.60000000e-04, 7.70000000e-04, 7.80000000e-04, 7.90000000e-04, 8.00000000e-04, 8.10000000e-04, 8.20000000e-04, 8.30000000e-04, 8.40000000e-04, 8.50000000e-04, 8.60000000e-04, 8.70000000e-04, 8.80000000e-04, 8.90000000e-04, 9.00000000e-04, 9.10000000e-04, 9.20000000e-04, 9.30000000e-04, 9.40000000e-04, 9.50000000e-04, 9.60000000e-04, 9.70000000e-04, 9.80000000e-04, 9.90000000e-04, 1.00000000e-03, 1.01000000e-03, 1.02000000e-03, 1.03000000e-03, 1.04000000e-03, 1.05000000e-03, 1.06000000e-03, 1.07000000e-03, 1.08000000e-03, 1.09000000e-03, 1.10000000e-03, 1.11000000e-03, 1.12000000e-03, 1.13000000e-03, 1.14000000e-03, 1.15000000e-03, 1.16000000e-03, 1.17000000e-03, 1.18000000e-03, 1.19000000e-03, 1.20000000e-03, 1.21000000e-03, 1.22000000e-03, 1.23000000e-03, 1.24000000e-03, 1.25000000e-03, 1.26000000e-03, 1.27000000e-03, 1.28000000e-03, 1.29000000e-03, 1.30000000e-03, 1.31000000e-03, 1.32000000e-03, 1.33000000e-03, 1.34000000e-03, 1.35000000e-03, 1.36000000e-03, 1.37000000e-03, 1.38000000e-03, 1.39000000e-03, 1.40000000e-03, 1.41000000e-03, 1.42000000e-03, 1.43000000e-03, 1.44000000e-03, 1.45000000e-03, 1.46000000e-03, 1.47000000e-03, 1.48000000e-03, 1.49000000e-03, 1.50000000e-03, 1.51000000e-03, 1.52000000e-03, 1.53000000e-03, 1.54000000e-03, 1.55000000e-03, 1.56000000e-03, 1.57000000e-03, 1.58000000e-03, 1.59000000e-03, 1.60000000e-03, 1.61000000e-03, 1.62000000e-03, 1.63000000e-03, 1.64000000e-03, 1.65000000e-03, 1.66000000e-03, 1.67000000e-03, 1.68000000e-03, 1.69000000e-03, 1.70000000e-03, 1.71000000e-03, 1.72000000e-03, 1.73000000e-03, 1.74000000e-03, 1.75000000e-03, 1.76000000e-03, 1.77000000e-03, 1.78000000e-03, 1.79000000e-03, 1.80000000e-03, 1.81000000e-03, 1.82000000e-03, 1.83000000e-03, 1.84000000e-03, 1.85000000e-03, 1.86000000e-03, 1.87000000e-03, 1.88000000e-03, 1.89000000e-03, 1.90000000e-03, 1.91000000e-03, 1.92000000e-03, 1.93000000e-03, 1.94000000e-03, 1.95000000e-03, 1.96000000e-03, 1.97000000e-03, 1.98000000e-03, 1.99000000e-03, 2.00000000e-03, 2.01000000e-03, 2.02000000e-03, 2.03000000e-03, 2.04000000e-03, 2.05000000e-03, 2.06000000e-03, 2.07000000e-03, 2.08000000e-03, 2.09000000e-03, 2.10000000e-03, 2.11000000e-03, 2.12000000e-03, 2.13000000e-03, 2.14000000e-03, 2.15000000e-03, 2.16000000e-03, 2.17000000e-03, 2.18000000e-03, 2.19000000e-03, 2.20000000e-03, 2.21000000e-03, 2.22000000e-03, 2.23000000e-03, 2.24000000e-03, 2.25000000e-03, 2.26000000e-03, 2.27000000e-03, 2.28000000e-03, 2.29000000e-03, 2.30000000e-03, 2.31000000e-03, 2.32000000e-03, 2.33000000e-03, 2.34000000e-03, 2.35000000e-03, 2.36000000e-03, 2.37000000e-03, 2.38000000e-03, 2.39000000e-03, 2.40000000e-03, 2.41000000e-03, 2.42000000e-03, 2.43000000e-03, 2.44000000e-03, 2.45000000e-03, 2.46000000e-03, 2.47000000e-03, 2.48000000e-03, 2.49000000e-03, 2.50000000e-03, 2.51000000e-03, 2.52000000e-03, 2.53000000e-03, 2.54000000e-03, 2.55000000e-03, 2.56000000e-03, 2.57000000e-03, 2.58000000e-03, 2.59000000e-03, 2.60000000e-03, 2.61000000e-03, 2.62000000e-03, 2.63000000e-03, 2.64000000e-03, 2.65000000e-03, 2.66000000e-03, 2.67000000e-03, 2.68000000e-03, 2.69000000e-03, 2.70000000e-03, 2.71000000e-03, 2.72000000e-03, 2.73000000e-03, 2.74000000e-03, 2.75000000e-03, 2.76000000e-03, 2.77000000e-03, 2.78000000e-03, 2.79000000e-03, 2.80000000e-03, 2.81000000e-03, 2.82000000e-03, 2.83000000e-03, 2.84000000e-03, 2.85000000e-03, 2.86000000e-03, 2.87000000e-03, 2.88000000e-03, 2.89000000e-03, 2.90000000e-03, 2.91000000e-03, 2.92000000e-03, 2.93000000e-03, 2.94000000e-03, 2.95000000e-03, 2.96000000e-03, 2.97000000e-03, 2.98000000e-03, 2.99000000e-03, 3.00000000e-03, 3.01000000e-03, 3.02000000e-03, 3.03000000e-03, 3.04000000e-03, 3.05000000e-03, 3.06000000e-03, 3.07000000e-03, 3.08000000e-03, 3.09000000e-03, 3.10000000e-03, 3.11000000e-03, 3.12000000e-03, 3.13000000e-03, 3.14000000e-03, 3.15000000e-03, 3.16000000e-03, 3.17000000e-03, 3.18000000e-03, 3.19000000e-03, 3.20000000e-03, 3.21000000e-03, 3.22000000e-03, 3.23000000e-03, 3.24000000e-03, 3.25000000e-03, 3.26000000e-03, 3.27000000e-03, 3.28000000e-03, 3.29000000e-03, 3.30000000e-03, 3.31000000e-03, 3.32000000e-03, 3.33000000e-03, 3.34000000e-03, 3.35000000e-03, 3.36000000e-03, 3.37000000e-03, 3.38000000e-03, 3.39000000e-03, 3.40000000e-03, 3.41000000e-03, 3.42000000e-03, 3.43000000e-03, 3.44000000e-03, 3.45000000e-03, 3.46000000e-03, 3.47000000e-03, 3.48000000e-03, 3.49000000e-03, 3.50000000e-03, 3.51000000e-03, 3.52000000e-03, 3.53000000e-03, 3.54000000e-03, 3.55000000e-03, 3.56000000e-03, 3.57000000e-03, 3.58000000e-03, 3.59000000e-03, 3.60000000e-03, 3.61000000e-03, 3.62000000e-03, 3.63000000e-03, 3.64000000e-03, 3.65000000e-03, 3.66000000e-03, 3.67000000e-03, 3.68000000e-03, 3.69000000e-03, 3.70000000e-03, 3.71000000e-03, 3.72000000e-03, 3.73000000e-03, 3.74000000e-03, 3.75000000e-03, 3.76000000e-03, 3.77000000e-03, 3.78000000e-03, 3.79000000e-03, 3.80000000e-03, 3.81000000e-03, 3.82000000e-03, 3.83000000e-03, 3.84000000e-03, 3.85000000e-03, 3.86000000e-03, 3.87000000e-03, 3.88000000e-03, 3.89000000e-03, 3.90000000e-03, 3.91000000e-03, 3.92000000e-03, 3.93000000e-03, 3.94000000e-03, 3.95000000e-03, 3.96000000e-03, 3.97000000e-03, 3.98000000e-03, 3.99000000e-03, 4.00000000e-03, 4.01000000e-03, 4.02000000e-03, 4.03000000e-03, 4.04000000e-03, 4.05000000e-03, 4.06000000e-03, 4.07000000e-03, 4.08000000e-03, 4.09000000e-03, 4.10000000e-03, 4.11000000e-03, 4.12000000e-03, 4.13000000e-03, 4.14000000e-03, 4.15000000e-03, 4.16000000e-03, 4.17000000e-03, 4.18000000e-03, 4.19000000e-03, 4.20000000e-03, 4.21000000e-03, 4.22000000e-03, 4.23000000e-03, 4.24000000e-03, 4.25000000e-03, 4.26000000e-03, 4.27000000e-03, 4.28000000e-03, 4.29000000e-03, 4.30000000e-03, 4.31000000e-03, 4.32000000e-03, 4.33000000e-03, 4.34000000e-03, 4.35000000e-03, 4.36000000e-03, 4.37000000e-03, 4.38000000e-03, 4.39000000e-03, 4.40000000e-03, 4.41000000e-03, 4.42000000e-03, 4.43000000e-03, 4.44000000e-03, 4.45000000e-03, 4.46000000e-03, 4.47000000e-03, 4.48000000e-03, 4.49000000e-03, 4.50000000e-03, 4.51000000e-03, 4.52000000e-03, 4.53000000e-03, 4.54000000e-03, 4.55000000e-03, 4.56000000e-03, 4.57000000e-03, 4.58000000e-03, 4.59000000e-03, 4.60000000e-03, 4.61000000e-03, 4.62000000e-03, 4.63000000e-03, 4.64000000e-03, 4.65000000e-03, 4.66000000e-03, 4.67000000e-03, 4.68000000e-03, 4.69000000e-03, 4.70000000e-03, 4.71000000e-03, 4.72000000e-03, 4.73000000e-03, 4.74000000e-03, 4.75000000e-03, 4.76000000e-03, 4.77000000e-03, 4.78000000e-03, 4.79000000e-03, 4.80000000e-03, 4.81000000e-03, 4.82000000e-03, 4.83000000e-03, 4.84000000e-03, 4.85000000e-03, 4.86000000e-03, 4.87000000e-03, 4.88000000e-03, 4.89000000e-03, 4.90000000e-03, 4.91000000e-03, 4.92000000e-03, 4.93000000e-03, 4.94000000e-03, 4.95000000e-03, 4.96000000e-03, 4.97000000e-03, 4.98000000e-03, 4.99000000e-03, 5.00000000e-03, 5.01000000e-03, 5.02000000e-03, 5.03000000e-03, 5.04000000e-03, 5.05000000e-03, 5.06000000e-03, 5.07000000e-03, 5.08000000e-03, 5.09000000e-03, 5.10000000e-03, 5.11000000e-03, 5.12000000e-03, 5.13000000e-03, 5.14000000e-03, 5.15000000e-03, 5.16000000e-03, 5.17000000e-03, 5.18000000e-03, 5.19000000e-03, 5.20000000e-03, 5.21000000e-03, 5.22000000e-03, 5.23000000e-03, 5.24000000e-03, 5.25000000e-03, 5.26000000e-03, 5.27000000e-03, 5.28000000e-03, 5.29000000e-03, 5.30000000e-03, 5.31000000e-03, 5.32000000e-03, 5.33000000e-03, 5.34000000e-03, 5.35000000e-03, 5.36000000e-03, 5.37000000e-03, 5.38000000e-03, 5.39000000e-03, 5.40000000e-03, 5.41000000e-03, 5.42000000e-03, 5.43000000e-03, 5.44000000e-03, 5.45000000e-03, 5.46000000e-03, 5.47000000e-03, 5.48000000e-03, 5.49000000e-03, 5.50000000e-03, 5.51000000e-03, 5.52000000e-03, 5.53000000e-03, 5.54000000e-03, 5.55000000e-03, 5.56000000e-03, 5.57000000e-03, 5.58000000e-03, 5.59000000e-03, 5.60000000e-03, 5.61000000e-03, 5.62000000e-03, 5.63000000e-03, 5.64000000e-03, 5.65000000e-03, 5.66000000e-03, 5.67000000e-03, 5.68000000e-03, 5.69000000e-03, 5.70000000e-03, 5.71000000e-03, 5.72000000e-03, 5.73000000e-03, 5.74000000e-03, 5.75000000e-03, 5.76000000e-03, 5.77000000e-03, 5.78000000e-03, 5.79000000e-03, 5.80000000e-03, 5.81000000e-03, 5.82000000e-03, 5.83000000e-03, 5.84000000e-03, 5.85000000e-03, 5.86000000e-03, 5.87000000e-03, 5.88000000e-03, 5.89000000e-03, 5.90000000e-03, 5.91000000e-03, 5.92000000e-03, 5.93000000e-03, 5.94000000e-03, 5.95000000e-03, 5.96000000e-03, 5.97000000e-03, 5.98000000e-03, 5.99000000e-03, 6.00000000e-03, 6.01000000e-03, 6.02000000e-03, 6.03000000e-03, 6.04000000e-03, 6.05000000e-03, 6.06000000e-03, 6.07000000e-03, 6.08000000e-03, 6.09000000e-03, 6.10000000e-03, 6.11000000e-03, 6.12000000e-03, 6.13000000e-03, 6.14000000e-03, 6.15000000e-03, 6.16000000e-03, 6.17000000e-03, 6.18000000e-03, 6.19000000e-03, 6.20000000e-03, 6.21000000e-03, 6.22000000e-03, 6.23000000e-03, 6.24000000e-03, 6.25000000e-03, 6.26000000e-03, 6.27000000e-03, 6.28000000e-03, 6.29000000e-03, 6.30000000e-03, 6.31000000e-03, 6.32000000e-03, 6.33000000e-03, 6.34000000e-03, 6.35000000e-03, 6.36000000e-03, 6.37000000e-03, 6.38000000e-03, 6.39000000e-03, 6.40000000e-03, 6.41000000e-03, 6.42000000e-03, 6.43000000e-03, 6.44000000e-03, 6.45000000e-03, 6.46000000e-03, 6.47000000e-03, 6.48000000e-03, 6.49000000e-03, 6.50000000e-03, 6.51000000e-03, 6.52000000e-03, 6.53000000e-03, 6.54000000e-03, 6.55000000e-03, 6.56000000e-03, 6.57000000e-03, 6.58000000e-03, 6.59000000e-03, 6.60000000e-03, 6.61000000e-03, 6.62000000e-03, 6.63000000e-03, 6.64000000e-03, 6.65000000e-03, 6.66000000e-03, 6.67000000e-03, 6.68000000e-03, 6.69000000e-03, 6.70000000e-03, 6.71000000e-03, 6.72000000e-03, 6.73000000e-03, 6.74000000e-03, 6.75000000e-03, 6.76000000e-03, 6.77000000e-03, 6.78000000e-03, 6.79000000e-03, 6.80000000e-03, 6.81000000e-03, 6.82000000e-03, 6.83000000e-03, 6.84000000e-03, 6.85000000e-03, 6.86000000e-03, 6.87000000e-03, 6.88000000e-03, 6.89000000e-03, 6.90000000e-03, 6.91000000e-03, 6.92000000e-03, 6.93000000e-03, 6.94000000e-03, 6.95000000e-03, 6.96000000e-03, 6.97000000e-03, 6.98000000e-03, 6.99000000e-03, 7.00000000e-03, 7.01000000e-03, 7.02000000e-03, 7.03000000e-03, 7.04000000e-03, 7.05000000e-03, 7.06000000e-03, 7.07000000e-03, 7.08000000e-03, 7.09000000e-03, 7.10000000e-03, 7.11000000e-03, 7.12000000e-03, 7.13000000e-03, 7.14000000e-03, 7.15000000e-03, 7.16000000e-03, 7.17000000e-03, 7.18000000e-03, 7.19000000e-03, 7.20000000e-03, 7.21000000e-03, 7.22000000e-03, 7.23000000e-03, 7.24000000e-03, 7.25000000e-03, 7.26000000e-03, 7.27000000e-03, 7.28000000e-03, 7.29000000e-03, 7.30000000e-03, 7.31000000e-03, 7.32000000e-03, 7.33000000e-03, 7.34000000e-03, 7.35000000e-03, 7.36000000e-03, 7.37000000e-03, 7.38000000e-03, 7.39000000e-03, 7.40000000e-03, 7.41000000e-03, 7.42000000e-03, 7.43000000e-03, 7.44000000e-03, 7.45000000e-03, 7.46000000e-03, 7.47000000e-03, 7.48000000e-03, 7.49000000e-03, 7.50000000e-03, 7.51000000e-03, 7.52000000e-03, 7.53000000e-03, 7.54000000e-03, 7.55000000e-03, 7.56000000e-03, 7.57000000e-03, 7.58000000e-03, 7.59000000e-03, 7.60000000e-03, 7.61000000e-03, 7.62000000e-03, 7.63000000e-03, 7.64000000e-03, 7.65000000e-03, 7.66000000e-03, 7.67000000e-03, 7.68000000e-03, 7.69000000e-03, 7.70000000e-03, 7.71000000e-03, 7.72000000e-03, 7.73000000e-03, 7.74000000e-03, 7.75000000e-03, 7.76000000e-03, 7.77000000e-03, 7.78000000e-03, 7.79000000e-03, 7.80000000e-03, 7.81000000e-03, 7.82000000e-03, 7.83000000e-03, 7.84000000e-03, 7.85000000e-03, 7.86000000e-03, 7.87000000e-03, 7.88000000e-03, 7.89000000e-03, 7.90000000e-03, 7.91000000e-03, 7.92000000e-03, 7.93000000e-03, 7.94000000e-03, 7.95000000e-03, 7.96000000e-03, 7.97000000e-03, 7.98000000e-03, 7.99000000e-03, 8.00000000e-03, 8.01000000e-03, 8.02000000e-03, 8.03000000e-03, 8.04000000e-03, 8.05000000e-03, 8.06000000e-03, 8.07000000e-03, 8.08000000e-03, 8.09000000e-03, 8.10000000e-03, 8.11000000e-03, 8.12000000e-03, 8.13000000e-03, 8.14000000e-03, 8.15000000e-03, 8.16000000e-03, 8.17000000e-03, 8.18000000e-03, 8.19000000e-03, 8.20000000e-03, 8.21000000e-03, 8.22000000e-03, 8.23000000e-03, 8.24000000e-03, 8.25000000e-03, 8.26000000e-03, 8.27000000e-03, 8.28000000e-03, 8.29000000e-03, 8.30000000e-03, 8.31000000e-03, 8.32000000e-03, 8.33000000e-03, 8.34000000e-03, 8.35000000e-03, 8.36000000e-03, 8.37000000e-03, 8.38000000e-03, 8.39000000e-03, 8.40000000e-03, 8.41000000e-03, 8.42000000e-03, 8.43000000e-03, 8.44000000e-03, 8.45000000e-03, 8.46000000e-03, 8.47000000e-03, 8.48000000e-03, 8.49000000e-03, 8.50000000e-03, 8.51000000e-03, 8.52000000e-03, 8.53000000e-03, 8.54000000e-03, 8.55000000e-03, 8.56000000e-03, 8.57000000e-03, 8.58000000e-03, 8.59000000e-03, 8.60000000e-03, 8.61000000e-03, 8.62000000e-03, 8.63000000e-03, 8.64000000e-03, 8.65000000e-03, 8.66000000e-03, 8.67000000e-03, 8.68000000e-03, 8.69000000e-03, 8.70000000e-03, 8.71000000e-03, 8.72000000e-03, 8.73000000e-03, 8.74000000e-03, 8.75000000e-03, 8.76000000e-03, 8.77000000e-03, 8.78000000e-03, 8.79000000e-03, 8.80000000e-03, 8.81000000e-03, 8.82000000e-03, 8.83000000e-03, 8.84000000e-03, 8.85000000e-03, 8.86000000e-03, 8.87000000e-03, 8.88000000e-03, 8.89000000e-03, 8.90000000e-03, 8.91000000e-03, 8.92000000e-03, 8.93000000e-03, 8.94000000e-03, 8.95000000e-03, 8.96000000e-03, 8.97000000e-03, 8.98000000e-03, 8.99000000e-03, 9.00000000e-03, 9.01000000e-03, 9.02000000e-03, 9.03000000e-03, 9.04000000e-03, 9.05000000e-03, 9.06000000e-03, 9.07000000e-03, 9.08000000e-03, 9.09000000e-03, 9.10000000e-03, 9.11000000e-03, 9.12000000e-03, 9.13000000e-03, 9.14000000e-03, 9.15000000e-03, 9.16000000e-03, 9.17000000e-03, 9.18000000e-03, 9.19000000e-03, 9.20000000e-03, 9.21000000e-03, 9.22000000e-03, 9.23000000e-03, 9.24000000e-03, 9.25000000e-03, 9.26000000e-03, 9.27000000e-03, 9.28000000e-03, 9.29000000e-03, 9.30000000e-03, 9.31000000e-03, 9.32000000e-03, 9.33000000e-03, 9.34000000e-03, 9.35000000e-03, 9.36000000e-03, 9.37000000e-03, 9.38000000e-03, 9.39000000e-03, 9.40000000e-03, 9.41000000e-03, 9.42000000e-03, 9.43000000e-03, 9.44000000e-03, 9.45000000e-03, 9.46000000e-03, 9.47000000e-03, 9.48000000e-03, 9.49000000e-03, 9.50000000e-03, 9.51000000e-03, 9.52000000e-03, 9.53000000e-03, 9.54000000e-03, 9.55000000e-03, 9.56000000e-03, 9.57000000e-03, 9.58000000e-03, 9.59000000e-03, 9.60000000e-03, 9.61000000e-03, 9.62000000e-03, 9.63000000e-03, 9.64000000e-03, 9.65000000e-03, 9.66000000e-03, 9.67000000e-03, 9.68000000e-03, 9.69000000e-03, 9.70000000e-03, 9.71000000e-03, 9.72000000e-03, 9.73000000e-03, 9.74000000e-03, 9.75000000e-03, 9.76000000e-03, 9.77000000e-03, 9.78000000e-03, 9.79000000e-03, 9.80000000e-03, 9.81000000e-03, 9.82000000e-03, 9.83000000e-03, 9.84000000e-03, 9.85000000e-03, 9.86000000e-03, 9.87000000e-03, 9.88000000e-03, 9.89000000e-03, 9.90000000e-03, 9.91000000e-03, 9.92000000e-03, 9.93000000e-03, 9.94000000e-03, 9.95000000e-03, 9.96000000e-03, 9.97000000e-03, 9.98000000e-03, 9.99000000e-03]), modes=[{'phase': 0.2, 'amp': 0.7, 'freq': 24000.0, 'mode_number': 3}, {'phase': 0.3, 'amp': 0.5, 'freq': 37000.0, 'mode_number': 4}], noise=0.2)

Generate synthetic multi-channel data for testing.

pyfusion.data.tests.get_n_channels(n_ch)

Return a list of n_ch channels.

timeseries Module

Timeseries data classes.

class pyfusion.data.timeseries.FlucStruc(svd_data, sv_list, timebase, min_dphase=-3.141592653589793, phase_pairs=None)

Bases: pyfusion.data.base.BaseData

fsplot_phase(input_data, *args, **kwargs)

plot the phase of a flucstruc, optionally inserting the first point at the end (if closed=True). Applies to closed arrays (e.g complete 2pi). Until Feb 2013, this version did not yet attempt to take into account angles, or check that adjacent channels are adjacent (i.e. ch2-ch1, ch2-c2 etc). Channel names are taken from the fs and plotted abbreviated

1/1/2011: TODO This appears to work only for database = None config 1/17/2011: bdb: May be fixed - I had used channel instead of channel.name

save()

svs()

class pyfusion.data.timeseries.SVDData(chrono_labels, topo_channels, svd_input)

Bases: pyfusion.data.base.BaseData

svd_input is a tuple as outputted by numpy.linalg.svd(data, 0)

plot_fs_groups(input_data, *args, **kwargs)

Show what would be generated by the supplied fs_grouper function.

An example of a grouper funcion is pyfusion.data.filters.fs_group_threshold

The grouper function takes an SVDData instance as its first input and returns an iterable where each element is a set of singular value indices which define a flucstruc.

self_cps()

svdplot(input_data, *args, **kwargs)

class pyfusion.data.timeseries.Signal

Bases: numpy.ndarray

Timeseries signal class with (not-yet-implemented) configurable digitisation.

see doc/subclassing.py in numpy code for details on subclassing ndarray

get_channel(channel_number, bounds=None)

allows us to use get_channel(0) no matter what ndim is

n_channels()

n_samples()

class pyfusion.data.timeseries.Timebase

Bases: numpy.ndarray

Timebase vector with parameterised internal representation.

see doc/subclassing.py in numpy code for details on subclassing ndarray

is_contiguous()

normalise_freq(input_freq)

Normalise input frequencies to [0,1] where 1 is omega = pi*sample_freq

class pyfusion.data.timeseries.TimeseriesData(timebase=None, signal=None, channels=None, bypass_length_check=False, **kwargs)

Bases: pyfusion.data.base.BaseData

copy(input_data, *args, **kwargs)

correlate(input_data, *args, **kwargs)

downsample(input_data, *args, **kwargs)

filter_fourier_bandpass(input_data, *args, **kwargs)

flucstruc(input_data, *args, **kwargs)

integrate(input_data, *args, **kwargs)

keys()

normalise(input_data, *args, **kwargs)

plot_signals(input_data, *args, **kwargs)

Plot a figure full of signals using n_columns[1],

sharey [=1] “gangs” y axes - sim for sharex - sharex=None stops this sharey: 2 gangs all but first (top) axis x axes are ganged by default: see Note:

fun, fun2: optionally plot a function of the signal. fun= refers to a function of one variable. fun2 is a function (t,x), such as one that returns a different timebase (diff should do this) if fun2 is given (a function of t and sig), then fun is ignored

labelfmt[“{short_name} {units}”] controls the channel labels.

The default ignores the shot and uses an abbreviated form of the channel name. If the short form is very short, it becomes a y label. A full version is “{name} {units}” and if > 8 chars, will become the x label. Even longer is “Shot={shot}, k_h={kh}, {name}”

labeleg: If ‘True’ put label in legend - else use this str as a legend lab linestyle: default of True means use ‘-‘ if no marker, and nothing if a

marker is given e.g. marker of ‘.’ and markersize<1 will produce “shaded” waveforms, good to see harmonic structure even without zooming in (need to adjust markersize or plot size for best results.

raw_names uses “digitiser” names, otherwise use names from the config file Note = sharex that to allow implicit overlay by using the same subplot specs the sharex must be the same between main and overlay - hence the use of explicit sharex = None

suptitle by default refers to the shot number

plot_spectrogram(input_data, *args, **kwargs)

title will be auto generated, if supplied, include ‘+’ to include the auto-generated part

reduce_time(input_data, *args, **kwargs)

remove_baseline(input_data, *args, **kwargs)

segment(input_data, *args, **kwargs)

sp_filter_butterworth_bandpass(input_data, *args, **kwargs)

subtract_mean(input_data, *args, **kwargs)

svd(input_data, *args, **kwargs)

pyfusion.data.timeseries.generate_timebase(t0=0.0, n_samples=10000.0, sample_freq=1000000.0)

pyfusion.data.timeseries.orm_load_flucstrucs(man)

pyfusion.data.timeseries.orm_load_svd_data(man)

pyfusion.data.timeseries.orm_load_timeseries_data(man)

utils Module

pyfusion.data.utils.bin2list(input_value)

pyfusion.data.utils.cps(a, b)

pyfusion.data.utils.find_peaks(arr, minratio=0.001, debug=0)

find the peaks in the data in arr, by selecting points See also Shauns using find_peaks running average where the slope changes sign, and the value is > minratio*max(arr)

pyfusion.data.utils.find_signal_spectral_peaks(timebase, signal, minratio=0.001, debug=0)

pyfusion.data.utils.get_axes_pixcells(ax)

return the pixcell coorindates of the axes ax This is useful in determining how many characters will fit

pyfusion.data.utils.list2bin(input_list)

pyfusion.data.utils.make_title(formatstr, input_data, channum=None, at_dict={}, min_length=3, raw_names=False)

Return a string describing the shot number, channel name etc using a formatstr which refers to items in a dictionary (at_dict), assembled in this routine, based on input_data and an optional dictionary which contains anything not otherwise available in input_data

pyfusion.data.utils.peak_freq(signal, timebase, minfreq=0, maxfreq=1e+18)

TODO: old code: needs review - since then bdb helped a bit... now also returns peaking factor this function only has a basic unittest to make sure it returns the correct freq in a simple case.

>>> tb = np.linspace(0,1,10000)
>>> int(peak_freq(np.sin(2*np.pi*567*tb), tb)[1])
567

pyfusion.data.utils.remap_periodic(input_array, min_val, period=6.283185307179586)

pyfusion.data.utils.split_names(names, pad=' ', min_length=3)

Given an array of strings, return an array of the part of the string (e.g. channel name) that varies, and optionally the prefix and suffix. The array of varying parts is first in the tuple in case others are not wanted. This is used to make the x labels of phase plots simpler and smaller. e.g.

>>> split_names(['MP01','MP10'],min_length=2)
(['01', '10'], 'MP', '')

The pad char is put on the end of shorter names - a better way would be to keep the end char the same, and pad in between the beginning and end the per channel part is at least min_length long. This is not really needed, as the routine chooses the lenght so that the results are not ambiguous (MP01,MP02 -> 1,2 but MP01,MP12 -> 01,12

pyfusion.data.utils.subdivide_interval(pts, overlap=None, debug=0)

return several intervals which straddle pts with overlap of ov The lowest x value is special - the point of division is much closer to that x then zero. overlap is a tuple (minoverlap, max), and describes the total overlap

pyfusion.data.utils.unique_id()

write_arff Module

for complex, had to unroll the loop to write a line of data - slows down about 7x!!

pyfusion.data.write_arff.debug_(debug, msg='', *args, **kwargs)

pyfusion.data.write_arff.split_vectors(dd, newfmts={})

pop any vectors and replace them with scalars consecutively numbered can supply a dictionary of formats for the new names - else use oldkey_0 etc Alters the dd to contain the items to be saved as scalars

pyfusion.data.write_arff.tomagphase(dd, key)

This would be happier as class method given a dictionary of arrays dd, and a key to complex data, add key and data for the magnitude and another for the phase, then remove the original key.

This destructive code has the advantage that the write is sped up by seven times.

pyfusion.data.write_arff.uniform(low=0.0, high=1.0, size=None)

Draw samples from a uniform distribution.

Samples are uniformly distributed over the half-open interval [low, high) (includes low, but excludes high). In other words, any value within the given interval is equally likely to be drawn by uniform.

Parameters:

• low (float, optional) – Lower boundary of the output interval. All values generated will be greater than or equal to low. The default value is 0.
• high (float) – Upper boundary of the output interval. All values generated will be less than high. The default value is 1.0.
• size (int or tuple of ints, optional) – Output shape. If the given shape is, e.g., (m, n, k), then m * n * k samples are drawn. Default is None, in which case a single value is returned.

Returns:

out – Drawn samples, with shape size.

Return type:

ndarray

See also

randint()

Discrete uniform distribution, yielding integers.

random_integers()

Discrete uniform distribution over the closed interval [low, high].

random_sample()

Floats uniformly distributed over [0, 1).

random()

Alias for random_sample.

rand()

Convenience function that accepts dimensions as input, e.g., rand(2,2) would generate a 2-by-2 array of floats, uniformly distributed over [0, 1).

Notes

The probability density function of the uniform distribution is

\[p(x) = \frac{1}{b - a}\]

anywhere within the interval [a, b), and zero elsewhere.

Examples

Draw samples from the distribution:

>>> s = np.random.uniform(-1,0,1000)

All values are within the given interval:

>>> np.all(s >= -1)
True
>>> np.all(s < 0)
True

Display the histogram of the samples, along with the probability density function:

>>> import matplotlib.pyplot as plt
>>> count, bins, ignored = plt.hist(s, 15, normed=True)
>>> plt.plot(bins, np.ones_like(bins), linewidth=2, color='r')
>>> plt.show()

pyfusion.data.write_arff.write_arff(da, filename='tmp.arff', use_keys=[])

use_keys - keys to save, empty list means all