binoculars package¶
Subpackages¶
- binoculars.backends package
- Submodules
- binoculars.backends.bm25 module
- binoculars.backends.example module
- binoculars.backends.id03 module
CylindricalQProjection
EH1
EH2
GammaDelta
GammaDeltaMu
GammaDeltaTheta
GisaxsDetector
HKLProjection
HKProjection
ID03Input
ID03Input.apply_mask()
ID03Input.dbg_pointno
ID03Input.dbg_scanno
ID03Input.find_edfs()
ID03Input.generate_jobs()
ID03Input.get_delayed_jobs()
ID03Input.get_delayed_scan()
ID03Input.get_destination_options()
ID03Input.get_images()
ID03Input.get_scan()
ID03Input.get_scan_params()
ID03Input.get_wavelength()
ID03Input.is_aborted()
ID03Input.is_zap()
ID03Input.parse_config()
ID03Input.process_job()
ID03Input.target()
ID03Input.wait_for_points()
QProjection
QTransformation
Qpp
SphericalQProjection
ThetaLProjection
TwoThetaProjection
load_matrix()
nrQProjection
pixels
specularangles
- binoculars.backends.id03_xu module
- binoculars.backends.sixs module
AnglesProjection
AnglesProjection2
DataFrame
Detector
Diffractometer
FLYMedVEiger
FlyMedH
FlyMedHS70
FlyMedV
FlyMedVS70
FlyScanUHV
FlyScanUHV2
FlyScanUHVS70
FlyScanUHVS70Andreazza
FlyScanUHVUfxc
GisaxUhvEiger
HKLProjection
HKProjection
M()
PDataFrame
Pixels
QIndex
QparQperIndexProjection
QparQperProjection
QxPolarProjection
QxQyIndexProjection
QxQyQzProjection
QxQzIndexProjection
QyPolarProjection
QyQzIndexProjection
QzPolarProjection
RealSpace
SBSFixedDetector
SBSMedH
SBSMedV
SBSMedVFixDetector
SIXS
Sample
Source
Stereo
SurfaceOrientation
dataframes()
get_detector()
get_diffractometer()
get_ki()
get_sample()
get_source()
hkl_matrix_to_numpy()
load_matrix()
normalized()
- Module contents
Submodules¶
binoculars.backend module¶
- class binoculars.backend.InputBase(config)[source]¶
Bases:
ConfigurableObject
Generate and process Job()s.
Note: there is no guarantee that generate_jobs() and process_jobs() will be called on the same instance, not even in the same process or on the same computer!
- class binoculars.backend.ProjectionBase(config)[source]¶
Bases:
ConfigurableObject
binoculars.dispatcher module¶
- class binoculars.dispatcher.Destination[source]¶
Bases:
object
- config = None¶
- filename = None¶
- limits = None¶
- opts: Dict[Any, Any] = {}¶
- overwrite = None¶
- type = None¶
- value = None¶
- class binoculars.dispatcher.DispatcherBase(config, main)[source]¶
Bases:
ConfigurableObject
- class binoculars.dispatcher.Local(config, main)[source]¶
Bases:
ReentrantBase
- actions: Tuple[str, ...] = ('user', 'job')¶
- class binoculars.dispatcher.Oar(config, main)[source]¶
Bases:
ReentrantBase
- actions: Tuple[str, ...] = ('user', 'process')¶
- class binoculars.dispatcher.ReentrantBase(config, main)[source]¶
Bases:
DispatcherBase
- actions: Tuple[str, ...] = ('user',)¶
binoculars.errors module¶
- exception binoculars.errors.BackendError[source]¶
Bases:
ExceptionBase
- exception binoculars.errors.CommunicationError[source]¶
Bases:
ExceptionBase
- exception binoculars.errors.ConfigError[source]¶
Bases:
ExceptionBase
- exception binoculars.errors.FileError[source]¶
Bases:
ExceptionBase
- exception binoculars.errors.SubprocessError[source]¶
Bases:
ExceptionBase
binoculars.fit module¶
- class binoculars.fit.FitBase(space, guess=None)[source]¶
Bases:
object
- fitdata = None¶
- guess = None¶
- parameters = None¶
- result = None¶
- summary = None¶
- class binoculars.fit.Gaussian1D(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.Lorentzian(space, guess=None)[source]¶
Bases:
AutoDimensionFit
- dimensions = {1: <class 'binoculars.fit.Lorentzian1D'>, 2: <class 'binoculars.fit.PolarLorentzian2D'>}¶
- class binoculars.fit.Lorentzian1D(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.Lorentzian1DNoBkg(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.Lorentzian2D(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.Lorentzian2Dnobkg(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.PolarLorentzian2D(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.PolarLorentzian2Dnobkg(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
- class binoculars.fit.Voigt1D(space, guess=None, loc=None)[source]¶
Bases:
PeakFitBase
binoculars.main module¶
binoculars.plot module¶
binoculars.space module¶
- class binoculars.space.Axes(axes)[source]¶
Bases:
object
Luxurious tuple of Axis objects.
- property dimension¶
- property memory_size¶
- property npoints¶
- class binoculars.space.Axis(min, max, res, label=None)[source]¶
Bases:
object
Represents a single dimension finite discrete grid centered at 0.
Important attributes: min lower bound max upper bound res step size / resolution label human-readable identifier
min, max and res are floats, but internally only integer operations are used. In particular min = imin * res, max = imax * res
- property max¶
- property min¶
- class binoculars.space.EmptySpace(config=None, metadata=None)[source]¶
Bases:
object
Convenience object for sum() and friends. Treated as zero for addition. Does not share a base class with Space for simplicity.
- class binoculars.space.EmptyVerse[source]¶
Bases:
object
Convenience object for sum() and friends. Treated as zero for addition.
- property dimension¶
- spaces = [EmptySpace]¶
- class binoculars.space.Multiverse(spaces)[source]¶
Bases:
object
A collection of spaces with basic support for addition. Only to be used when processing data. This makes it possible to process multiple limit sets in a combination of scans
- property dimension¶
- class binoculars.space.Space(axes, config=None, metadata=None)[source]¶
Bases:
object
Main data-storing object in BINoculars. Data is represented on an n-dimensional rectangular grid. Per grid point, the number of photons (~ intensity) and the number of original data points (pixels) contribution is stored.
- Important attributes:
axes Axes instances describing range and stepsizes of each of the dimensions photons n-dimension numpy float array, total intensity per grid point contribitions n-dimensional numpy integer array, number of original datapoints (pixels) per grid point dimension n
- property config¶
util.ConfigFile instance describing configuration file used to create this Space instance
- property dimension¶
- classmethod from_image(resolutions, labels, coordinates, intensity, weights, limits=None)[source]¶
Create Space from image data.
resolutions n-tuple of axis resolutions labels n-tuple of axis labels coordinates n-tuple of data coordinate arrays intensity data intensity array
- classmethod fromfile(file: str, key: Optional[Sequence[slice]] = None) Space [source]¶
Load Space from HDF5 file.
file filename string or h5py.Group instance key sliced (subset) loading, should be an n-tuple of slice()s in data coordinates
- get_grid() Tuple[ndarray] [source]¶
Returns the data coordinates of each grid point, as n-tuple of n-dimensinonal arrays. Basically numpy.mgrid() in data coordinates.
- get_key(key: Union[float, slice, Tuple[Union[float, slice, ndarray], ...], List[Union[float, slice, ndarray]]]) Tuple[Union[int, slice, ndarray], ...] [source]¶
Convert the n-dimensional interval described by key (as used by e.g. __getitem__()) from data coordinates to indices.
- get_masked() MaskedArray [source]¶
Returns photons/contributions, but with divide-by-zero’s masked out.
- property memory_size¶
Returns approximate memory consumption of this Space. Only considers size of .photons and .contributions, does not take into account the overhead.
- property metadata¶
util.MetaData instance describing metadata used to create this Space instance
- property npoints¶
- process_image(coordinates, intensity, weights)[source]¶
Load image data into Space.
coordinates n-tuple of data coordinate arrays intensity data intensity array weights weights array, supply numpy.ones_like(intensity) for equal weights
- project(axis: Union[str, int], *more_axes) Space [source]¶
Reduce dimensionality of Space by projecting onto ‘axis’. All data (photons, contributions) is summed along this axis.
axis the label of the axis or the index *more_axis also project on these axes
- rebin(resolutions)[source]¶
Change bin size.
resolution n-tuple of floats, new resolution of each axis
- binoculars.space.chunked_sum(verses: Iterable[Multiverse], chunksize=10) Union[EmptyVerse, Multiverse] [source]¶
Calculate sum of iterable of Multiverse instances. Creates intermediate sums to avoid growing a large space at every summation.
verses iterable of Multiverse instances chunksize number of Multiverse instances in each intermediate sum
- binoculars.space.silence_numpy_errors()[source]¶
Silence numpy warnings about zero division. Normal usage of Space() will trigger these warnings.
- binoculars.space.sum_onto(a, axis)[source]¶
Numpy convenience. Project all dimensions of an array onto an axis, i.e. apply sum() to all axes except the one given.
- binoculars.space.verse_sum(verses: Generator[Multiverse, None, None]) Multiverse [source]¶
binoculars.util module¶
- class binoculars.util.OrderedOperation(option_strings, dest, nargs=None, const=None, default=None, type=None, choices=None, required=False, help=None, metavar=None)[source]¶
Bases:
Action
- binoculars.util.atomic_write(filename)[source]¶
Atomically write data into ‘filename’ using a temporary file and os.rename()
Rename on success, clean up on failure (any exception).
Example: with atomic_write(filename) as tmpfile
- with open(tmpfile, ‘w’) as fp:
fp.write(…)
- binoculars.util.chunk_slicer(count, chunksize)[source]¶
yields slice() objects that split an array of length ‘count’ into equal sized chunks of at most ‘chunksize’
- binoculars.util.cluster_jobs2(jobs, target_weight)[source]¶
Taking the first n jobs that together add up to target_weight. Here as opposed to cluster_jobs the total number of jobs does not have to be known beforehand
- binoculars.util.loop_delayer(delay)[source]¶
Delay a loop such that it runs at most once every ‘delay’ seconds. Usage example: delay = loop_delayer(5) while some_condition:
next(delay) do_other_tasks
- binoculars.util.status(line, eol=False)[source]¶
Prints a status line to sys.stdout, overwriting the previous one. Set eol to True to append a newline to the end of the line
- binoculars.util.wait_for_files(filelist, timeout=None)[source]¶
Wait until the files in ‘filelist’ have appeared, for a maximum of ‘timeout’ seconds. Returns True on success, False on timeout.
- binoculars.util.yield_when_exists(filelist, timeout=None)[source]¶
Wait for files in ‘filelist’ to appear, for a maximum of ‘timeout’ seconds, yielding them in arbitrary order as soon as they appear. If ‘filelist’ is a set, it will be modified in place, and on timeout it will contain the files that have not appeared yet.
Module contents¶
- binoculars.fitspace(space, function, guess=None)[source]¶
fit the space data.
Parameters space: binoculars space function: list
a string with the name of the desired function. supported are: lorentzian (automatically selects 1d or 2d), gaussian1d and voigt1d
- guess: list
a list of length N with the resolution per label
Returns A binoculars fit object.
Examples: >>> fit = binoculars.fitspace(space, ‘lorentzian’) >>> print(fit.summary)
I: 1.081e-07 +/- inf loc: 0.3703 +/- inf gamma: 0.02383 +/- inf slope: 0.004559 +/- inf offset: -0.001888 +/- inf
>>> parameters = fit.parameters >>> data = fit.fitdata >>> binoculars.plotspace(space, fit = data)
- binoculars.info(filename)[source]¶
Explore the file without loading the file, or after loading the file
Parameters filename: filename or space
Examples: >>> print binoculars.info(‘test.hdf5’) Axes (3 dimensions, 46466628 points, 531.0 MB) {
Axis H (min=-0.1184, max=0.0632, res=0.0008, count=228) Axis K (min=-1.1184, max=-0.9136, res=0.0008, count=257) Axis L (min=0.125, max=4.085, res=0.005, count=793)
} ConfigFile{
[dispatcher] [projection] [input]
} origin = test.hdf5 >>> space = binoculars.load(‘test.hdf5’) >>> print binoculars.info(space) Axes (3 dimensions, 46466628 points, 531.0 MB) {
Axis H (min=-0.1184, max=0.0632, res=0.0008, count=228) Axis K (min=-1.1184, max=-0.9136, res=0.0008, count=257) Axis L (min=0.125, max=4.085, res=0.005, count=793)
} ConfigFile{
[dispatcher] [projection] [input]
} origin = test.hdf5
- binoculars.load(filename, key=None)[source]¶
Parameters filename: string
Only hdf5 files are acceptable
key: a tuple with slices in as much dimensions as the space is
Returns A binoculars space
Examples: >>> space = binoculars.load(‘test.hdf5’) >>> space Axes (3 dimensions, 2848 points, 33.0 kB) {
Axis qx (min=-0.01, max=0.0, res=0.01, count=2) Axis qy (min=-0.04, max=-0.01, res=0.01, count=4) Axis qz (min=0.48, max=4.03, res=0.01, count=356)
}
- binoculars.plotspace(space, log=True, clipping=0.0, fit=None, norm=None, colorbar=True, labels=True, **plotopts)[source]¶
plots a space with the correct axes. The space can be either one or two dimensinal.
Parameters space: binoculars space
the space containing the data that needs to be plotted
- log: bool
axes or colorscale logarithmic
- clipping: 0 < float < 1
cuts a lowest and highst value on the color scale
fit: numpy.ndarray
same shape and the space. If one dimensional the fit will be overlayed.
- norm: matplotlib.colors
object defining the colorscale
- colorbar: bool
show or not show the colorbar
- labels: bool
show or not show the labels
- plotopts: keyword arguments
keywords that will be accepted by matplotlib.pyplot.plot or matplotlib.pyplot.imshow
Examples: >>> space Axes (3 dimensions, 2848 points, 33.0 kB) {
Axis qx (min=-0.01, max=0.0, res=0.01, count=2) Axis qy (min=-0.04, max=-0.01, res=0.01, count=4) Axis qz (min=0.48, max=4.03, res=0.01, count=356)
} >>> binoculars.plotspace(‘test.hdf5’)
- binoculars.run(args)[source]¶
Parameters args: string
String as if typed in terminal. The string must consist of the location of the configuration file and the command for specifying the jobs that need to be processed. All additonal configuration file overides can be included
Returns A tuple of binoculars spaces
Examples: >>> space = binoculars.run(‘config.txt 10’) >>> space[0] Axes (3 dimensions, 2848 points, 33.0 kB) {
Axis qx (min=-0.01, max=0.0, res=0.01, count=2) Axis qy (min=-0.04, max=-0.01, res=0.01, count=4) Axis qz (min=0.48, max=4.03, res=0.01, count=356)
}
- binoculars.save(filename, space)[source]¶
Save a space to file
Parameters filename: string
filename to which the data is saved. ‘.txt’, ‘.hdf5’ are supported.
- space: binoculars space
the space containing the data that needs to be saved
Examples: >>> space Axes (3 dimensions, 2848 points, 33.0 kB) {
Axis qx (min=-0.01, max=0.0, res=0.01, count=2) Axis qy (min=-0.04, max=-0.01, res=0.01, count=4) Axis qz (min=0.48, max=4.03, res=0.01, count=356)
} >>> binoculars.save(‘test.hdf5’, space)
- binoculars.transform(space, labels, resolutions, exprs)[source]¶
transformation of the coordinates.
Parameters space: binoculars space labels: list
a list of length N with the labels
- resolutions: list
a list of length N with the resolution per label
- exprs: list
a list of length N with strings containing the expressions that will be evaluated. all numpy funtions can be called without adding ‘numpy.’ to the functions.
Returns A binoculars space of dimension N with labels and resolutions specified in the input
Examples: >>> space = binoculars.load(‘test.hdf5’) >>> space Axes (3 dimensions, 2848 points, 33.0 kB) {
Axis qx (min=-0.01, max=0.0, res=0.01, count=2) Axis qy (min=-0.04, max=-0.01, res=0.01, count=4) Axis qz (min=0.48, max=4.03, res=0.01, count=356)
} >>> newspace = binoculars.transform(space, [‘twotheta’], [0.003], [‘2 * arcsin(0.51 * (sqrt(qx**2 + qy**2 + qz**2) / (4 * pi)) / (pi * 180))’]) # noqa >>> newspace Axes (1 dimensions, 152 points, 1.0 kB) {
Axis twotheta (min=0.066, max=0.519, res=0.003, count=152)
}