This note is about n-d array processing algorithms implemented in ImgAlgos.PyAlgos. Algorithms can be called from python but low level implementation is done on C++ with boost/python wrapper. All examples are shown for python level interface.

Content

Common features of algorithms

n-d arrays

LCLS detector data come from DAQ as n-d arrays (ndarray in C++ or numpy.array in Python). In simple case camera data is an image presented by the 2-d array. For composite detectors like CSPAD, CSPAD2X2, EPIX, PNCCD, etc. data comes from a set of sensors as 3-d or 4-d arrays. If relative sensors' positions are known, then sensors can be composed in 2-d image. But this image contains significant portion of "fake" empty pixels, that may be up to ~20-25% in case of CSPAD. Most efficient data processing algorithms should be able to work with n-d arrays.

Windows

In some experiments not all sensors contain useful data. It might be more efficient to select Region of Interest (ROI) on sensors, where data need to be processed. To support this feature a tuple (or list) of windows is passed as a constructor parameter. Each window is presented by the tuple of 5 parameters (segnum, rowmin, rowmax, colmin, colmax), where segnum is a sensor index in the n-d array, other parameters constrain window 2-d matrix rows and columns. Several windows can be defined for the same sensor using the same segnum. For 2-d arrays segnumparameter is not used, but still needs to be presented in the window tuple by any integer number. To increase algorithm efficiency only pixels in windows are processed. If windows=None, all sensors will be processed.

The array of windows can be converted in 3-d or 2-d array of mask using method pyimgalgos.GlobalUtils.mask_from_windows.

Mask

Alternatively ROI can be defined by the mask of good/bad (1/0) pixels. For 2-d image mask can easily be defined in user's code. In case of ≥3-d arrays the Mask Editor helps to produce ROI mask. Entire procedure includes

...

In addition mask accounts for bad pixels which should be discarded in processing. Total mask may be a product of ROI and other masks representing good/bad pixels.

Make object and set parameters

Any algorithm object can be created as shown below.

import numpy as np
from ImgAlgos.PyAlgos import PyAlgos 

# create object:
alg = PyAlgos(windows=winds, mask=mask, pbits=0)

Define ROI using windows and/or mask

Region Of Interest (ROI)is defined by the set of rectangular windows on segments and mask, as shown in example below.

# List of windows
winds = None # entire size of all segments will be used for peak finding
winds = (( 0, 0, 185, 0, 388),              
         ( 1, 20,160, 30,300),
         ( 7, 0, 185, 0, 388))

# Mask
mask = None                   # (default) all pixels in windows will be used for peak finding
mask = det.mask()             # see class Detector.PyDetector
mask = np.loadtxt(fname_mask) # 
mask.shape = <should be the same as shape of data n-d array>

Hit finders

Hit finders return simple values for decision on event selection. Two algorithms are implemented in ImgAlgos.PyAlgos. They count number of pixels and intensity above threshold in the Region Of Interest (ROI) defined by windows and mask parameters in object constructor.

Both hit-finders receive input n-d drray data and threshold thr parameters and return a single value in accordance with method name.

Number of pixels above threshold

number_of_pix_above_thr

npix = alg.number_of_pix_above_thr(data, thr=10)

Total intensity above threshold

intensity_of_pix_above_thr

intensity = alg.intensity_of_pix_above_thr(data, thr=12)

Peak finders

Peak selection parameters

Internal peak selection is done at the end of each peak finder, but all peak selection parameters need to be defined right after algorithm object is created.  These peak selection parameters are set for all peak-finders:

...

• nda - calibrated n-d array of data, pedestals and background should be subtracted, common mode - corrected
• thr* - different type of thresholds
• rank - peak rank as explained in  section below.
• r0, dr - ring internal radius and width to evaluate background and noise rms as explained in  section below.

Two threshold "Droplet finder"

two-threshold peak-finding algorithm in restricted region around pixel with maximal intensity. Two threshold allows to speed-up this algorithms. It is assumed that only pixels with intensity above thr_high are pretending to be peak candidate centers. Candidates are considered as a peak if their intensity is maximal in the (square) region of radius around them. Low threshold in the same region is used to account for contributing to peak pixels.

peak_finder_v1

peaks = alg.peak_finder_v1(nda, thr_low=10, thr_high=150, radius=5, dr=0.05)

...

• defines (square) region to search for local maximum with intensity above thr_high and contributing pixels with intensity above  thr_lo,
• is used as a r0 parameter to evaluate background and noise rms as explained in section below.

peak_finder_v4

peaks = alg.peak_finder_v4(nda, thr_low=10, thr_high=150, rank=4, r0=5, dr=0.05)

The same algorithm as  peak_finder_v1 , but parameter radius is split for two (unsigned) rank and (float)r0 with the same meaning as in  peak_finder_v3 .

Flood filling algorithm

define peaks for regions of connected pixels above threshold

peak_finder_v2

peaks = alg.peak_finder_v2(nda, thr=10, r0=5, dr=0.05)

Two neighbor pixels are assumed connected if have common side. Pixels with intensity above threshold thr are considered only.

Local maximums search algorithm

define peaks in local maximums of specified rank (radius), for example rank=2 means 5x5 pixel region around central pixel.

peak_finder_v3

peaks = alg.peak_finder_v3(nda, rank=2, r0=5, dr=0.05)

...

Demonstration for local maximum map

Test for 100x100 image with random normal distribution of intensities

...

Anchor
bkg_rms bkg_rms

Evaluation of the background level, rms, and S/N ratio

