5. smallanglescattering (sas)¶
This module allows smearing/desmearing of SAXS/SANS data.
Smearing is for line collimation as in Kratky SAXS cameras and for point collimation for SANS. For SANS the resolution smearing a la Pedersen is realised (resFunct).
For desmearing the Lake algorithm is an iterative procedure to desmear smeared data. We follow here the improvements according to Vad using a convergence criterion and smoothing.
2D SAXS dat can be read into sasImage, background subtracted , transmission corrected….
As references the waterXray scattering and a AgBe reference spectrum are available.
For form factors and structure factors see the respective modules.
5.1. Smear/desmear 1D¶
smear (data, beamProfile, **kwargs) |
Smearing data for line-collimated SAXS (Kratky camera) or as point collimation SANS/SAXS. |
desmear (Ios, beamProfile[, NIterations, …]) |
Desmearing according to Lake algorithm with possibility to stop recursion at best desmearing. |
resFunct (S[, collDist, collAperture, …]) |
Resolution smearing of small angle scattering for SANS or SAXS according to Pedersen for radial averaged data. |
prepareBeamProfile ([data]) |
Prepare beam profile from Beam Profile measurement or according to given parameters. |
getBeamWidth (empty[, minmax, show]) |
Extract primary beam of empty cell or buffer measurement for semitransparent beam stops. |
plotBeamProfile (beam[, p]) |
Plots beam profile and weight function according to parameters in beam. |
transmissionCorrection (data, dark[, …]) |
Subtract dark current, find primary beam, get transmission and normalize by transmission and exposure time. |
waterXrayScattering ([composition, T, units]) |
Absolute scattering of water with components (salt, buffer) at Q=0 as reference for X-ray. |
AgBeReference (q, wavelength[, n, ampn, …]) |
The scattering intensity expected from AgBe as a reference for wavelength calibration. |
5.2. 2D sasImage¶
Read 2D image files (TIFF) from SAXS cameras and extract the corresponding data.
The sasImage is a 2D array that allows direct subtraction and multiplication (e.g. transmission) respecting given masks in operations. E.g.
sample=js.sas.sasImage('sample.tiff')
solvent=js.sas.sasImage('solvent.tiff')
corrected = sample/sampletransmission - solvent/solventtransmission
Calibration of detector distance, radial average, size reduction and more. .pickBeamcenter allows sensitive detection of the beamcenter.
An example is shown in sasImage
.
sasImage |
Creates sasImage as maskedArray from a detector image for evaluation. |
5.2.1. sasImage methods¶
maskReset
()Reset the mask. maskRegion
(xmin, xmax, ymin, ymax)Mask rectangular region. maskFromImage
(image)Use/copy mask from image. maskRegion
(xmin, xmax, ymin, ymax)Mask rectangular region. maskRegions
(regions)Mask several regions. maskbelowLine
(p1, p2)Mask points at one side of line. maskTriangle
(p1, p2, p3[, invert])Mask inside triangle. mask4Polygon
(p1, p2, p3, p4[, invert])Mask inside polygon of 4 points. maskCircle
(center, radius[, invert])Mask points inside circle. maskSectors
(angles, width[, radialmax, …])Mask sector around beamcenter. asImage
([scale, colormap, inverse, …])Returns the sasImage as 8bit RGB image using PIL. saveAsTIF
(filename[, fill])Save the sasImage as float32 tif without loosing information. pickBeamcenter
([levels, symmetry])Open image to pick the beamcenter from a calibration sample as AgBe. setDetectorDistance
(detector_distance[, offset])Set detector distance from calibration . findCenterOfIntensity
([beamcenter, size])Find beam center as center of intensity around beamcenter. radialAverage
([beamcenter, number, kind])Radial average of image and conversion to wavevector q. lineAverage
([beamcenter, number, minmax, show])Line average of image and conversion to wavevector q for line collimation cameras. recalibrateDetDistance
([beamcenter, number, …])Recalibration of detectorDistance by AgBe reference for point collimation. gaussianFilter
([sigma])Gaussian filter in place. showPolar
([beamcenter, scaleR, offset])Show image transformed to polar coordinates around beamcenter. reduceSize
([pixelsize, center, border])Reduce size of image using uniform average in box. show
(**kwargs)Show sasImage as matplotlib figure. getfromcomment
(name[, replace, newname])Extract name from .artist or .imageDescription with attribute name in front. array
Strip of all attributes and return a simple array without mask. asdataArray
([masked])Return representation of sasImage as dataArray representing wavevectors (qx,qy) against intensity. interpolateMaskedRadial
([radial])Interpolate masked values from radial averaged image or function.
5.3. 2D sasImage convenience¶
createLogPNG (filenames[, center, size, …]) |
Create .png files from grayscale images with log scale conversion to values between [1,255]. |
createImageDescriptions (images) |
Create text file with image descriptions as list of content. |
readImages (filenames) |
Read a list of images returning sasImage`s. |
5.4. Housekeeping¶
readpdh (pdhFileName) |
Opens and reads a SAXS data file in the .pdh (Primary Data Handling) format. |
autoscaleYinoverlapX (dataa[, key, keep]) |
Scales elements of data to have same mean .Y value in the overlap region of .X . |
removeSpikesMinmaxMethod (dataa[, order, …]) |
Takes a dataset and removes single spikes from data by substitution with spline. |
removeSpikes (dataa[, xmin, xmax, medwindow, …]) |
Takes a dataset and removes single spikes. |
locateFiles (pattern[, root]) |
Locate all files matching supplied filename pattern in and below supplied root directory. |
copyFiles (pattern[, root, destination, link]) |
Copies all files matching pattern in tree below root to destination directory |
addXMLParameter (data) |
Adds the parameters stored in xml part of a .pdh file as eg. |
moveSAXSPACE (pattern[, root, destination, …]) |
Read SAXSPACE .pdh files and removes spikes by removeSpikes. |
This module allows smearing/desmearing of SAXS/SANS data.
Smearing is for line collimation as in Kratky SAXS cameras and for point collimation for SANS. For SANS the resolution smearing a la Pedersen is realised (resFunct).
For desmearing the Lake algorithm is an iterative procedure to desmear smeared data. We follow here the improvements according to Vad using a convergence criterion and smoothing.
2D SAXS dat can be read into sasImage, background subtracted , transmission corrected….
As references the waterXray scattering and a AgBe reference spectrum are available.
For form factors and structure factors see the respective modules.
-
jscatter.smallanglescattering.
AgBeReference
(q, wavelength, n=array([1, 2, 3, 4, 5, 6, 7, 8, 9]), ampn=None, domainsize=100, udw=0.1, asym=0, lg=1)[source]¶ The scattering intensity expected from AgBe as a reference for wavelength calibration.
The intensities assume a d-spacing of 5.8378 nm and a reduction of the intensity as q**-2. The domain size determines the width according to Scherrer equation [2]. The first peak is at 1.076 1/nm. The result needs to be convoluted with the instrument resolution by resFunct or smear.
Parameters: - q : array
Wavevector
- wavelength : float
Wavelength
- n : array of int
Order of the peaks.
- ampn : list of float
Amplitudes of the peaks
- domainsize : float
Domainsize of AgBe crystals in nm. default 100 nm as is given in [1].
- udw : float
Displacement u in Debye Waller factor exp(-u**2*q**2/3)
- asym : float
Factor asymmetry in Voigt function describing the peaks.
- lg : float
Lorenzian/gaussian fraction of both FWHM, describes the contributions of gaussian and lorenzian shape. See Voigt for details.
Returns: - dataArray
References
[1] (1, 2) T. C. Huang, H. Toraya, T. N. Blanton and Y. Wu X-ray Powder Diffraction Analysis of Silver Behenate, a Possible Low-Angle Diffraction Standard J. Appl.Cryst.(1993).26,180-184 [2] (1, 2) Patterson, A. The Scherrer Formula for X-Ray Particle Size Determination Phys. Rev. 56 (10): 978–982 (1939) doi:10.1103/PhysRev.56.978.
-
jscatter.smallanglescattering.
addXMLParameter
(data)[source]¶ Adds the parameters stored in xml part of a .pdh file as eg. in SAXSPACE .pdh files.
Parameters: - data : dataArray
Already read pdh file. XML content is found in comments of the read files and starts with ‘<’.
-
jscatter.smallanglescattering.
autoscaleYinoverlapX
(dataa, key=None, keep='lowest')[source]¶ Scales elements of data to have same mean .Y value in the overlap region of .X .
Parameters: - dataa : dataList
data to scale
- key : string
Data are grouped into unique values of attribute key before scaling. E.g. to do it for a series of concentrations for each concentration individually.
- keep : default ‘l’
If ‘l’ the lowest X are kept and higher X are scaled successively to next lower X. Anything else highest X are kept and other are scaled to next higher.
Returns: - dataList
new scaled dataList
Notes
First data are sorted along .X.mean() scaling value is stored in .autoscalefactor
-
jscatter.smallanglescattering.
copyFiles
(pattern, root='.', destination='copy', link=False)[source]¶ Copies all files matching pattern in tree below root to destination directory
Parameters: - pattern : file pattern
Pattern used in fnmatch.filter
- root : directory, default is os.curdir
Directory where to start
- destination : dirname
Destination
- link : bool
If True links are created.
-
jscatter.smallanglescattering.
desmear
(Ios, beamProfile, NIterations=-15, windowsize=4, qmax=4, output=True, **kwargs)[source]¶ Desmearing according to Lake algorithm with possibility to stop recursion at best desmearing.
For negative NIterations the iterations are stopped if a convergence criterion reaches a minimum as described by Vad [2]. In each step a smoothing based on the ratio desmeared/observed as described in [2] is used (point average with windowsize).
Parameters: - Ios : dataArray
Original smeared data
- beamProfile : dataArray
Beam profile as prepared with prepareBeamProfile
- NIterations : int, default=-15
Number of iterations to stop. Negative values indicate to stop when convergence criterion increases again (as described by Vad [2]) and abs(NIterations) is maximum number of iterations.
- qmax : float, default=4
Maximum in scattering vector q up to where the convergence criterion is evaluated. This reduces the influence of the noise at larger a.
- windowsize : odd int , default=4
Window size for smoothing in each step of desmearing (running average).
- output : bool
Print output.
Returns: - dataArray
References
[1] Lake, J. A. (1967). Acta Cryst. 23, 191–194. [2] (1, 2, 3, 4) Comparison of iterative desmearing procedures for one-dimensional small-angle scattering data Vad and Sager, J. Appl. Cryst. (2011) 44,32-42
-
jscatter.smallanglescattering.
fitDarkCurrent
(darkfile)[source]¶ Fits dark current with 5th order polynom + cosine.
This is dangerous as the darkcurrent has a noise that is not random but depends on used detector and count time.
-
jscatter.smallanglescattering.
getBeamWidth
(empty, minmax='auto', show=False)[source]¶ Extract primary beam of empty cell or buffer measurement for semitransparent beam stops.
The primary beam is searched and cut between the next minima found, then normalized. Additionally a Gaussian fit is done and hwhm is included in result profile.
Parameters: - empty : dataArray
Empty cell measurement with the transmitted beam included.
- minmax : ‘auto’,[float,float]
Automatic or interval for search of primary beam. E.g. [-0.03,0.03] allow for explicitly setting the interval.
- show : bool
Show the fit result
Returns: - dataArray with beam width profile
.sigma sigma of fit with Gaussian .hwhm half width half maximum
-
jscatter.smallanglescattering.
locateFiles
(pattern, root='.')[source]¶ Locate all files matching supplied filename pattern in and below supplied root directory.
Parameters: - pattern : file pattern
Pattern used in fnmatch.filter
- root : directory, default is os.curdir
Directory where to start
Returns: - file list
-
jscatter.smallanglescattering.
moveSAXSPACE
(pattern, root='./', destination='./despiked', medwindow=5, SGwindow=5, sigma=0.2, order=2)[source]¶ Read SAXSPACE .pdh files and removes spikes by removeSpikes.
This is mainly for use at JCNS SAXSPACE with CCD camera as detector :-))))
Parameters: - pattern : string
Search pattern for filenames
- root : string
Root path
- destination : string
Where to save the files
- medwindow : odd integer
Window size of scipy.signal.medfilt
- SGwindow : odd int, None
Savitzky-Golay filter see scipy.signal.savgol_filter
- order : int
Polynominal order of scipy.signal.savgol_filter
- sigma : float
Deviation factor of eY If datapoint-median> sigma*std its a spike
Notes
Default values are adjusted to typical SAXSPACE measurement.
-
jscatter.smallanglescattering.
plotBeamProfile
(beam, p=None)[source]¶ Plots beam profile and weight function according to parameters in beam.
Parameters: - beam
beam with parameters
- p : GracePlot instance
Reuse the given plot
-
jscatter.smallanglescattering.
prepareBeamProfile
(data=None, **kwargs)[source]¶ Prepare beam profile from Beam Profile measurement or according to given parameters.
Parameters: - data : dataArray,’trapez’,’SANS’
- Contains measured beam profile, explicit Gaussian width list or type ‘SANS’, ‘trapz’.
- dataArray Line collimation as measured can be given and will be smoothed and made symmetric.
- dataArray with explicit given Gaussian width for each Q values, missing values will be interpolated.
- ‘trapez’ : Line collimation with trapezoidal parameters a, b, bxw, dIW.
- ‘SANS’ : Smearing a la Pedersen; see resFunct for parameters
- collDist,collAperture,detDist,sampleAperture : float
Parameters as described in resFunct for SANS These are determined from the experimental setup.
- wavelength,wavespread,dpixelWidth,dringwidth,extrapolfunc : float
Parameters as described in resFunct for SANS These are determined from the experimental setup.
- a : float
Larger full length of trapezoidal profile in detector q units. Ignored for measured profile.
- b : float
Shorter full length of trapezoidal profile in detector q units If a=b ==> a=a*(1+1e-7), b=b*(1-1e-7) Ignored for measured profile.
- bxw : float,dataArray
Beam width profile. Use getBeamWidth to cut the primary beam and fit a Gaussian. A float describes the beam half-width at half maximum (hwhm of Gaussian). If bxw is the profile prepared by getBeamWidth the experimental profile is used.
- dIW : float
Detector slit width in detector q units. Length on detector to integrate parallel to beam length for line collimation. On my SAXSpace this is 1.332 as given in the header of the file (20 mm in real coordinates).
- wavelength : float,
Wavelength in nm default 0.155418 nm for SAXS 0.5 nm for SANS
- detDist : float, default 305.3558
Detector distance in units mm Default is detector distance of SAXSpace
- explicit : int
For explicit given Gaussian width the index of the column with the width. For merged dataFiles of KWS2@MLZ this is the forth column with index 3.
Returns: - beam profile as dataArray
Notes
- For measured beam profiles parameters a,b are determined from the flanks for trapezoidal profile.
- Detector q units are equivalent to the pixel distance as expressed in a corrected measurement.
- For ‘explicit’ Gaussian width a SANS measurement as on KWS2 can be used which has sigma as 4th column. Missing values are interpolated.
Examples
# use as # prepare measured beamprofile mbeam = js.sas.prepareBeamProfile(beam, bxw=0.01, dIW=1.) # prepare profile with trapezoidal shape (a,b are fitted above) tbeam = js.sas.prepareBeamProfile('trapz', a=mbeam.a, b=mbeam.b, bxw=0.01, dIW=1) # prepare profile SANS (missing parameters get defaults) Sbeam = js.sas.prepareBeamProfile('SANS', detDist=2000,wavelength=0.4,wavespread=0.15) # prepare profile with explicit given Gaussian width in column 3 as e.g. KWS2@MLZ Gbeam = js.sas.prepareBeamProfile(measurement,explicit=3)
-
jscatter.smallanglescattering.
readpdh
(pdhFileName)[source]¶ Opens and reads a SAXS data file in the .pdh (Primary Data Handling) format.
If data contain X values <0 it is assumed that the primary beam is included as for SAXSpace instruments. Use js.sas.transmissionCorrection to subtract dark counts and do transmission correction.
Parameters: - pdhFileName : string
file name
Returns: - dataArray
Notes
Alternatively the files can be read ignoring the information in the header (it is stored in the comments)
data=dA(pdhFileName,lines2parameter=[2,3,4])
PDH format used in the PCG SAXS software suite developed by the Glatter group at the University of Graz, Austria. This format is described in the appendix of the PCG manual (below from version 4.05.12 page 159).
In the PDH format, lines 1-5 contain header information, followed by the SAXS data.
line 1: format A80 -> description line 2: format 16(A4,1X) -> description in 16x4 character groups (1X = space separated) line 3: format 8(I9,1X) -> 8 integers (1X = space separated) line 4: format 5(E14.6,1X) -> 5 float (1X = space separated) line 5: format 5(E14.6,1X -> 5 float(1X = space separated) line 6+ format 3(E14.6,1X) - SAXS data x,y,error (1X = space separated)
- with:
- line 3 field 0 : number of points
- line 4 field 4 : normalization constant (default 1, never zero!!)
Anything else can have different meanings.
- The SAXSpace and SAXSess (AntonPaar) format add:
- line 4 field 2 : detector distance in mm
- line 4 field 5 : wavelength
- line 5 field 2 : detector slit length (equivalent to width of integration area) in q units
Additional xml parameter as in the SAXSPACE format appended can be extracted to attributes by addXMLParameter. Mainly this is “Exposure” time.
Example data for SAXSpace
<emptyline> SAXS BOX 2057 0 0 0 0 0 0 0 0.000000E+00 3.052516E+02 0.000000E+00 1.000000E+00 1.541800E-01 0.000000E+00 1.332843E+00 0.000000E+00 0.000000E+00 0.000000E+00 -1.033335E-01 2.241656E+03 1.024389E+00 -1.001430E-01 2.199972E+03 1.052537E+00 ...
-
jscatter.smallanglescattering.
removeSpikes
(dataa, xmin=None, xmax=None, medwindow=5, SGwindow=None, sigma=0.2, SGorder=2)[source]¶ Takes a dataset and removes single spikes.
A median filter is used to find single spikes. If abs(data.Y-medianY)/data.eY>sigma then the medianY value is used. If SGwindow!=None Savitzky-Golay filtered values are used. If sigma is 0 then new values (median or Savitzky-Golay filtered) are used everywhere.
Parameters: - xmin,xmax : float
Minimum and maximum X values
- dataa : dataArray
dataset with eY data
- medwindow : odd integer
window size of scipy.signal.medfilt
- SGwindow : odd int, None
Savitzky-Golay filter see scipy.signal.savgol_filter without the spikes; window should be smaller than instrument resolution
- SGorder : int
Polynomial order of scipy.signal.savgol_filter needs to be smaller than SGwindow
- sigma : float
relative deviation from eY if datapoint-median> sigma*eY its a spike
Returns: - Filtered and smoothed dataArray
-
jscatter.smallanglescattering.
removeSpikesMinmaxMethod
(dataa, order=7, sigma=2, nrepeat=1, removePoints=None)[source]¶ Takes a dataset and removes single spikes from data by substitution with spline.
Find minima and maxima of data including double point spikes; no 3 point spikes scipy.signal.argrelextrema with np.greater and np.less are used to find extrema
Parameters: - dataa : dataArray
Dataset with eY data.
- order : int
Number of points see scipy.signal.argrelextrema. Distance between extrema.
- sigma : float
Deviation factor from std dev; from eY. If datapoint -spline> sigma*std its a spike.
- nrepeat : int
Repeat the procedure nrepeat times.
- removePoints : list of integer
Instrument related points to remove because of dead pixels. ‘JCNS’ results in a list for SAXSPACE at JCNS Jülich.
Returns: - data with spikes removed
-
jscatter.smallanglescattering.
resFunct
(S, collDist=8000.0, collAperture=10, detDist=8000.0, sampleAperture=10, wavelength=0.5, wavespread=0.2, dpixelWidth=10, dringwidth=1, extrapolfunc=None)[source]¶ Resolution smearing of small angle scattering for SANS or SAXS according to Pedersen for radial averaged data.
I(q0)= Integral{(R(q,q0)*S(q)}dq with Kernel R(q,q0) of equ. 33 in [1] including wavelength spread, finite collimation and detector resolution. Default parameters are typical for a SANS machine like KWS2@JCNS with rectangular apertures. Low Q can be extrapolated as power law or Guinier like or constant. Best practice is to calculate the used model for large Q range, smear it and prune to the needed data range. This is demonstrated in example 2.
Parameters: - S : array like
dataArray with X. and .Y theoretical Scattering function q in nm^-1 .Y can be arbitrary unit
- collDist : float, default 8000
collimation distance in mm
- collAperture : float, default 10
collimation rectangular aperture size in mm
- detDist : float, default 8000
detector distance in mm
- sampleAperture : float, default 10
sample rectangular aperture size in mm
- wavelength : float, default 0.5
wavelength in nm
- wavespread : float, default 0.1
FWHM wavelengthspread dlambda/lambda
- dpixelWidth : float, default 10
Detector pixel width in mm
- dringwidth : integer, default 1
number of pixel for averaging
- extrapolfunc : list , default None
Type of extrapolation at low and high X edges as list for [low edge,high edge]. If a singe value is given this is used for both.
- Low edge towards zero
- float : Power law extrapolation
- None : A constant value as Y(X.min()).
- anything else : Low X data are log scaled, then X**2 extrapolated as Guinier like extrapolation.
- High edge towards infinity
- float : Power law extrapolation of low X e.g. a=-4 for q**a for Porod scaling.
- None : A constant value as Y(X.max()).
Returns: - dataArray
columns [‘wavevector; smeared scattering; unsmeared scattering; half width smearing function’]
Notes
- HalfWidthSmearingFunction is the FWHM the Gaussian used for smearing including all effects.
- The resolution is assumed to be equal in direction parallel and perpendicular to q on a 2D detector as described in chap. 2.5 in [1].
- We neglect additional smearing due to radial averaging (last paragraph in chap 2.5 of [1]).
- Defaults correspond to a typical medium resolution measurement.
- extrapolfunc allows extrapolation at both edges to reduce edge effects. The best values depends on the measured signal shape at the edge. The optimal way is to calculate the used model for the whole Q range, smear it and prune to the needed range. This is demonstrated in example 2. If at low q only a power law is observed use the corresponding extrapolation. Same for high q.
- An example for SANS fitting with resFunc is given in example_SANSsmearing.py.
References
[1] (1, 2, 3, 4, 5) Analytical Treatment of the Resolution Function for Small-Angle Scattering JAN SKOV PEDERSEN, DORTHE POSSELTAND KELL MORTENSEN J. Appl. Cryst. (1990). 23, 321-333 Examples
Reproducing Table 1 of [1]
import jscatter as js q=js.loglist(0.1,10,500) S=js.ff.sphere(q,6) # this is the direct call of resFunc, use smear instead as shown in next example Sr=js.sas.resFunct(S, collDist=2000.0, collAperture=20, detDist=2000.0, sampleAperture=10, wavelength=0.5, wavespread=0.2,dpixelWidth=0,dringwidth=0) # plot it p=js.grace() p.plot( S,sy=[1,0.3],li=1,legend='sphere') p.plot( Sr,sy=[2,0.3,2],li=2,legend='smeared sphere') p.plot(Sr.X,Sr[-1],li=4,sy=0,legend='FWHM in nm\S-1 ') p.yaxis(min=1e-3,scale='l',charsize=1.5,label='I(q) / a.u.',tick=[10,9]) p.yaxis(min=1e-1,tick=[10,9]) p.xaxis(scale='l',charsize=1.5,label='q / nm\S-1',tick=[10,9]) p.legend(x=0.8,y=5e5)
Example 2:
# smear model over full range and interpolate to needed data # this is the best way to smear a model for fitting, but is not possible for desmearing meas=js.dA('measureddata.dat') # load data # define profile resol2m = js.sas.prepareBeamProfile('SANS', detDist=2000,collDist=2000.,wavelength=0.4,wavespread=0.15) q=np.r_[0.01:5:0.01] # or np.r_[0:meas.X.min():0.01,meas.X,meas.X.max():meas.X.max()*2:0.1] # calc model temp=js.ff.ellipsoid(q,2,3) # smear it smearedmodel=js.sas.smear(temp,resol2m).interpolate(X=meas.X)
-
jscatter.smallanglescattering.
resFunctExplicit
(S, beamprofile, extrapolfunc=None)[source]¶ Resolution smearing of small angle scattering for SANS or SAXS according to explict given Gaussian width.
I(q0)= Integral{(R(q,q0)*S(q)}dq with Gaussian kernel R(q,q0). E.g. for merged dataFiles of KWS2@MLZ the explicit width is given in the 4th column.
Parameters: - S : array like
dataArray with X. and .Y theoretical Scattering function q in nm^-1 .Y can be arbitrary unit
- beamprofile : beamProfile ‘explicit’
Beam profile as prepared from prepareBeamProfile ‘explicit’
- extrapolfunc : list , default None
Type of extrapolation at low and high X edges for handling of the border as list for [low edge,high edge]. If a singe value is given this is used for both.
- Low edge
- float : Power law extrapolation of low X e.g. a=-4 for q**a for Porod scaling.
- None : A constant value as Y(X.min()).
- anything else : Low X data are log scaled, then X**2 extrapolated as Guinier like extrapolation.
- High edge
- float : Power law extrapolation of low X e.g. a=-4 for q**a for Porod scaling.
- None : A constant value as Y(X.max()).
Returns: - dataArray
columns [‘wavevector; smeared scattering; unsmeared scattering; half width smearing function’]
Notes
- HalfWidthSmearingFunction is the FWHM the Gaussian used for smearing including all effects.
-
jscatter.smallanglescattering.
smear
(data, beamProfile, **kwargs)[source]¶ Smearing data for line-collimated SAXS (Kratky camera) or as point collimation SANS/SAXS.
The full resolution for point collimation SAXS/SANS is described in resFunct.
Parameters: - data : dataArray
Data to be smeared.
- beamProfile : beamProfile or ‘trap’, ‘SANS’, ‘explicit’, dataArray
Beam profile as prepared from prepareBeamProfile or type as ‘trapezoidal’, ‘SANS’,’explicit’ or a measured beam profile as dataArray for line collimation. Measured Profile is treated by prepareBeamProfile.
- kwargs :
See prepareBeamProfile for kwargs.
Returns: - dataArray
Notes
- If data has attributes a, b, dIW, bxw, detDist these are used, if not given in function call.
- If wavelength is missing in data a default of 0.155418 nm for Xray K_{\alpha} line is assumed. For SANS 0.6 A.
- During smearing for Kratky camera an integration over the beam width and beam length are performed. In this integration q_{w,l}= ((q+q_{w})^2)+q_{l}^2)^{1/2} is used with q_{w} along the beam width and q_{l} along the beam length. in regions q_{w,l} > max(q_{data}) we estimate the measured scattering intensity by the mean of the last 10 points of the measured spectra to allow for a maximum in q range. The strictly valid q range can be estimated by calculating q_{x,y} < max(q) with 2 times the used beam width and beam length. As the smearing for larger q has no real effect the estimate might be still ok.
Examples
# use as # prepare measured line collimation beamprofile mbeam = js.sas.prepareBeamProfile(beam, bxw=0.01, dIW=1.) # prepare profile with trapezoidal shape (a,b can be fitted above) tbeam = js.sas.prepareBeamProfile('trapz', a=mbeam.a, b=mbeam.b, bxw=0.01, dIW=1) # prepare profile SANS (missing parameters get defaults, see resFunct) Sbeam = js.sas.prepareBeamProfile('SANS', detDist=2000,wavelength=0.4,wavespread=0.15) # prepare profile with explicit given Gaussian width in column 3 as e.g. KWS2@JCNS Gbeam = js.sas.prepareBeamProfile(measurement,explicit=3) # smear datasm= js.sas.smear(data,mbeam) datast= js.sas.smear(data,tbeam) datasS= js.sas.smear(data,Sbeam) datasG= js.sas.smear(data,Gbeam)
-
jscatter.smallanglescattering.
transmissionCorrection
(data, dark, emptybeam=None, edge=0.03, exposure=None)[source]¶ Subtract dark current, find primary beam, get transmission and normalize by transmission and exposure time. IN PLACE.
For measurements including the primary beam from a semitransparent beamstop as from SAXSpace. Transmission is the maximum of the primary beam peak after dark subtraction. Allows easier comparison of measurements with aged source (primary intensity change).
Parameters: - data : dataArray
A measurement from a SAXSpace instrument read by js.dA(‘filename.pdh’,lines2parameter=[2,3,4])
- dark : dataArray
Dark current measurement.
- emptybeam : dataArray,default=None
Empty beam measurement to subtract from measurement. dark will be subtracted of empty beam.
- edge : float, default 0.03
Wavevector value below beam stop edge. The primary beam is searched below this value.
- exposure : float, default None
Exposure time in unit ‘s’. If not given the xml description at the end of the file is examined.
Returns: - None
Notes
Files from SAXSpace instrument (.pdh) can be read by
js.dA('filename.pdh',lines2parameter=[2,3,4])
Correction
Brulet at al [1] describe the data correction for SANS, which is in principle also valid for SAXS, if incoherent contributions are neglected.
The difference is, that SAXS has typical transmission around ~0.3 for 1mm water sample in quartz cell due to absorption, while in SANS typical values are around ~0.9 for D2O. Larger volume fractions in the sample play a more important rule for SANS as hydrogenated ingredients reduce the transmission significantly, while in SAXS still the water and the cell (quartz) dominate.
One finds for a sample inside of a container with thicknesses (z) for sample, buffer (solvent), empty cell and empty beam measurement (omitting the overall q dependence):
I_s = \frac{1}{z_S}\big((\frac{I_S-I_{dark}}{T_S}-I_{b}T_S\big) - \big(\frac{I_{EC}-I_{dark}}{T_{EC}}-I_{b}T_{EC})\big) - \frac{1}{z_B}\big((\frac{I_B-I_{dark}}{T_B}-I_{b}T_B\big) - \big(\frac{I_{EC}-I_{dark}}{T_{EC}}-I_{b}T_{EC})\big)
- where
- I_s is the interesting species
- I_S is the sample of species in solvent (buffer)
- I_B is the pure solvent (describing a constant background)
- I_{dark} is the dark current measurement
- I_b is the empty beam measurement
- I_{EC} is the empty cell measurement
- z_x corresponding sample thickness
- T_x corresponding transmission
The recurring pattern \big((\frac{I-I_{dark}}{T}-I_{b}T\big) shows that the the beam tail (border of primary beam not absorbed by the beam stop) is attenuated by the corresponding sample.
For equal sample thickness z the empty beam is included in subtraction of I_B:
I_s = \frac{1}{z} \big((\frac{I_S-I_{dark}}{T_S}-I_{b}T_S) - (\frac{I_B-I_{dark}}{T_B}-I_{b}T_B)\big)
The simple case
If the transmissions are nearly equal as for e.g. protein samples with low concentration (T_S \approx T_B) we only need to subtract the transmission and dark current corrected buffer measurement from the sample.
I_s = \frac{1}{z} \big((\frac{I_S-I_{dark}}{T_S}) - (\frac{I_B-I_{dark}}{T_B}\big)
Higher accuracy for large volume fractions
For larger volume fractions \Phi the transmission might be different and we have to take into account that only 1-\Phi of solvent contributes to I_S. We may incorporate this in the sense of an optical density changing the effective thickness \frac{1}{z_B}\rightarrow\frac{1-\Phi}{z_B} resulting in different thicknesses z_S \neq z_B
Transmission
The transmission is measured as the ratio T=\frac{I(q=0)_{sample}}{I(q=0)_{emptybeam}} with I(q=0) as the primary beam intensity.
If the primary beam tail is neglected in the above equation I(q=0)_{emptybeam} only gives a common scaling factor and can be omitted if arbitrary units are used. Alternatively one can scale to the EC transmission with T_{EC}=1 For absolute calibration the same needs to be done. One finds T_{sample in cell}=T_{empty cell}T_{sample without cell}.
References
[1] (1, 2) Improvement of data treatment in small-angle neutron scattering Brûlet et al Journal of Applied Crystallography 40, 165-177 (2007)
-
jscatter.smallanglescattering.
waterXrayScattering
(composition='h2o1', T=293, units='mol')[source]¶ Absolute scattering of water with components (salt, buffer) at Q=0 as reference for X-ray.
According to [2] a buffer of water with components might be used. Ions need to be given separatly as [‘55.51h2o1’,‘0.15Na’,‘0.15Cl’] for 0.15 M NaCl solution. It is accounted for the temperature dependence of water density and compressibility.
Parameters: - composition : string
Buffer composition as in scatteringLengthDensityCalc give dissociated ions separatly as [‘1Na’,‘1Cl’] with concentration in mol prepended the additional scattering as ionic liquid of the ions in water is taken into account see [2] mass in g; 1000g water are 55.508 mol
- T : float
temperature in °K
- units : ‘mol’
anything except ‘mol’ prepended unit is mass fraction ‘mol’ prepended units is mol and mass fraction is calculated as mass=[mol] mass_{molecule} e.g. 1l Water with 123 mmol NaCl [‘55.508H2O1’,‘0.123Na1Cl1’]
Returns: - float absolute scattering length in Units 1/cm
Notes
I(0)=(\sigma_{water}^2f_e^2 n_{ew}^2 k_B T \chi + \sum_{ci} n_i N_A 1000 n_{ei}^2 f_e^2 )/100
with
\sigma_{water} water density \chi compressibility n_{ew} number of electrons per water molecule f_e cross section of electron in nm k_B Boltzmann constant n_i concentration component i n_{ei} number of electrons per molecule component i in Mol \sum_{ci} is done for all ions separately if givenReferences
[1] SAXS experiments on absolute scale with Kratky systems using water as a secondary standard Doris Orthaber et al. J. Appl. Cryst. (2000). 33, 218±225 [2] (1, 2, 3) A high sensitivity pinhole camera for soft condensed matter T. Zemb, O. Tache, F. Né, and O. Spalla, J. Appl. Crystallogr. 36, 800 (2003).
Read 2D image files (TIFF) from SAXS cameras and extract the corresponding data.
The sasImage is a 2D array that allows direct subtraction and multiplication (e.g. transmission) respecting given masks in operations. E.g.
sample=js.sas.sasImage('sample.tiff')
solvent=js.sas.sasImage('solvent.tiff')
corrected = sample/sampletransmission - solvent/solventtransmission
Calibration of detector distance, radial average, size reduction and more. .pickBeamcenter allows sensitive detection of the beamcenter.
An example is shown in sasImage
.
-
jscatter.sasimagelib.
AgBepeaks
= [1.0753, 2.1521, 3.2286, 4.3049, 5.3813, 6.4576, 7.5339, 8.6102, 9.6865, 10.7628]¶ AgBe peak positions
-
class
jscatter.sasimagelib.
SubArray
[source]¶ Bases:
numpy.ndarray
-
array
¶
-
attr
¶ Show specific attribute names as sorted list of attribute names.
-
setattr
(objekt, prepend='', keyadd='_')[source]¶ Set (copy) attributes from objekt.
Parameters: - objekt : objekt with attr or dictionary
Can be a dictionary of names:value pairs like {‘name’:[1,2,3,7,9]} If object has property attr the returned attribut names are copied.
- prepend : string, default ‘’
Prepend this string to all attribute names.
- keyadd : char, default=’_’
If reserved attributes (T, mean, ..) are found the name is ‘T’+keyadd
-
-
jscatter.sasimagelib.
createImageDescriptions
(images)[source]¶ Create text file with image descriptions as list of content.
Parameters: - images : list of sasImages or glob pattern
List of images
-
jscatter.sasimagelib.
createLogPNG
(filenames, center=None, size=None, colormap='jet', equalize=False, contrast=None)[source]¶ Create .png files from grayscale images with log scale conversion to values between [1,255].
This generates images viewable in simple image viewers as overview. The new files are stored in the same folder as the original files.
Parameters: - filenames : string
Filename with glob pattern as ‘file*.tif’
- center : [int,int]
Center of crop region.
- size : int
- Size of crop region.
- If center is given a box with 2*size around center is used.
- If center is None the border is cut by size.
- colormap : string, None
Colormap from matplotlib or None for grayscale. For standard colormap names look in mpl.showColors().
- equalize : bool
Equalize the images.
- contrast : None, float
Autocontrast for the image. The value (0.1=10%) determines how much percent are cut from the intensity histogram before linear spread of intensities.
-
jscatter.sasimagelib.
readImages
(filenames)[source]¶ Read a list of images returning sasImage`s.
Parameters: - filenames : string
Glob pattern to read
Returns: - list of sasImage`s
Notes
To get a list of image descriptions:
images=js.sas.readImages(path+'/latest*.tiff') [i.description for i in images]
-
class
jscatter.sasimagelib.
sasImage
[source]¶ Bases:
jscatter.sasimagelib.SubArray
,numpy.ma.core.MaskedArray
Creates sasImage as maskedArray from a detector image for evaluation.
- Reads a .tif file including the information in the EXIF tag.
- All methods of maskedArrays including masking of invalid areas work.
- Masked areas are automatically masked for all math operations.
- Arithmetic operations for sasImages work as for numpy arrays e.g. to subtract background image or multiplying with transmission. Use the numpy.ma methods.
- Coordinates are [height,width]
Parameters: - file : string
Filename to open.
- detector_distance : float, sasImage
Detector distance from calibration measurement or calibrated image. Overwrites value in the file EXIF tag.
- beamcenter : None, list 2xfloat, sasImage
Beamcenter is [height, width] position of primary beam. If sasImage is given the corresponding beamcenter is copied. Overwrites value given in the file EXIF tag.
- maskbelow : float, default =0
Mask values below this value.
Returns: - image : sasImage with attributes
- .beamcenter : beam center
- .iX : Height pixel positions
- .iY : Width pixel positions
- .filename
- .artist : Additional attributes from EXIF Tag Artist
- .imageDescription : Additional attributes from EXIF Tag ImageDescription
Notes
Unmasked data can be accessed as .data
The mask is .mask and initial set to all negative values.
Masking of a pixel is done as image[i,j]=np.ma.masked. Use mask methods as implemented.
Geometry mask methods can be used and additional masking methods from numpy masked Arrays.
import jscatter as js from numpy import ma cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal.mask = ma.nomask # reset mask cal[cal<0]= ma.mask # mask negative values cal[(cal>30) & (cal<100)] = ma.mask # mask region of values
TIFF tags with index above 700 are ignored.
Tested for reading tiff image files from Pilatus detectors as given from our metal jet SAXS machines Ganesha and Galaxi at JCNS, Jülich.
Additional SAXSpace TIFF files are supported which show frames per pixel on the Y axis. This allows to examine the time evolution of the measurement on these line collimation cameras (Kratky camera). Instead of the old PIL the newer fork Pillow is needed for the multi page TIFFs. Additional the pixel_size is set to 0.024 (µm) as for the JCNS CCD camera.
- Beamcenter & orientation:
The x,y orientation is not well defined and dependent on the implementation on the specific camera setup. The default used here corresponds to our in house Ganesha which needs to revert the EXIF beamcenter. We use the lower left image corner as zero with X as lower axis. Please check if your beamcenter corresponds to this. If not just change it.
References
[1] Everything SAXS: small-angle scattering pattern collection and correction Brian Richard Pauw J. Phys.: Condens. Matter 25, 383201 (2013) DOI https://doi.org/10.1088/0953-8984/25/38/383201 Examples
import jscatter as js # # Look at calibration measurement calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') # Check beamcenter # For correct beamcenter it should show straight lines (change beamcenter to see change) calibration.showPolar(beamcenter=[254,122],scaleR=3) # or use pickBeamcenter which seems to be more accurate calibration.pickBeamcenter() # Recalibrate with previous found beamcenter (calibration sets it already) calibration.recalibrateDetDistance(showfits=True) iqcal=calibration.radialAverage() # This might be used to calibrate detector distance for following measurements as # empty.setDetectorDistance(calibration) # empty = js.sas.sasImage(js.examples.datapath+'/emptycell.tiff') # Mask beamstop (not the same as calibration, unluckily) empty.mask4Polygon([185,92],[190,92],[233,0],[228,0]) empty.maskCircle(empty.beamcenter, 9) empty.show() buffer = js.sas.sasImage(js.examples.datapath+'/buffer.tiff') buffer.maskFromImage(empty) buffer.show() bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') bsa.maskFromImage(empty) bsa.show() # by default a log scaled image # # subtract buffer (transmission factor is just a guess here, sorry) new=bsa-buffer*0.2 new.show() # iqempty=empty.radialAverage() iqbuffer=buffer.radialAverage() iqbsa=bsa.radialAverage() # p=js.grace(1,1) p.plot(iqempty,le='empty cell') p.plot(iqbuffer,le='buffer') p.plot(iqbsa,le='bsa 11 mg/ml') p.title('raw data, no transmission correction') p.yaxis(min=1,max=1e3,scale='l',label='I(q) / a.u.') p.xaxis(scale='l',label='q / nm\S-1') p.legend()
-
array
¶ Strip of all attributes and return a simple array without mask.
-
asImage
(scale='log', colormap='jet', inverse=False, linthresh=1.0, linscale=1.0)[source]¶ Returns the sasImage as 8bit RGB image using PIL.
See PIL(Pillow) for more info about PIL images and image manipulation possibilities as e.g. in notes. Conversion to 8bit RGB looses floating point information but is for presenting and publication.
Parameters: - scale : ‘log’, ‘symlog’, default = ‘norm’
Scale for intensities.
- ‘norm’ Linear scale.
- ‘log’ Logarithmic scale
- ‘symlog’ Symmetrical logarithmic scale is logarithmic in both the positive and negative directions from the origin. Use linthresh, linscale to adjust.
- colormap : string, None
Colormap from matplotlib or None for grayscale. For standard colormap names look in js.mpl.showColors().
- inverse : bool
Inverse colormap.
- linthresh : float, default = 1
Only used for scale ‘sym’. The range within which the plot is linear (-linthresh to linthresh).
- linscale : float, default = 1
Only used for scale ‘sym’. Its value is the number of decades to use for each half of the linear range.
Returns: - PIL image
Notes
Pillow (fork of PIL) allows image manipulation. As a prerequisite of Jscatter it is installed on your system and can be imported as
import PIL
image=mysasimage.asImage() image.show() # show in system default viewer image.save('test.pdf', format=None, **params) # save the image in different formats image.save('test.jpg',subsampling=0, quality=100) # use these for best jpg quality image.save('test.png',transparency=(0,0,0)) # png image with black as transparent image.crop((10,10,200,200)) # cut border import PIL.ImageOps as PIO nimage=PIO.equalize(image, mask=None) # Equalize the image histogram. nimage=PIO.autocontrast(image, cutoff=0, ignore=None) # Automatic contrast nimage=PIO.expand(image, border=20, fill=(255,255,255)) # add border to image (here white)
Examples
import jscatter as js import PIL.ImageOps as PIO cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') # create image for later usage image=cal.asImage(colormap='inferno',scale='log',inverse=1) # create image and just show it cal.asImage(colormap='inferno',scale='log').show() # expand image and show it or save it PIO.expand(image, border=20, fill=(255,255,255)).show() PIO.expand(image, border=20, fill=(255,255,255)).save('myimageas.pdf')
-
asdataArray
(masked=0)[source]¶ Return representation of sasImage as dataArray representing wavevectors (qx,qy) against intensity.
Parameters: - masked : float, None, string, default=0
- How to deal with masked values.
- float : Set masked pixels to this value
- None : Remove from dataArray.
- To recover the image the masked pixels need to be interpolated on a regular grid.
- ‘linear’, ‘cubic’, ‘nearest’ : interpolate masked points by scipy.interpolate.griddata
- using specified order of interpolation on 2D image.
- ‘radial’ Use the radial averaged data to interpolate.
Returns: - dataArray with [qx,qy,I(qx,qy) ]
- .qx, .qz : original qx values to recover the image
Examples
This demo will show the interpolation in the masked regions of an artificial intensity distribution.
# manipulate data (not the mask) calibration.data[:150,30:60]=100 calibration.data[:150,60:90]=300 calibration.data[:150,90:]=600 # mask a circle calibration.maskCircle([100,100], 60) cal=calibration.asdataArray('linear') cal.Y[cal.Y<=0.1]=1.1 js.mpl.surface(cal.X, cal.Z, cal.Y) cal2=calibration.asdataArray(None) # this is reduced in size due to the mask
Use radial averaged data to interpolate
import jscatter as js import numpy as np bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') bsa.maskCircle(bsa.beamcenter,20) bsar=bsa.asdataArray('radial') js.mpl.surface(bsar.X, bsar.Z, bsar.Y)
-
findCenterOfIntensity
(beamcenter=None, size=100)[source]¶ Find beam center as center of intensity around beamcenter.
Only values above the mean value are used to calc center of intensity. Use an image with a clear symmetric and strong scattering sample as AgBe. Use .showPolar([600,699],scaleR=5) to see if peak is symmetric.
Parameters: - beamcenter : list 2x int
First estimate of beamcenter as [height, width] position. If not given preliminary beamcenter is estimated as center of intensity of full image.
- size : int
Defines size of rectangular region of interest (ROI) around the beamcenter to look at.
Returns: - Adds (replaces) .beamcenter as attribute.
Notes
If ROI is to large the result may be biased due to asymmetry of the intensity distribution inside of ROI.
-
gaussianFilter
(sigma=2)[source]¶ Gaussian filter in place.
Uses ndimage.filters.gaussian_filter with default parameters except sigma.
Parameters: - sigma : float
Gaussian kernel sigma.
-
getPixelQ
()[source]¶ Get scattering vector along pixel dimension around beamcenter.
Returns: - qx,qy with image x and y dimension
-
getfromcomment
(name, replace=None, newname=None)[source]¶ Extract name from .artist or .imageDescription with attribute name in front.
If multiple names start with parname first one is used. Used line is deleted from .artist or .imageDescription.
Parameters: - name : string
Name of the parameter in first place.
- replace : dict
Dictionary with pairs to replace in all lines.
- newname : string
New attribute name
-
iX
¶ X pixel coordinates
-
iY
¶ Y pixel coordinates
-
interpolateMaskedRadial
(radial=None)[source]¶ Interpolate masked values from radial averaged image or function.
This can be used to “extrapolate” over masked regions if e.g a background was measured at wrong distance.
Parameters: - radial : dataArray, function, default = None
- Determines how to determine masked values based on radial q values from beamcenter.
- Function accepting array to calculate masked data.
- dataArray for linear interpolating masked points.
- None uses the radialAverage image.
Returns: - sasImage including original parameters.
Examples
Use radial averaged data to interpolate
import jscatter as js import numpy as np calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal=calibration.interpolateMaskedRadial() # or # cal=calibration.interpolateMaskedRadial(calibration.radialAverage()) cal.show()
Generate image for different detector distance
cal.setDetectorDistance(0.3) # mask whole image cal.mask=np.ma.masked # recover image with radial average from original cal2=cal.interpolateMaskedRadial(calibration.radialAverage()) cal2.show()
-
lineAverage
(beamcenter=None, number=None, minmax='auto', show=False)[source]¶ Line average of image and conversion to wavevector q for line collimation cameras.
Remember to set .detector_distance to calibrated value.
Parameters: - beamcenter : float
Sets beam center in data and uses this. If not given the beam center is determined from semitransparent beam.
- number : int, default None
Number of intervals on new X scale. None means all pixels.
- minmax : [int,int], ‘auto’
Interval for determination of beamcenter.
- show : bool
Show the fit of the primary beam.
Returns: - dataArray
- .filename
- .detector_distance
- .description
- .beamcenter
Notes
- Detector distance in attributes is used.
- The primary beam is automatically detected.
- Correction for flat detector projected to Ewald sphere included.
-
mask4Polygon
(p1, p2, p3, p4, invert=False)[source]¶ Mask inside polygon of 4 points.
Points need to be given in right hand order.
Parameters: - p1,p2,p3,p4 : list of 2x float
Edge points.
- invert : bool
Invert region. Mask outside circle.
-
maskCircle
(center, radius, invert=False)[source]¶ Mask points inside circle.
Parameters: - center : list of 2x float
Center point.
- radius : float
Radius in pixel units
- invert : bool
Invert region. Mask outside circle.
-
maskFromImage
(image)[source]¶ Use/copy mask from image.
Parameters: - image : sasImage
sasImage to use mask for resetting mask. image needs to have same dimension.
-
maskRegion
(xmin, xmax, ymin, ymax)[source]¶ Mask rectangular region.
Parameters: - xmin,xmax,ymin,ymax : int
Corners of the region to mask
-
maskRegions
(regions)[source]¶ Mask several regions.
Parameters: - regions : list
List of regions as in maskRegion.
-
maskReset
()[source]¶ Reset the mask.
By default values smaller 0 are automatically masked again as is also default for reading
-
maskSectors
(angles, width, radialmax=None, beamcenter=None, invert=False)[source]¶ Mask sector around beamcenter.
Zero angle is
Parameters: - angles : list of float
Center angles of sectors in grad.
- width : float or list of float
Width of the sectors in grad. If single value all sectors are equal.
- radialmax : float
Maximum radius in pixels.
- beamcenter : 2x float
Center if different from stored beamcenter.
- invert : bool
Invert mask or not.
Examples
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal.maskSectors([0,90,180],20,radialmax=100,invert=True) cal.show()
-
maskTriangle
(p1, p2, p3, invert=False)[source]¶ Mask inside triangle.
Parameters: - p1,p2,p3 : list of 2x float
Edge points of triangle.
- invert : bool
Invert region. Mask outside circle.
-
maskbelowLine
(p1, p2)[source]¶ Mask points at one side of line.
The masked side is left looking from p1 to p2.
Parameters: - p1, p2 : list of 2x float
Points in pixel coordinates defining line.
-
pickBeamcenter
(levels=8, symmetry=6)[source]¶ Open image to pick the beamcenter from a calibration sample as AgBe.
Radial averaged sectors allow to find the optimal beamcenter with best overlap of peaks. Closing the image accepts the actual selected beamcenter.
Parameters: - levels : int
Number of levels in contour image.
- symmetry : int
Number of sectors around beamcenter for radial averages.
Returns: - After closing the selected beamcenter is saved in the sasImage.
Notes
A figure with the AgBe picture (right) and a radial average over sectors is shown (left, symmetry defines number of sectors) .
- A circle is shown around the mouse point after clicking. By default the radius corresponds to an AgBe reflex. The radius can be changed by mouse wheel or +/- keys.
- Radius of the circle (center of left plot data) can be increased/decreased by +/-.
- Width around radius can be increased/decrease by ctrl++/ctrl+-.
- The beamcenter can be changed by the mouse pointer (click).
- The beamcenter can be moved by arrow keys (+-1) or ctrl+arrow (+-0.1)
- In the sectors a radial average (after some smoothing) is calculated and shown in the left axes.
- The beamcenter is OK if the peaks show maximum overlap.
Examples
import jscatter as js # calibration = js.sas.sasImage(js.examples.datapath+’/calibration.tiff’) # use pickBeamcenter calibration.pickBeamcenter()
-
radialAverage
(beamcenter=None, number=300, kind='log')[source]¶ Radial average of image and conversion to wavevector q.
Remember to set .detector_distance to calibrated value.
Parameters: - beamcenter : list 2x float
Sets beam center or radial center in data and uses this. If not given the attribut beamcenter in the data is used.
- number : int, default 500
Number of intervals on new X scale.
- kind : ‘lin’, default ‘log’
Determines how points are distributed.
Returns: - dataArray
Notes
- Correction of pixel size for flat detector projected to Ewald sphere included.
- The value in a q binning is the average count rate c(q)=(\sum c_i)/N with counts in pixel i c_i and number of pixels N
- The error (standard deviation) is calculated in a q binning as e=(\sum c_i)^{1/2}/N
- The error is valid for single photon counting detectors showing Poisson statistics
- as the today typical Pilatus detectors from DECTRIS.
Examples
Mask and do radial average over sectors.
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') p=js.grace() calc=cal.copy() calc.maskSectors([0,180],20,radialmax=100,invert=True) calc.show() icalc=calc.radialAverage() p.plot(icalc,le='horizontal') calc=cal.copy() calc.maskSectors([90+0,90+180],20,radialmax=100,invert=True) calc.show() icalc=calc.radialAverage() p.plot(icalc,le='vertical') p.yaxis(scale='l') p.legend() p.title('The AgBe is not isotropically ordered')
-
recalibrateDetDistance
(beamcenter=None, number=500, showfits=False)[source]¶ Recalibration of detectorDistance by AgBe reference for point collimation.
Use only for AgBe reference measurements to determine the correction factor. For non AgBe measurements set during reading or .detector_distance to the new value. May not work if the detector distance is totally wrong.
Parameters: - beamcenter : list 2x float
Sets beam center or radial center in data and uses this. If not given the attribut beamcenter in the data is used.
- number : int, default 1000
number of intervals on new X scale.
- showfits : bool
Show the AgBe peak fits.
Notes
- .distanceCorrection will contain factor for correction. Repeating this results in a .distanceCorrection close to 1.
-
reduceSize
(pixelsize=2, center=None, border=None)[source]¶ Reduce size of image using uniform average in box.
XResolution,YResolution,beamcenter are scaled correspondingly.
Parameters: - pixelsize : int
Pixel size of the box to average. Also factor for reduction.
- center : [int,int]
Center of crop region.
- border : int
- Size of crop region.
- If center is given a box with 2*size around center is used.
- If center is None the border is cut by size.
Returns: - sasImage
-
saveAsTIF
(filename, fill=None, **params)[source]¶ Save the sasImage as float32 tif without loosing information.
Conversion from float64 to float32 is necessary. To save colored images use asImage.save() (see
asImage()
)Parameters: - filename : string
Filename to save to.
- fill : float, ‘min’ default None
Fill value for masked values. By default this is -1.
- ‘min’ uses the minimal value of the respective data type
- np.iinfo(np.int32).min = -2147483648 for int32
- np.finfo(np.float32).min = -3.4028235e+38 for float32
- params : kwargs
Additional kwargs for PIL.Image.save if needed.
Examples
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal2=cal/2. cal2.saveAsTIF('mycal',fill=-100) mycal=js.sas.sasImage('mycal.tif',maskbelow=-200) mycal.show()
-
setBeamcenter
(beamcenter)[source]¶ Set beamcenter .
Parameters: - beamcenter : float, sasImage
New value for beamcenter as [height, width]. If sasImage the beamcenter is copied.
-
setDetectorDistance
(detector_distance, offset=0)[source]¶ Set detector distance from calibration .
Parameters: - detector_distance : float, sasImage
New value for detector distance. If sasImage the detector_distance is copied.
- offset : float
Offset for sample compared to calibration sample.
-
show
(**kwargs)[source]¶ Show sasImage as matplotlib figure.
Parameters: - scale : ‘log’, ‘symlog’, default = ‘norm’
Scale for intensities.
- ‘norm’ Linear scale.
- ‘log’ Logarithmic scale
- ‘symlog’ Symmetrical logarithmic scale is logarithmic in both the positive and negative directions from the origin. This works also for only positive data. Use linthresh, linscale to adjust.
- levels : int, None
Number of contour levels.
- axis : ‘pixel’, None
Force pixel as coordinates, otherwise wavevectors if possible.
- colorMap : string
Get a colormap instance from name. Standard mpl colormap name (see showColors).
- linthresh : float, default = 1
Only used for scale ‘symlog’. The range within which the plot is linear (-linthresh to linthresh).
- linscale : float, default = 1
Only used for scale ‘symlog’. Its value is the number of decades to use for each half of the linear range. E.g. 10 uses 1 decade.
- lineMap : string
Label color Colormap name as in colorMap, otherwise as cs in in Axes.clabel * if None, the color of each label matches the color of the corresponding contour * if one string color, e.g., colors = ‘r’ or colors = ‘red’, all labels will be plotted in this color * if a tuple of matplotlib color args (string, float, rgb, etc),
different labels will be plotted in different colors in the order specified
- fontsize : int, default 10
Size of line labels in pixel
- title : None, string
Title of the plot. May be set by fig.axes[0].set_title(‘title’)
- axis : None, ‘pixel’
If coordinates should be forced to pixel.
- invert_yaxis,invert_xaxis : bool
Invert corresponding axis.
- block : bool
Open in blocking or non-blocking mode
- origin : ‘lower’,’upper’
Origin of the plot. See matplotlib imshow.
Returns: - image handle
Examples
Use radial averaged data to interpolate
import jscatter as js import numpy as np calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') calibration.show(colorMap='ocean')
Use
scale='symlog'
for mixed lin=log scaling to pronounce low scattering.import jscatter as js import numpy as np # sets negative values to zero bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') fig=js.mpl.contourImage(bsa,scale='sym',linthresh=30, linscale=10)
-
showPolar
(beamcenter=None, scaleR=1, offset=0)[source]¶ Show image transformed to polar coordinates around beamcenter.
Azimuth corresponds: center line to beamcenter upwards, upper quarter beamcenter to right upper/lower edge = beamcenter downwards, lower quarter beamcenter to left
Parameters: - beamcenter : [int,int]
Beamcenter
- scaleR : float
Scaling factor for radial component to zoom the beamcenter.
- offset : float
Offset to remove beamcenter from polar image.
Returns: - Handle to figure
-
jscatter.sasimagelib.
shortprint
(values, threshold=6, edgeitems=2)[source]¶ Creates a short handy representation string for array values.
Parameters: - values : object
Values to print.
- threshold: int default 6
Number of elements to switch to reduced form.
- edgeitems : int default 2
Items at the edge.
-
jscatter.sasimagelib.
subarray
¶ alias of
jscatter.sasimagelib.SubArray
-
class
jscatter.sasimagelib.
sasImage
[source] Creates sasImage as maskedArray from a detector image for evaluation.
- Reads a .tif file including the information in the EXIF tag.
- All methods of maskedArrays including masking of invalid areas work.
- Masked areas are automatically masked for all math operations.
- Arithmetic operations for sasImages work as for numpy arrays e.g. to subtract background image or multiplying with transmission. Use the numpy.ma methods.
- Coordinates are [height,width]
Parameters: - file : string
Filename to open.
- detector_distance : float, sasImage
Detector distance from calibration measurement or calibrated image. Overwrites value in the file EXIF tag.
- beamcenter : None, list 2xfloat, sasImage
Beamcenter is [height, width] position of primary beam. If sasImage is given the corresponding beamcenter is copied. Overwrites value given in the file EXIF tag.
- maskbelow : float, default =0
Mask values below this value.
Returns: - image : sasImage with attributes
- .beamcenter : beam center
- .iX : Height pixel positions
- .iY : Width pixel positions
- .filename
- .artist : Additional attributes from EXIF Tag Artist
- .imageDescription : Additional attributes from EXIF Tag ImageDescription
Notes
Unmasked data can be accessed as .data
The mask is .mask and initial set to all negative values.
Masking of a pixel is done as image[i,j]=np.ma.masked. Use mask methods as implemented.
Geometry mask methods can be used and additional masking methods from numpy masked Arrays.
import jscatter as js from numpy import ma cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal.mask = ma.nomask # reset mask cal[cal<0]= ma.mask # mask negative values cal[(cal>30) & (cal<100)] = ma.mask # mask region of values
TIFF tags with index above 700 are ignored.
Tested for reading tiff image files from Pilatus detectors as given from our metal jet SAXS machines Ganesha and Galaxi at JCNS, Jülich.
Additional SAXSpace TIFF files are supported which show frames per pixel on the Y axis. This allows to examine the time evolution of the measurement on these line collimation cameras (Kratky camera). Instead of the old PIL the newer fork Pillow is needed for the multi page TIFFs. Additional the pixel_size is set to 0.024 (µm) as for the JCNS CCD camera.
- Beamcenter & orientation:
The x,y orientation is not well defined and dependent on the implementation on the specific camera setup. The default used here corresponds to our in house Ganesha which needs to revert the EXIF beamcenter. We use the lower left image corner as zero with X as lower axis. Please check if your beamcenter corresponds to this. If not just change it.
References
[1] Everything SAXS: small-angle scattering pattern collection and correction Brian Richard Pauw J. Phys.: Condens. Matter 25, 383201 (2013) DOI https://doi.org/10.1088/0953-8984/25/38/383201 Examples
import jscatter as js # # Look at calibration measurement calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') # Check beamcenter # For correct beamcenter it should show straight lines (change beamcenter to see change) calibration.showPolar(beamcenter=[254,122],scaleR=3) # or use pickBeamcenter which seems to be more accurate calibration.pickBeamcenter() # Recalibrate with previous found beamcenter (calibration sets it already) calibration.recalibrateDetDistance(showfits=True) iqcal=calibration.radialAverage() # This might be used to calibrate detector distance for following measurements as # empty.setDetectorDistance(calibration) # empty = js.sas.sasImage(js.examples.datapath+'/emptycell.tiff') # Mask beamstop (not the same as calibration, unluckily) empty.mask4Polygon([185,92],[190,92],[233,0],[228,0]) empty.maskCircle(empty.beamcenter, 9) empty.show() buffer = js.sas.sasImage(js.examples.datapath+'/buffer.tiff') buffer.maskFromImage(empty) buffer.show() bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') bsa.maskFromImage(empty) bsa.show() # by default a log scaled image # # subtract buffer (transmission factor is just a guess here, sorry) new=bsa-buffer*0.2 new.show() # iqempty=empty.radialAverage() iqbuffer=buffer.radialAverage() iqbsa=bsa.radialAverage() # p=js.grace(1,1) p.plot(iqempty,le='empty cell') p.plot(iqbuffer,le='buffer') p.plot(iqbsa,le='bsa 11 mg/ml') p.title('raw data, no transmission correction') p.yaxis(min=1,max=1e3,scale='l',label='I(q) / a.u.') p.xaxis(scale='l',label='q / nm\S-1') p.legend()
-
array
Strip of all attributes and return a simple array without mask.
-
asImage
(scale='log', colormap='jet', inverse=False, linthresh=1.0, linscale=1.0)[source] Returns the sasImage as 8bit RGB image using PIL.
See PIL(Pillow) for more info about PIL images and image manipulation possibilities as e.g. in notes. Conversion to 8bit RGB looses floating point information but is for presenting and publication.
Parameters: - scale : ‘log’, ‘symlog’, default = ‘norm’
Scale for intensities.
- ‘norm’ Linear scale.
- ‘log’ Logarithmic scale
- ‘symlog’ Symmetrical logarithmic scale is logarithmic in both the positive and negative directions from the origin. Use linthresh, linscale to adjust.
- colormap : string, None
Colormap from matplotlib or None for grayscale. For standard colormap names look in js.mpl.showColors().
- inverse : bool
Inverse colormap.
- linthresh : float, default = 1
Only used for scale ‘sym’. The range within which the plot is linear (-linthresh to linthresh).
- linscale : float, default = 1
Only used for scale ‘sym’. Its value is the number of decades to use for each half of the linear range.
Returns: - PIL image
Notes
Pillow (fork of PIL) allows image manipulation. As a prerequisite of Jscatter it is installed on your system and can be imported as
import PIL
image=mysasimage.asImage() image.show() # show in system default viewer image.save('test.pdf', format=None, **params) # save the image in different formats image.save('test.jpg',subsampling=0, quality=100) # use these for best jpg quality image.save('test.png',transparency=(0,0,0)) # png image with black as transparent image.crop((10,10,200,200)) # cut border import PIL.ImageOps as PIO nimage=PIO.equalize(image, mask=None) # Equalize the image histogram. nimage=PIO.autocontrast(image, cutoff=0, ignore=None) # Automatic contrast nimage=PIO.expand(image, border=20, fill=(255,255,255)) # add border to image (here white)
Examples
import jscatter as js import PIL.ImageOps as PIO cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') # create image for later usage image=cal.asImage(colormap='inferno',scale='log',inverse=1) # create image and just show it cal.asImage(colormap='inferno',scale='log').show() # expand image and show it or save it PIO.expand(image, border=20, fill=(255,255,255)).show() PIO.expand(image, border=20, fill=(255,255,255)).save('myimageas.pdf')
-
asdataArray
(masked=0)[source] Return representation of sasImage as dataArray representing wavevectors (qx,qy) against intensity.
Parameters: - masked : float, None, string, default=0
- How to deal with masked values.
- float : Set masked pixels to this value
- None : Remove from dataArray.
- To recover the image the masked pixels need to be interpolated on a regular grid.
- ‘linear’, ‘cubic’, ‘nearest’ : interpolate masked points by scipy.interpolate.griddata
- using specified order of interpolation on 2D image.
- ‘radial’ Use the radial averaged data to interpolate.
Returns: - dataArray with [qx,qy,I(qx,qy) ]
- .qx, .qz : original qx values to recover the image
Examples
This demo will show the interpolation in the masked regions of an artificial intensity distribution.
# manipulate data (not the mask) calibration.data[:150,30:60]=100 calibration.data[:150,60:90]=300 calibration.data[:150,90:]=600 # mask a circle calibration.maskCircle([100,100], 60) cal=calibration.asdataArray('linear') cal.Y[cal.Y<=0.1]=1.1 js.mpl.surface(cal.X, cal.Z, cal.Y) cal2=calibration.asdataArray(None) # this is reduced in size due to the mask
Use radial averaged data to interpolate
import jscatter as js import numpy as np bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') bsa.maskCircle(bsa.beamcenter,20) bsar=bsa.asdataArray('radial') js.mpl.surface(bsar.X, bsar.Z, bsar.Y)
-
findCenterOfIntensity
(beamcenter=None, size=100)[source] Find beam center as center of intensity around beamcenter.
Only values above the mean value are used to calc center of intensity. Use an image with a clear symmetric and strong scattering sample as AgBe. Use .showPolar([600,699],scaleR=5) to see if peak is symmetric.
Parameters: - beamcenter : list 2x int
First estimate of beamcenter as [height, width] position. If not given preliminary beamcenter is estimated as center of intensity of full image.
- size : int
Defines size of rectangular region of interest (ROI) around the beamcenter to look at.
Returns: - Adds (replaces) .beamcenter as attribute.
Notes
If ROI is to large the result may be biased due to asymmetry of the intensity distribution inside of ROI.
-
gaussianFilter
(sigma=2)[source] Gaussian filter in place.
Uses ndimage.filters.gaussian_filter with default parameters except sigma.
Parameters: - sigma : float
Gaussian kernel sigma.
-
getPixelQ
()[source] Get scattering vector along pixel dimension around beamcenter.
Returns: - qx,qy with image x and y dimension
-
getfromcomment
(name, replace=None, newname=None)[source] Extract name from .artist or .imageDescription with attribute name in front.
If multiple names start with parname first one is used. Used line is deleted from .artist or .imageDescription.
Parameters: - name : string
Name of the parameter in first place.
- replace : dict
Dictionary with pairs to replace in all lines.
- newname : string
New attribute name
-
iX
X pixel coordinates
-
iY
Y pixel coordinates
-
interpolateMaskedRadial
(radial=None)[source] Interpolate masked values from radial averaged image or function.
This can be used to “extrapolate” over masked regions if e.g a background was measured at wrong distance.
Parameters: - radial : dataArray, function, default = None
- Determines how to determine masked values based on radial q values from beamcenter.
- Function accepting array to calculate masked data.
- dataArray for linear interpolating masked points.
- None uses the radialAverage image.
Returns: - sasImage including original parameters.
Examples
Use radial averaged data to interpolate
import jscatter as js import numpy as np calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal=calibration.interpolateMaskedRadial() # or # cal=calibration.interpolateMaskedRadial(calibration.radialAverage()) cal.show()
Generate image for different detector distance
cal.setDetectorDistance(0.3) # mask whole image cal.mask=np.ma.masked # recover image with radial average from original cal2=cal.interpolateMaskedRadial(calibration.radialAverage()) cal2.show()
-
lineAverage
(beamcenter=None, number=None, minmax='auto', show=False)[source] Line average of image and conversion to wavevector q for line collimation cameras.
Remember to set .detector_distance to calibrated value.
Parameters: - beamcenter : float
Sets beam center in data and uses this. If not given the beam center is determined from semitransparent beam.
- number : int, default None
Number of intervals on new X scale. None means all pixels.
- minmax : [int,int], ‘auto’
Interval for determination of beamcenter.
- show : bool
Show the fit of the primary beam.
Returns: - dataArray
- .filename
- .detector_distance
- .description
- .beamcenter
Notes
- Detector distance in attributes is used.
- The primary beam is automatically detected.
- Correction for flat detector projected to Ewald sphere included.
-
mask4Polygon
(p1, p2, p3, p4, invert=False)[source] Mask inside polygon of 4 points.
Points need to be given in right hand order.
Parameters: - p1,p2,p3,p4 : list of 2x float
Edge points.
- invert : bool
Invert region. Mask outside circle.
-
maskCircle
(center, radius, invert=False)[source] Mask points inside circle.
Parameters: - center : list of 2x float
Center point.
- radius : float
Radius in pixel units
- invert : bool
Invert region. Mask outside circle.
-
maskFromImage
(image)[source] Use/copy mask from image.
Parameters: - image : sasImage
sasImage to use mask for resetting mask. image needs to have same dimension.
-
maskRegion
(xmin, xmax, ymin, ymax)[source] Mask rectangular region.
Parameters: - xmin,xmax,ymin,ymax : int
Corners of the region to mask
-
maskRegions
(regions)[source] Mask several regions.
Parameters: - regions : list
List of regions as in maskRegion.
-
maskReset
()[source] Reset the mask.
By default values smaller 0 are automatically masked again as is also default for reading
-
maskSectors
(angles, width, radialmax=None, beamcenter=None, invert=False)[source] Mask sector around beamcenter.
Zero angle is
Parameters: - angles : list of float
Center angles of sectors in grad.
- width : float or list of float
Width of the sectors in grad. If single value all sectors are equal.
- radialmax : float
Maximum radius in pixels.
- beamcenter : 2x float
Center if different from stored beamcenter.
- invert : bool
Invert mask or not.
Examples
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal.maskSectors([0,90,180],20,radialmax=100,invert=True) cal.show()
-
maskTriangle
(p1, p2, p3, invert=False)[source] Mask inside triangle.
Parameters: - p1,p2,p3 : list of 2x float
Edge points of triangle.
- invert : bool
Invert region. Mask outside circle.
-
maskbelowLine
(p1, p2)[source] Mask points at one side of line.
The masked side is left looking from p1 to p2.
Parameters: - p1, p2 : list of 2x float
Points in pixel coordinates defining line.
-
pickBeamcenter
(levels=8, symmetry=6)[source] Open image to pick the beamcenter from a calibration sample as AgBe.
Radial averaged sectors allow to find the optimal beamcenter with best overlap of peaks. Closing the image accepts the actual selected beamcenter.
Parameters: - levels : int
Number of levels in contour image.
- symmetry : int
Number of sectors around beamcenter for radial averages.
Returns: - After closing the selected beamcenter is saved in the sasImage.
Notes
A figure with the AgBe picture (right) and a radial average over sectors is shown (left, symmetry defines number of sectors) .
- A circle is shown around the mouse point after clicking. By default the radius corresponds to an AgBe reflex. The radius can be changed by mouse wheel or +/- keys.
- Radius of the circle (center of left plot data) can be increased/decreased by +/-.
- Width around radius can be increased/decrease by ctrl++/ctrl+-.
- The beamcenter can be changed by the mouse pointer (click).
- The beamcenter can be moved by arrow keys (+-1) or ctrl+arrow (+-0.1)
- In the sectors a radial average (after some smoothing) is calculated and shown in the left axes.
- The beamcenter is OK if the peaks show maximum overlap.
Examples
import jscatter as js # calibration = js.sas.sasImage(js.examples.datapath+’/calibration.tiff’) # use pickBeamcenter calibration.pickBeamcenter()
-
radialAverage
(beamcenter=None, number=300, kind='log')[source] Radial average of image and conversion to wavevector q.
Remember to set .detector_distance to calibrated value.
Parameters: - beamcenter : list 2x float
Sets beam center or radial center in data and uses this. If not given the attribut beamcenter in the data is used.
- number : int, default 500
Number of intervals on new X scale.
- kind : ‘lin’, default ‘log’
Determines how points are distributed.
Returns: - dataArray
Notes
- Correction of pixel size for flat detector projected to Ewald sphere included.
- The value in a q binning is the average count rate c(q)=(\sum c_i)/N with counts in pixel i c_i and number of pixels N
- The error (standard deviation) is calculated in a q binning as e=(\sum c_i)^{1/2}/N
- The error is valid for single photon counting detectors showing Poisson statistics
- as the today typical Pilatus detectors from DECTRIS.
Examples
Mask and do radial average over sectors.
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') p=js.grace() calc=cal.copy() calc.maskSectors([0,180],20,radialmax=100,invert=True) calc.show() icalc=calc.radialAverage() p.plot(icalc,le='horizontal') calc=cal.copy() calc.maskSectors([90+0,90+180],20,radialmax=100,invert=True) calc.show() icalc=calc.radialAverage() p.plot(icalc,le='vertical') p.yaxis(scale='l') p.legend() p.title('The AgBe is not isotropically ordered')
-
recalibrateDetDistance
(beamcenter=None, number=500, showfits=False)[source] Recalibration of detectorDistance by AgBe reference for point collimation.
Use only for AgBe reference measurements to determine the correction factor. For non AgBe measurements set during reading or .detector_distance to the new value. May not work if the detector distance is totally wrong.
Parameters: - beamcenter : list 2x float
Sets beam center or radial center in data and uses this. If not given the attribut beamcenter in the data is used.
- number : int, default 1000
number of intervals on new X scale.
- showfits : bool
Show the AgBe peak fits.
Notes
- .distanceCorrection will contain factor for correction. Repeating this results in a .distanceCorrection close to 1.
-
reduceSize
(pixelsize=2, center=None, border=None)[source] Reduce size of image using uniform average in box.
XResolution,YResolution,beamcenter are scaled correspondingly.
Parameters: - pixelsize : int
Pixel size of the box to average. Also factor for reduction.
- center : [int,int]
Center of crop region.
- border : int
- Size of crop region.
- If center is given a box with 2*size around center is used.
- If center is None the border is cut by size.
Returns: - sasImage
-
saveAsTIF
(filename, fill=None, **params)[source] Save the sasImage as float32 tif without loosing information.
Conversion from float64 to float32 is necessary. To save colored images use asImage.save() (see
asImage()
)Parameters: - filename : string
Filename to save to.
- fill : float, ‘min’ default None
Fill value for masked values. By default this is -1.
- ‘min’ uses the minimal value of the respective data type
- np.iinfo(np.int32).min = -2147483648 for int32
- np.finfo(np.float32).min = -3.4028235e+38 for float32
- params : kwargs
Additional kwargs for PIL.Image.save if needed.
Examples
import jscatter as js cal = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') cal2=cal/2. cal2.saveAsTIF('mycal',fill=-100) mycal=js.sas.sasImage('mycal.tif',maskbelow=-200) mycal.show()
-
setBeamcenter
(beamcenter)[source] Set beamcenter .
Parameters: - beamcenter : float, sasImage
New value for beamcenter as [height, width]. If sasImage the beamcenter is copied.
-
setDetectorDistance
(detector_distance, offset=0)[source] Set detector distance from calibration .
Parameters: - detector_distance : float, sasImage
New value for detector distance. If sasImage the detector_distance is copied.
- offset : float
Offset for sample compared to calibration sample.
-
show
(**kwargs)[source] Show sasImage as matplotlib figure.
Parameters: - scale : ‘log’, ‘symlog’, default = ‘norm’
Scale for intensities.
- ‘norm’ Linear scale.
- ‘log’ Logarithmic scale
- ‘symlog’ Symmetrical logarithmic scale is logarithmic in both the positive and negative directions from the origin. This works also for only positive data. Use linthresh, linscale to adjust.
- levels : int, None
Number of contour levels.
- axis : ‘pixel’, None
Force pixel as coordinates, otherwise wavevectors if possible.
- colorMap : string
Get a colormap instance from name. Standard mpl colormap name (see showColors).
- linthresh : float, default = 1
Only used for scale ‘symlog’. The range within which the plot is linear (-linthresh to linthresh).
- linscale : float, default = 1
Only used for scale ‘symlog’. Its value is the number of decades to use for each half of the linear range. E.g. 10 uses 1 decade.
- lineMap : string
Label color Colormap name as in colorMap, otherwise as cs in in Axes.clabel * if None, the color of each label matches the color of the corresponding contour * if one string color, e.g., colors = ‘r’ or colors = ‘red’, all labels will be plotted in this color * if a tuple of matplotlib color args (string, float, rgb, etc),
different labels will be plotted in different colors in the order specified
- fontsize : int, default 10
Size of line labels in pixel
- title : None, string
Title of the plot. May be set by fig.axes[0].set_title(‘title’)
- axis : None, ‘pixel’
If coordinates should be forced to pixel.
- invert_yaxis,invert_xaxis : bool
Invert corresponding axis.
- block : bool
Open in blocking or non-blocking mode
- origin : ‘lower’,’upper’
Origin of the plot. See matplotlib imshow.
Returns: - image handle
Examples
Use radial averaged data to interpolate
import jscatter as js import numpy as np calibration = js.sas.sasImage(js.examples.datapath+'/calibration.tiff') calibration.show(colorMap='ocean')
Use
scale='symlog'
for mixed lin=log scaling to pronounce low scattering.import jscatter as js import numpy as np # sets negative values to zero bsa = js.sas.sasImage(js.examples.datapath+'/BSA11mg.tiff') fig=js.mpl.contourImage(bsa,scale='sym',linthresh=30, linscale=10)
-
showPolar
(beamcenter=None, scaleR=1, offset=0)[source] Show image transformed to polar coordinates around beamcenter.
Azimuth corresponds: center line to beamcenter upwards, upper quarter beamcenter to right upper/lower edge = beamcenter downwards, lower quarter beamcenter to left
Parameters: - beamcenter : [int,int]
Beamcenter
- scaleR : float
Scaling factor for radial component to zoom the beamcenter.
- offset : float
Offset to remove beamcenter from polar image.
Returns: - Handle to figure