Page History
| Composition Setup |
|---|
cloak.toggle.type = text
cloak.toggle.open = + show
cloak.toggle.close = - hide
|
A quick walk-through of the tools that exist for analysis of xtc files with python.
The main focus is on pyana, and the examples are from and for XPP primarily,
but may be useful examples to other experiments too.
| Toggle Cloak | ||
|---|---|---|
|
...
| Table of Contents |
|---|
...
|
...
|
The Basics
Python
http://docs.python.org/tutorial/
Pyana
Analysis Workbook. Python-based Analysis
...
| Panel | ||
|---|---|---|
| ||
|
...
Loops through the xtc datagrams and dumps info to screen. I recommend piping it to 'less'.
| Panel | ||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| ||||||||||||||
Try the same with different verbosity levels:
|
...
| Code Block | ||||
|---|---|---|---|---|
| ||||
Scanning....
Start parsing files:
['/reg/d/psdm/xpp/xppi0310/xtc/e81-r0098-s00-c00.xtc', '/reg/d/psdm/xpp/xppi0310/xtc/e81-r0098-s01-c00.xtc']
14826 datagrams read in 4.120000 s . . . . . . .
-------------------------------------------------------------
XtcScanner information:
- 61 calibration cycles.
- Events per calib cycle:
[240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240]
Information from 1 control channels found:
fs2:ramp_angsft_target
Information from 11 devices found
BldInfo:EBeam: EBeamBld_V1 (14641)
BldInfo:FEEGasDetEnergy: FEEGasDetEnergy (14563) Any (78)
BldInfo:NH2-SB1-IPM-01: SharedIpimb (14641)
BldInfo:PhaseCavity: PhaseCavity (14641)
DetInfo:EpicsArch-0|NoDevice-0: Epics_V1 (107580)
DetInfo:NoDetector-0|Evr-0: EvrConfig_V4 (62) EvrData_V3 (14640)
DetInfo:XppSb2Ipm-1|Ipimb-0: IpimbConfig_V1 (1) IpmFexConfig_V1 (1) IpimbData_V1 (14640) IpmFex_V1 (14640)
DetInfo:XppSb3Ipm-1|Ipimb-0: IpimbConfig_V1 (1) IpmFexConfig_V1 (1) IpimbData_V1 (14640) IpmFex_V1 (14640)
DetInfo:XppSb3Pim-1|Ipimb-0: IpimbConfig_V1 (1) IpmFexConfig_V1 (1) IpimbData_V1 (14640) IpmFex_V1 (14640)
DetInfo:XppSb4Pim-1|Ipimb-0: IpimbConfig_V1 (1) IpmFexConfig_V1 (1) IpimbData_V1 (14640) IpmFex_V1 (14640)
ProcInfo:: RunControlConfig_V1 (62)
XtcScanner is done!
-------------------------------------------------------------
|
xtcexplorer
XTC Explorer - Old - GUI interface that builds pyana modules for you.
...
| Panel | |||||||
|---|---|---|---|---|---|---|---|
| |||||||
(' Now you have a local version of the XtcExplorer package in your directory. That allows you to edit the Exercise for later: |
...
For some more useful analysis examples, in the following we'll stick to writing customized pyana modules and running pyana from the command line.
Extracting the data with pyana, some examples
Outline of a pyana module
Like the other frameworks, pyana is an executable that loops through the XTC file and calls all
requested user modules at certain transitions. All the analysts need to do is to fill in the
relevant functions in their user analysis module:
| Toggle Cloak | ||
|---|---|---|
|
...
| id | kappe |
|---|
...
But before getting to the pyana modules, I'll briefly touch on a few items general to python that may be useful: saving files, matplotlib for plotting, and IPython for interactive work.
NumPy, SciPy and MatPlotLib
These are packages that you may want to look into. Pretty much all our examples here are using them:
- http://numpy.scipy.org/ - Numerical package for python (arrays etc)
- http://www.scipy.org/ - Scientific tools
- http://matplotlib.sourceforge.net/ - plotting package
Other useful links:
- http://www.scipy.org/NumPy_for_Matlab_Users
- http://www.scipy.org/
- http://www.sagemath.org/
- http://code.google.com/p/spyderlib/
Saving data arrays
Here are a few examples of how you can save data arrays in python.
| Panel | ||
|---|---|---|
| ||
|
| Panel | ||||
|---|---|---|---|---|
| ||||
|
| Panel | ||||
|---|---|---|---|---|
|
For the following two examples, check out the latest version of the pyana_examples package:
...
addpkg pyana_examples HEAD
scons
Point detector delay scan
| Panel | |||||||
|---|---|---|---|---|---|---|---|
| Code Block | none | none | |||||
[pyana]
modules pyana_examples.xppt_delayscan
[pyana_examples.xppt_delayscan]
controlpv = fs2:ramp_angsft_target
ipimb_norm = XppSb3Ipm-1|Ipimb-0
ipimb_sig = XppSb3Pim-1|Ipimb-0
threshold = 0
outputfile = point_scan_delay.npy
| |||||||
|
| Toggle Cloak | ||
|---|---|---|
|
xppt_delayscan.py:...
| id | delayscan |
|---|
...
Image peak finding
Here are a collection of useful algorithms for image analysis: http://docs.scipy.org/doc/scipy/reference/ndimage.html
This particular example is done with a CSPad image, but only a single section is available. For more typical CSPad module, see next section.
...
def event( self, evt, env ) :
elements = evt.getCsPadQuads(self.source, env)
image = elements[0].data().reshape(185, 388)
...
For more examples, see How to access HDF5 data from Python and http://code.google.com/p/h5py/ |
Plotting with MatPlotLib
One of the most commonly used tools for plotting in python: matplotlib. Other alternatives exist too.
Matplotlib:
- The plotting can be done directly in the pyana module, but be aware that you need to disable plotting for the
module to run successfully in a batch job.Code Block import matplotlib.pyplot as plt plt.plot(array) plt.show()
- Or you can load arrays from a file and interactively plot them in iPython. The same ('recommended') syntax as above can be used, or if you use 'import *' you don't need to prepend the commands with the package name, which is handy when plotting interactively:
Code Block from matplotlib.pyplot import * ion() plot(array) draw()
Interactive analysis with IPython
The LCLS offline analysis group does have plans for a real interactive pyana, but currently this is not available.
2011-11-04 iPsana Interactive Analysis Framework.pdf
The version available in our offline release system is
IPython 0.9.1 – An enhanced Interactive Python.
so this is the one I've been using in these examples.
Not a whole lot more than a python shell.
However, the latest IPython has loads of new and interesting features...
| Panel | ||
|---|---|---|
| ||
This example reads in a file produced by the "point detector delay scan" example below.
Sometimes you need to issue the draw() command twice, for some reason. After drawing you can keep working on the arrays and plot more... |
Extracting the data with pyana, some examples
Outline of a pyana module
Like the other frameworks, pyana is an executable that loops through the XTC file and calls all
requested user modules at certain transitions. All the analysts need to do is to fill in the
relevant functions in their user analysis module:
| Toggle Cloak | ||
|---|---|---|
|
...
| Cloak | ||
|---|---|---|
| ||
| Code Block |
|---|
# useful imports
import numpy as np
import matplotlib.pyplot as plt
from pypdsdata.xtc import TypeId
class mymodule (object) :
"""Class whose instance will be used as a user analysis module. """
def __init__ ( self, |
...
# Region of Interest (RoI)
if self.roi is None:
self.roi = [ 0, image.shape[1], 0, image.shape[0] ] # [x1,x2,y1,y2]
print "ROI [x1, x2, y1, y2] = ", self.roi
...
roi_array = image[self.roi[2]:self.roi[3],self.roi[0]:self.roi[1]]
cms = scipy.ndimage.measurements.center_of_mass(roi_array)
print "Center-of-mass of the ROI: (x, y) = (%.2f, %.2f)" %(self.roi[0]+cms[1],self.roi[2]+cms[0])
...
|
...
|
...
|
...
|
...
source = "" |
...
|
...
|
...
|
...
|
...
threshold = "" ): """Class constructor. The parameters |
...
to |
...
the |
...
constructor |
...
are passed from pyana configuration file. |
...
|
...
If parameters do not have default values here then the must be defined in |
...
pyana.cfg. All parameters are passed as strings, convert to correct type before use. |
...
|
...
|
...
|
...
|
...
@param |
...
source |
...
|
...
name of device, format 'Det-ID|Dev-ID' @param threshold threshold value (remember to convert from string) """ |
...
self.source = source self.threshold = float(threshold) def beginjob( self, evt, |
...
env |
...
) : """This method is called |
...
once at the beginning of the job. It should do a one-time initialization possible extracting |
...
values |
...
from event data (which is a Configure object) or |
...
environment. @param evt event data object |
...
|
...
@param env environment object |
...
|
...
""" pass def beginrun( |
...
self, evt, env ) : """This optional method is called if present at the beginning of the new run. @param evt event data object @param env environment object |
...
CSPad images and tile arangements
CSPad data structure
CSPad data in xtc is a list of elements. In pyana get the list from the evt (event) object (notice the need for the env (environment) object too!):
...
elements = evt.getCsPadQuads(self.source, env)
elements here is a python list of ElementV1 or ElementV2 (or later versions) objects, each representing one quadrant. The list is not ordered, so to know which quadrant you have, you have to check with element.quad(). To store a local array of the whole CSPad detector, you can do the following.
...
""" pass def begincalibcycle( self, evt, env ) : """This optional method is called if present at the beginning of the new calibration cycle. @param evt event data object |
...
|
...
@param env |
...
environment object |
...
|
...
|
...
""" pass |
...
def event( self, evt, env ) : """This method is called for every L1Accept |
...
transition. @param evt event data |
...
object |
...
@param env environment object |
...
|
...
|
...
|
...
|
...
""" |
...
|
...
|
...
pass |
...
|
...
def |
...
endcalibcycle( self, env ) : """This optional method is called if present at the end of the |
...
|
...
calibration cycle. @param env |
...
|
...
|
...
environment object """ |
...
pass |
...
|
...
|
...
|
...
def |
...
endrun( |
...
self, |
...
env |
...
) |
...
: |
...
"""This optional method is called if present at the end of the run. @param env environment object |
...
""" pass |
...
def |
...
endjob( |
...
self, |
...
env |
...
) : """This method is called at the end of |
...
the |
...
job. |
...
It |
...
should do final cleanup, e.g. close all open files. |
...
|
...
@param env environment object """ |
...
|
...
|
...
pass |
...
| Cloak |
|---|
| Panel |
|---|
For the following two examples, check out the latest version of the |
:
|
...
|
...
(Note, if you don't already have a Kerberos ticket, you need to issue a ' |
Datatypes, and how to find data from your detector/device in the xtc file.
Pyana and psana has follows this naming scheme for labeling the datatypes from various devices. You can find the
names by investigating the xtc file with the above-mentioned tools (pyxtcreader, xtcscanner, xtcexplorer).
To see some examples of how to fetch the various data types in pyana (or psana), look at Devices and Datatypes.
Point detector delay scan
The python code for this pyana module resides in pyana_examples/src/xppt_delayscan.py. In this example, we do a point detector delay scan, where we get the time as scan points via a control PV, and where time rebinning based on phase cavity measurement is used to improve the time resolution. One IPIMB device (a.k.a. IPM3) is used for normalization (i0, I Zero) (parameter name ipimb_norm) and another IPIMB device (a.k.a. PIM3) channel 1 is used as the signal (parameter name ipimb_sig).
| Panel | ||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Open an editor and save the following in a file named pyana.cfg:
If you look at the code (pyana_examples/src/xppt_delayscan.py) you'll notice there are no detector names in there. The names of the detectors in the XTC file are passed as parameters from the configuration file above. The Run pyana (start with 200 events):
|
| Toggle Cloak | ||
|---|---|---|
|
...
| Cloak | ||
|---|---|---|
| ||
- Fetching the ControlPV information:
ControlPV is available from theenvobject, and since it only changes at the beginning
of each calibration cycle, thebegincalibcyclefunction is the appropriate place to get it:Code Block none none def begincalibcycle( self, evt, env ) :
The ControlConfig object may contain several pvControl and pvMonitor objects. In this case
there's only one, but make sure the name matches anyway:Code Block none none ctrl_config = env.getConfig(TypeId.Type.Id_ControlConfig
CSPad tile arrangement
To get a rough picture of the full detector, here's an example of how XtcExplorer/src/cspad.py does it:
- For each Quadrant (cspad_layout.txt):
Code Block none none def get_quad_image( self, data3d, qn) : """get_quad_image Get an image for this quad (qn) @param data3d 3d data array (row vs. col vs. section) @param qn quad number """ pairs = [] for i in range (8) : # 1) insert gap between asics in the 2x1 asics = np.hsplit( data3d[i], 2) gap = np.zeros( (185,3), dtype=data3d.dtype ) # # gap should be 3 pixels wide pair = np.hstack( (asics[0], gap, asics[1]) ) for ic in range # all sections are originally 185 (rows) x 388 (columns) (0, ctrl_config.npvControls() ): cpv = ctrl_config.pvControl(ic) # Re-orient each section in the quadif cpv.name()=="fs2:ramp_angsft_target": if i==0 or i==1 : # store the value in a class variable (visible in every class method) pair = pair[:,::-1].T # reverse columns, switch columns to rows. if i==4 or i==5 :self.current_pv_value = cpv.value() )
- Fetching the IPIMB and PhaseCavity information:
All the other information that we need, is available through theevtobject, and
eventmember function is the place to get it:Code Block none none def event( self, evt, env pair = pair[::-1,:].T # reverse rows, switch rows to columns ) :
Use "XppSb3Ipm-1|Ipimb-0" (a.k.a. IPM3) sum of all channels for normalization and filteringCode Block none none ipmN_raw = evt.get(TypeId.Type.Id_IpimbData, "XppSb3Ipm-1|Ipimb-0") ipmN_fex = pairsevt.append( pair get(TypeId.Type.Id_IpmFex, "XppSb3Ipm-1|Ipimb-0") ipmN_norm if= ipmN_fex.sum
Use "XppSb3Pim-1|Ipimb-0" (a.k.a. PIM3) channel 1 as signalCode Block none none self.small_angle_tilt : pairipmS_raw = scipyevt.ndimageget(TypeId.interpolation.rotate(pair,self.tilt_array[qn][i]) # make the array for this quadrantType.Id_IpimbData, "XppSb3Pim-1|Ipimb-0" ) quadrantipmS_fex = npevt.zeros( (850, 850), dtype=data3d.dtypeget(TypeId.Type.Id_IpmFex, "XppSb3Pim-1|Ipimb-0" ) #ipm_sig = ipmS_fex.channel[1]
Get the phase cavity:Code Block none none insert the 2x1 sections according to forpc sec in range (8):= evt.getPhaseCavity() nrows, ncolsphasecav1 = pairs[sec].shape pc.fFitTime1 phasecav2 = pc.fFitTime2 # colp,rowp are where the top-left cornercharge1 of a section should be placed= pc.fCharge1 rowpcharge2 = 850 - self.sec_offset[0] - (self.section_centers[0][qn][sec] + nrows/2)pc.fCharge2
Compute delay time and fill histogramsCode Block none none colpdelaytime = 850 - self.sec_offset[1] - (self.section_centers[1][qn][sec] + ncols/2)current_pv_value + phasecav1*1e3 # The "histograms" are nothing but python lists. Append to them, and turn them into quadrant[rowp:rowp+nrows, colp:colp+ncols] = pairs[sec][0:nrows,0:ncols] arrays at the end. return quadrant
Then combine all four quadrant images into the full detector image:
...
self.h_ipm_rsig.append( ipm_sig ) self.
...
h_ipm_nsig.append( ipm_sig/ipm_norm )
...
...
self.
...
h_delaytime.append( delaytime )
| Cloak |
|---|
Image peak finding
Here are a collection of useful algorithms for image analysis: http://docs.scipy.org/doc/scipy/reference/ndimage.html
The python code for this pyana module example resides in pyana_examples/src/xppt_image_analysis.py.
This particular example is done with a CSPad image, but only a single section is available. For more typical CSPad module, see next section.
| Panel | |||||||||||||||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Edit
Then run the xppt_image_analysis pyana module on xppi0112 run 55:
*Edit
|
Here are some code snippet highlights from the xppt_image_analysis.py module:
- For each event, fetch the CsPad information, and get the image array:
Code Block none none def event( self, evt, env ) : elements = evt.getCsPadQuads(self.source, env) image = elements[0].data().reshape(185, 388) - Select a region of interest. If none is given (optional module parameter), set RoI to be the whole image.
Code Block none none # Region of Interest (RoI) if self.roi is None: self.roi = [ 0, image.shape[1], 0, image.shape[0] ] # [x1,x2,y1,y2] print "ROI [x1, x2, y1, y2] = ", self.roi - Using only the RoI subset of the image, compute center-of-mass using one of the SciPy.ndimamge algorithms. Add to it the position of the RoI to get the CMS in global pixel coordinates:
Code Block none none roi_array = image[self.roi[2]:self.roi[3],self.roi[0]:self.roi[1]] cms = scipy.ndimage.measurements.center_of_mass(roi_array) print "Center-of-mass of the ROI: (x, y) = (%.2f, %.2f)" %(self.roi[0]+cms[1],self.roi[2]+cms[0]) - Here's an example how you can make an interactive plot and select the Region of Interest with the mouse. Here we plot the image in two axes (subpads on the canvas). The first will always show the full image. In the second axes, you can select a rectangular region in "Zoom" mode (click on the Toolbar's Zoom button). The selected region will be drawn on top of the full image to the left, while the right plot will zoom into the selected region:
Code Block none none fig = plt.figure(1,figsize=(16,5)) axes1 = fig.add_subplot(121) axes2 = fig.add_subplot(122) axim1 = axes1.imshow(image) axes1.set_title("Full image") axim2 = axes2.imshow(roi_array, extent=(self.roi[0],self.roi[1],self.roi[3],self.roi[2])) axes2.set_title("Region of Interest") # rectangular ROI selector rect = UpdatingRect([0, 0], 0, 0, facecolor='None', edgecolor='red', picker=10) rect.set_bounds(*axes2.viewLim.bounds) axes1.add_patch(rect) # Connect for changing the view limits axes2.callbacks.connect('xlim_changed', rect) axes2.callbacks.connect('ylim_changed', rect) - To compute the center-of-mass of the selected region, revert back to non-zoom mode (hit the 'zoom' button again) and click on the rectangle. The rectangle is connected to the 'onpick' function which updates
self.roiand computes the center-of-mass:Code Block none none def onpick(event): xrange = axes2.get_xbound() yrange = axes2.get_ybound() self.roi = [ xrange[0], xrange[1], yrange[0], yrange[1]] roi_array = image[self.roi[2]:self.roi[3],self.roi[0]:self.roi[1]] cms = scipy.ndimage.measurements.center_of_mass(roi_array) print "Center-of-mass of the ROI: (x, y) = (%.2f, %.2f)" % (self.roi[0]+cms[1],self.roi[2]+cms[0]) fig.canvas.mpl_connect('pick_event', onpick) plt.draw()
CSPad images and tile arangements
The python code for this pyana module example resides in XtcExplorer/src/pyana_image.py.
| Panel | |||||||
|---|---|---|---|---|---|---|---|
Try some plotting of CSPad data using xtcexplorer. Launch the explorer and load xpp48712 run 66 (a dark run):
|
CSPad data structure
CSPad data in xtc is a list of elements. In pyana get the list from the evt (event) object (notice the need for the env (environment) object too!):
| Code Block | ||||
|---|---|---|---|---|
| ||||
elements = evt.getCsPadQuads(self.source, env)
|
elements here is a python list of ElementV1 or ElementV2 (or later versions) objects, each representing one quadrant. The list is not ordered, so to know which quadrant you have, you have to check with element.quad(). To store a local array of the whole CSPad detector, you can do the following.
- In beginjob, find out from the configuration object what part of the CSPad was in use (sometimes sections are missing):
Code Block none none def beginjob ( self, evt, env ) : config = env.getConfig(xtc.TypeId.Type.Id_CspadConfig, self.source) if not config: print '*** cspad config object is missing ***' return quads = range(4) # memorize this list of sections for later self.sections = map(config.sections, quads) - In each event, get the current CSPad data:
Code Block none none def event(self, evt, env): elements = evt.getCsPadQuads(self.source,env) pixel_array = np.zeros((4,8,185,388), dtype="uint16") for element in elements: data = element.data() # the 3-dimensional data array (list of 2d sections) quad = element.quad() # current quadrant number (integer value) # if any sections are missing, insert zeros if len( data ) < 8 : zsec = np.zeros( (185,388), dtype=data.dtype) for i in range (8) : if i not in self.sections[quad] : data = np.insert( data, i, zsec, axis=0 ) pixel_array[quad] = data
What we have so far gives you a 4d numpy array of all pixels. And if you want to store it in e.g. a numpy array, you can reshape it down to 2 dimensions (this is the format of the official pedestal files made by the translator):
| Code Block | ||||
|---|---|---|---|---|
| ||||
pixels = pixel_array.reshape(1480,388)
np.save("pixel_pedestal_file.npy", pixels )
|
CSPad tile arrangement
To get a rough picture of the full detector, here's an example of how XtcExplorer/src/cspad.py does it:
- For each Quadrant (cspad_layout.txt):
Code Block none none def get_quad_image( self, data3d, qn) : """get_quad_image Get an image for this quad (qn) @param data3d 3d data array (row vs. col vs. section) @param qn quad number """ pairs = [] for i in range (8) : # 1) insert gap between asics in the 2x1 asics = np.hsplit( data3d[i], 2) gap = np.zeros( (185,3), dtype=data3d.dtype ) # # gap should be 3 pixels wide pair = np.hstack( (asics[0], gap, asics[1]) ) # all sections are originally 185 (rows) x 388 (columns) # Re-orient each section in the quad if i==0 or i==1 : pair = pair[:,::-1].T # reverse columns, switch columns to rows. if i==4 or i==5 : pair = pair[::-1,:].T # reverse rows, switch rows to columns pairs.append( pair ) if self.small_angle_tilt : pair = scipy.ndimage.interpolation.rotate(pair,self.tilt_array[qn][i]) # make the array for this quadrant quadrant = np.zeros( (850, 850), dtype=data3d.dtype ) # insert the 2x1 sections according to for sec in range (8): nrows, ncols = pairs[sec].shape # colp,rowp are where the top-left corner of a section should be placed rowp = 850 - self.sec_offset[0] - (self.section_centers[0][qn][sec] + nrows/2) colp = 850 - self.sec_offset[1] - (self.section_centers[1][qn][sec] + ncols/2) quadrant[rowp:rowp+nrows, colp:colp+ncols] = pairs[sec][0:nrows,0:ncols] return quadrant
Then combine all four quadrant images into the full detector image:
| Code Block | ||||
|---|---|---|---|---|
| ||||
self.image = np.zeros((2*850+100, 2*850+100 ), dtype="float64")
for quad in xrange (4):
quad_image = self.get_quad_image( self.pixels[quad], quad )
self.qimages[quad] = quad_image
if quad>0:
# reorient the quad_image as needed
quad_image = np.rot90( quad_image, 4-quad)
qoff_x = self.quad_offset[0,quad]
qoff_y = self.quad_offset[1,quad]
self.image[qoff_x:qoff_x+850, qoff_y:qoff_y+850]=quad_image
return self.image
|
Fine tuning
Notice that the code snippets above make use of some predefined quantities which it reads in from "calibration files". The files contains calibrated numerical values for individual sections' and quads' rotations and shifts. All of these files are located in the experiment's 'calib' folder, but is not generated automatically. The XtcExplorer currently has a local version which is not correct but which is close enough to display a reasonable image. For how the files have been read in, you can take a look at XtcExplorer/src/cspad.py's read_alignment function.
For how to find the correct constants for each experiment, look at the CSPAD Alignment page.
Fine tuning
Notice that the code snippets above make use of some predefined quantities which it reads in from "calibration files". The files contains calibrated numerical values for individual sections' and quads' rotations and shifts. All of these files are located in the experiment's 'calib' folder, but is not generated automatically. The XtcExplorer currently has a local version which is not correct but which is close enough to display a reasonable image. For how the files have been read in, you can take a look at XtcExplorer/src/cspad.py's read_alignment function.
For how to find the correct constants for each experiment, look at the CSPad alignment page.
Saving data arrays
saving numpy arrays to numpy file
| Code Block |
|---|
import numpy as np
myarray = np.arange(0,100)
np.save( "output_file_name.npy", myarray)
np.savetxt( "output_file_name.txt", myarray)
|
Both of these work with arrays of maximum 2 dimensions. And only one array per file.
saving to MatLab file
| Code Block |
|---|
import scipy.io
N_array = np.arange(0,100)
x_array = np.random(100)
y_array = np.random(100)
scipy.io.savemat( "output_file_name.mat", mdict={'N': N_array, 'x' : x_array, 'y' : y_array } )
|
saving to HDF5 file
...
import h5py
ofile = h5py.File("output_file_name.h5",'w')
group = ofile.create_group("MyGroup")
group.create_dataset('delaytime',data=np.array(self.h_delaytime))
group.create_dataset('rawsignal',data=np.array(self.h_ipm_rsig))
group.create_dataset('normsignal',data=np.array(self.h_ipm_nsig))
ofile.close()
For more examples, see How to access HDF5 data from Python and http://code.google.com/p/h5py/
Interactive analysis with IPython
The version available in our offline release system is
IPython 0.9.1 – An enhanced Interactive Python.
so this is the one I've been using in these examples.
Not a whole lot more than a python shell.
However, the latest IPython has loads of new and interesting features...
| Panel | ||
|---|---|---|
| ||
|
Plotting with MatPlotLib
Matplotlib:
- The plotting can be done directly in the pyana module, but be aware that you need to disable plotting for the
module to run successfully in a batch job.Code Block import matplotlib.pyplot as plt plt.plot(array) plt.show()
- Or you can load arrays from a file and interactively plot them in iPython. The same ('recommended') syntax as above can be used, or if you use 'import *' you don't need to prepend the commands with the package name, which is handy when plotting interactively:
Code Block from matplotlib.pyplot import * ion() plot(array) draw()
Related useful tools and links
...
Non-interactive batch analysis
...
Overview
Content Tools