When peak is found, its parameters can be precised for background level, noise rms, and signal over background ratio (S/N) can be estimated. All these values can be evaluated using pixels surrounding the peak on some distance. For all peak-finders we use the same algorithm. Surrounding pixels are defined by the ring with internal radial parameter r0 and ring width dr (both in pixels). The number of surrounding pixels depends on r0 and dr parameters as shown in matrices below. We use notation

• + central pixel with maximal intensity,
• 1 pixels counted in calculation of averaged background level and noise rms,
• 0 pixels not counted.

Matrices of pixels for r0=3 and 4 and different dr values

r0=3  dr=0.1 (4 pixels)   r0=3  dr=0.5 (12 pixels)   r0=3  dr=1 (24 pixels) 
0 0 0 0 0 0 0 0 0         0 0 0 0 0 0 0 0 0          0 0 0 0 1 0 0 0 0       
0 0 0 0 1 0 0 0 0         0 0 0 1 1 1 0 0 0          0 0 1 1 1 1 1 0 0       
0 0 0 0 0 0 0 0 0         0 0 0 0 0 0 0 0 0          0 1 0 0 0 0 0 1 0       
0 0 0 0 0 0 0 0 0         0 1 0 0 0 0 0 1 0          0 1 0 0 0 0 0 1 0       
0 1 0 0 + 0 0 1 0         0 1 0 0 + 0 0 1 0          1 1 0 0 + 0 0 1 1       
0 0 0 0 0 0 0 0 0         0 1 0 0 0 0 0 1 0          0 1 0 0 0 0 0 1 0       
0 0 0 0 0 0 0 0 0         0 0 0 0 0 0 0 0 0          0 1 0 0 0 0 0 1 0       
0 0 0 0 1 0 0 0 0         0 0 0 1 1 1 0 0 0          0 0 1 1 1 1 1 0 0       
0 0 0 0 0 0 0 0 0         0 0 0 0 0 0 0 0 0          0 0 0 0 1 0 0 0 0

r0=4  dr=0.2 (12 pixels)   r0=4  dr=0.3 (16 pixels)    r0=4  dr=0.5 (24 pixels)
0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0       0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 1 1 0 0 0 0      0 0 0 0 1 1 1 0 0 0 0       0 0 0 1 1 1 1 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0      0 0 1 0 0 0 0 0 1 0 0       0 0 1 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0       0 1 0 0 0 0 0 0 0 1 0
0 1 0 0 0 0 0 0 0 1 0      0 1 0 0 0 0 0 0 0 1 0       0 1 0 0 0 0 0 0 0 1 0
0 1 0 0 0 + 0 0 0 1 0      0 1 0 0 0 + 0 0 0 1 0       0 1 0 0 0 + 0 0 0 1 0
0 1 0 0 0 0 0 0 0 1 0      0 1 0 0 0 0 0 0 0 1 0       0 1 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0       0 1 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0      0 0 1 0 0 0 0 0 1 0 0       0 0 1 0 0 0 0 0 1 0 0
0 0 0 0 1 1 1 0 0 0 0      0 0 0 0 1 1 1 0 0 0 0       0 0 0 1 1 1 1 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0       0 0 0 0 0 0 0 0 0 0 0 

Matrices of pixels for r0=5 and 6 and different dr values

r0=5  dr=0.05 (12 pixels)      r0=5  dr=0.5  (28 pixels)            
0 0 0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 0 0 0 0 0 0      0 0 0 0 1 1 1 1 1 0 0 0 0
0 0 0 1 0 0 0 0 0 1 0 0 0      0 0 0 1 0 0 0 0 0 1 0 0 0
0 0 1 0 0 0 0 0 0 0 1 0 0      0 0 1 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0      0 1 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0      0 1 0 0 0 0 0 0 0 0 0 1 0
0 1 0 0 0 0 + 0 0 0 0 1 0      0 1 0 0 0 0 + 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0      0 1 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0      0 1 0 0 0 0 0 0 0 0 0 1 0
0 0 1 0 0 0 0 0 0 0 1 0 0      0 0 1 0 0 0 0 0 0 0 1 0 0
0 0 0 1 0 0 0 0 0 1 0 0 0      0 0 0 1 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 1 0 0 0 0 0 0      0 0 0 0 1 1 1 1 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0      0 0 0 0 0 0 0 0 0 0 0 0 0

r0=6  dr=0.2 (12 pixels)        r0=6  dr=0.5 (28 pixels)      
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 1 1 1 0 0 0 0 0 0     0 0 0 0 0 1 1 1 1 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 1 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 1 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 1 0 0 0 0 0 0 0 0 0 0 0 1 0
0 1 0 0 0 0 0 0 0 0 0 0 0 1 0     0 1 0 0 0 0 0 0 0 0 0 0 0 1 0
0 1 0 0 0 0 0 + 0 0 0 0 0 1 0     0 1 0 0 0 0 0 + 0 0 0 0 0 1 0
0 1 0 0 0 0 0 0 0 0 0 0 0 1 0     0 1 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 1 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 1 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 1 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 1 1 1 0 0 0 0 0 0     0 0 0 0 0 1 1 1 1 1 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0     0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

References

• ImgAlgos.PyAlgos - code example in Sphinx documentation
• Peak Finding - short announcement about peak finders
• Hit and Peak Finders - examples in Chris' tutorial
• Peak Finding Module - (depricated) psana module, it demonstaration examples and results
• Psana Module Catalog - (depricated) peak finding  psana modules
• Psana Module Examples - (depricated) peak finding examples in psana modules
• GUI for tuning peak finding - Chun's page in development
  
  

...