LSST Applications g0f08755f38+9c285cab97,g1635faa6d4+13f3999e92,g1653933729+a8ce1bb630,g1a0ca8cf93+bf6eb00ceb,g28da252d5a+0829b12dee,g29321ee8c0+5700dc9eac,g2bbee38e9b+9634bc57db,g2bc492864f+9634bc57db,g2cdde0e794+c2c89b37c4,g3156d2b45e+41e33cbcdc,g347aa1857d+9634bc57db,g35bb328faa+a8ce1bb630,g3a166c0a6a+9634bc57db,g3e281a1b8c+9f2c4e2fc3,g414038480c+077ccc18e7,g41af890bb2+fde0dd39b6,g5fbc88fb19+17cd334064,g781aacb6e4+a8ce1bb630,g80478fca09+55a9465950,g82479be7b0+d730eedb7d,g858d7b2824+9c285cab97,g9125e01d80+a8ce1bb630,g9726552aa6+10f999ec6a,ga5288a1d22+2a84bb7594,gacf8899fa4+c69c5206e8,gae0086650b+a8ce1bb630,gb58c049af0+d64f4d3760,gc28159a63d+9634bc57db,gcf0d15dbbd+4b7d09cae4,gda3e153d99+9c285cab97,gda6a2b7d83+4b7d09cae4,gdaeeff99f8+1711a396fd,ge2409df99d+5e831397f4,ge79ae78c31+9634bc57db,gf0baf85859+147a0692ba,gf3967379c6+41c94011de,gf3fb38a9a8+8f07a9901b,gfb92a5be7c+9c285cab97,w.2024.46
LSST Data Management Base Package
Loading...
Searching...
No Matches
Classes | Functions
lsst.scarlet.lite.fft Namespace Reference

Classes

class  Fourier
 

Functions

np.ndarray centered (np.ndarray arr, Sequence[int] newshape)
 
np.ndarray fast_zero_pad (np.ndarray arr, Sequence[Sequence[int]] pad_width)
 
np.ndarray _pad (np.ndarray arr, Sequence[int] newshape, int|Sequence[int]|None axes=None, str mode="constant", float constant_values=0)
 
tuple get_fft_shape (np.ndarray|Sequence[int] im_or_shape1, np.ndarray|Sequence[int] im_or_shape2, int padding=3, int|Sequence[int]|None axes=None, bool use_max=False)
 
Fourier _kspace_operation (Fourier image1, Fourier image2, int padding, Callable op, Sequence[int] shape, int|Sequence[int] axes)
 
Fourier|np.ndarray match_kernel (np.ndarray|Fourier kernel1, np.ndarray|Fourier kernel2, int padding=3, int|Sequence[int] axes=(-2, -1), bool return_fourier=True, bool normalize=False)
 
np.ndarray|Fourier convolve (np.ndarray|Fourier image, np.ndarray|Fourier kernel, int padding=3, int|Sequence[int] axes=(-2, -1), bool return_fourier=True, bool normalize=False)
 

Function Documentation

◆ _kspace_operation()

Fourier lsst.scarlet.lite.fft._kspace_operation ( Fourier image1,
Fourier image2,
int padding,
Callable op,
Sequence[int] shape,
int | Sequence[int] axes )
protected 
Combine two images in k-space using a given `operator`

Parameters
----------
image1:
    The LHS of the equation.
image2:
    The RHS of the equation.
padding:
    The amount of padding to add before transforming into k-space.
op:
    The operator used to combine the two images.
    This is either ``operator.mul`` for a convolution
    or ``operator.truediv`` for deconvolution.
shape:
    The shape of the output image.
axes:
    The dimension(s) of the array that will be transformed.

Definition at line 374 of file fft.py.

381) -> Fourier:
382 """Combine two images in k-space using a given `operator`
383
384 Parameters
385 ----------
386 image1:
387 The LHS of the equation.
388 image2:
389 The RHS of the equation.
390 padding:
391 The amount of padding to add before transforming into k-space.
392 op:
393 The operator used to combine the two images.
394 This is either ``operator.mul`` for a convolution
395 or ``operator.truediv`` for deconvolution.
396 shape:
397 The shape of the output image.
398 axes:
399 The dimension(s) of the array that will be transformed.
400 """
401 if len(image1.shape) != len(image2.shape):
402 msg = (
403 "Both images must have the same number of axes, "
404 f"got {len(image1.shape)} and {len(image2.shape)}"
405 )
406 raise ValueError(msg)
407
408 fft_shape = get_fft_shape(image1.image, image2.image, padding, axes)
409 if (
410 op == operator.truediv
411 or op == operator.floordiv
412 or op == operator.itruediv
413 or op == operator.ifloordiv
414 ):
415 # prevent divide by zero
416 lhs = image1.fft(fft_shape, axes)
417 rhs = image2.fft(fft_shape, axes)
418
419 # Broadcast, if necessary
420 if rhs.shape[0] == 1 and lhs.shape[0] != rhs.shape[0]:
421 rhs = np.tile(rhs, (lhs.shape[0],) + (1,) * len(rhs.shape[1:]))
422 if lhs.shape[0] == 1 and lhs.shape[0] != rhs.shape[0]:
423 lhs = np.tile(lhs, (rhs.shape[0],) + (1,) * len(lhs.shape[1:]))
424 # only select non-zero elements for the denominator
425 cuts = rhs != 0
426 transformed_fft = np.zeros(lhs.shape, dtype=lhs.dtype)
427 transformed_fft[cuts] = op(lhs[cuts], rhs[cuts])
428 else:
429 transformed_fft = op(image1.fft(fft_shape, axes), image2.fft(fft_shape, axes))
430 return Fourier.from_fft(transformed_fft, fft_shape, shape, axes, image1.image.dtype)
431
432

◆ _pad()

np.ndarray lsst.scarlet.lite.fft._pad ( np.ndarray arr,
Sequence[int] newshape,
int | Sequence[int] | None axes = None,
str mode = "constant",
float constant_values = 0 )
protected 
Pad an array to fit into newshape

Pad `arr` with zeros to fit into newshape,
which uses the `np.fft.fftshift` convention of moving
the center pixel of `arr` (if `arr.shape` is odd) to
the center-right pixel in an even shaped `newshape`.

Parameters
----------
arr:
    The arrray to pad.
newshape:
    The new shape of the array.
axes:
    The axes that are being reshaped.
mode:
    The numpy mode used to pad the array.
    In other words, how to fill the new padded elements.
    See ``numpy.pad`` for details.
constant_values:
    If `mode` == "constant" then this is the value to set all of
    the new padded elements to.

Definition at line 98 of file fft.py.

104) -> np.ndarray:
105 """Pad an array to fit into newshape
106
107 Pad `arr` with zeros to fit into newshape,
108 which uses the `np.fft.fftshift` convention of moving
109 the center pixel of `arr` (if `arr.shape` is odd) to
110 the center-right pixel in an even shaped `newshape`.
111
112 Parameters
113 ----------
114 arr:
115 The arrray to pad.
116 newshape:
117 The new shape of the array.
118 axes:
119 The axes that are being reshaped.
120 mode:
121 The numpy mode used to pad the array.
122 In other words, how to fill the new padded elements.
123 See ``numpy.pad`` for details.
124 constant_values:
125 If `mode` == "constant" then this is the value to set all of
126 the new padded elements to.
127 """
128 _newshape = np.asarray(newshape)
129 if axes is None:
130 currshape = np.array(arr.shape)
131 diff = _newshape - currshape
132 startind = (diff + 1) // 2
133 endind = diff - startind
134 pad_width = list(zip(startind, endind))
135 else:
136 # only pad the axes that will be transformed
137 pad_width = [(0, 0) for _ in arr.shape]
138 if isinstance(axes, int):
139 axes = [axes]
140 for a, axis in enumerate(axes):
141 diff = _newshape[a] - arr.shape[axis]
142 startind = (diff + 1) // 2
143 endind = diff - startind
144 pad_width[axis] = (startind, endind)
145 if mode == "constant" and constant_values == 0:
146 result = fast_zero_pad(arr, pad_width)
147 else:
148 result = np.pad(arr, tuple(pad_width), mode=mode) # type: ignore
149 return result
150
151

◆ centered()

np.ndarray lsst.scarlet.lite.fft.centered ( np.ndarray arr,
Sequence[int] newshape ) 
Return the central newshape portion of the array.

Parameters
----------
arr:
    The array to center.
newshape:
    The new shape of the array.

Notes
-----
If the array shape is odd and the target is even,
the center of `arr` is shifted to the center-right
pixel position.
This is slightly different than the scipy implementation,
which uses the center-left pixel for the array center.
The reason for the difference is that we have
adopted the convention of `np.fft.fftshift` in order
to make sure that changing back and forth from
fft standard order (0 frequency and position is
in the bottom left) to 0 position in the center.

Definition at line 34 of file fft.py.

34def centered(arr: np.ndarray, newshape: Sequence[int]) -> np.ndarray:
35 """Return the central newshape portion of the array.
36
37 Parameters
38 ----------
39 arr:
40 The array to center.
41 newshape:
42 The new shape of the array.
43
44 Notes
45 -----
46 If the array shape is odd and the target is even,
47 the center of `arr` is shifted to the center-right
48 pixel position.
49 This is slightly different than the scipy implementation,
50 which uses the center-left pixel for the array center.
51 The reason for the difference is that we have
52 adopted the convention of `np.fft.fftshift` in order
53 to make sure that changing back and forth from
54 fft standard order (0 frequency and position is
55 in the bottom left) to 0 position in the center.
56 """
57 _newshape = np.array(newshape)
58 currshape = np.array(arr.shape)
59
60 if not np.all(_newshape <= currshape):
61 msg = f"arr must be larger than newshape in both dimensions, received {arr.shape}, and {_newshape}"
62 raise ValueError(msg)
63
64 startind = (currshape - _newshape + 1) // 2
65 endind = startind + _newshape
66 myslice = [slice(startind[k], endind[k]) for k in range(len(endind))]
67
68 return arr[tuple(myslice)]
69
70

◆ convolve()

np.ndarray | Fourier lsst.scarlet.lite.fft.convolve ( np.ndarray | Fourier image,
np.ndarray | Fourier kernel,
int padding = 3,
int | Sequence[int] axes = (-2, -1),
bool return_fourier = True,
bool normalize = False ) 
Convolve image with a kernel

Parameters
----------
image:
    Image either as array or as `Fourier` object
kernel:
    Convolution kernel either as array or as `Fourier` object
padding:
    Additional padding to use when generating the FFT
    to suppress artifacts.
axes:
    Axes that contain the spatial information for the PSFs.
return_fourier:
    Whether to return `Fourier` or array
normalize:
    Whether or not to normalize the input kernels.

Returns
-------
result:
    The convolution of the image with the kernel.

Definition at line 481 of file fft.py.

488) -> np.ndarray | Fourier:
489 """Convolve image with a kernel
490
491 Parameters
492 ----------
493 image:
494 Image either as array or as `Fourier` object
495 kernel:
496 Convolution kernel either as array or as `Fourier` object
497 padding:
498 Additional padding to use when generating the FFT
499 to suppress artifacts.
500 axes:
501 Axes that contain the spatial information for the PSFs.
502 return_fourier:
503 Whether to return `Fourier` or array
504 normalize:
505 Whether or not to normalize the input kernels.
506
507 Returns
508 -------
509 result:
510 The convolution of the image with the kernel.
511 """
512 if not isinstance(image, Fourier):
513 image = Fourier(image)
514 if not isinstance(kernel, Fourier):
515 kernel = Fourier(kernel)
516
517 convolved = _kspace_operation(image, kernel, padding, operator.mul, image.shape, axes=axes)
518 if return_fourier:
519 return convolved
520 else:
521 return np.real(convolved.image)

◆ fast_zero_pad()

np.ndarray lsst.scarlet.lite.fft.fast_zero_pad ( np.ndarray arr,
Sequence[Sequence[int]] pad_width ) 
Fast version of numpy.pad when `mode="constant"`

Executing `numpy.pad` with zeros is ~1000 times slower
because it doesn't make use of the `zeros` method for padding.

Parameters
---------
arr:
    The array to pad
pad_width:
    Number of values padded to the edges of each axis.
    See numpy.pad docs for more.

Returns
-------
result: np.ndarray
    The array padded with `constant_values`

Definition at line 71 of file fft.py.

71def fast_zero_pad(arr: np.ndarray, pad_width: Sequence[Sequence[int]]) -> np.ndarray:
72 """Fast version of numpy.pad when `mode="constant"`
73
74 Executing `numpy.pad` with zeros is ~1000 times slower
75 because it doesn't make use of the `zeros` method for padding.
76
77 Parameters
78 ---------
79 arr:
80 The array to pad
81 pad_width:
82 Number of values padded to the edges of each axis.
83 See numpy.pad docs for more.
84
85 Returns
86 -------
87 result: np.ndarray
88 The array padded with `constant_values`
89 """
90 newshape = tuple([a + ps[0] + ps[1] for a, ps in zip(arr.shape, pad_width)])
91
92 result = np.zeros(newshape, dtype=arr.dtype)
93 slices = tuple([slice(start, s - end) for s, (start, end) in zip(result.shape, pad_width)])
94 result[slices] = arr
95 return result
96
97

◆ get_fft_shape()

tuple lsst.scarlet.lite.fft.get_fft_shape ( np.ndarray | Sequence[int] im_or_shape1,
np.ndarray | Sequence[int] im_or_shape2,
int padding = 3,
int | Sequence[int] | None axes = None,
bool use_max = False ) 
Return the fast fft shapes for each spatial axis

Calculate the fast fft shape for each dimension in
axes.

Parameters
----------
im_or_shape1:
    The left image or shape of an image.
im_or_shape2:
    The right image or shape of an image.
padding:
    Any additional padding to add to the final shape.
axes:
    The axes that are being transformed.
use_max:
    Whether or not to use the maximum of the two shapes,
    or the sum of the two shapes.

Returns
-------
shape:
    Tuple of the shape to use when the two images are transformed
    into k-space.

Definition at line 152 of file fft.py.

158) -> tuple:
159 """Return the fast fft shapes for each spatial axis
160
161 Calculate the fast fft shape for each dimension in
162 axes.
163
164 Parameters
165 ----------
166 im_or_shape1:
167 The left image or shape of an image.
168 im_or_shape2:
169 The right image or shape of an image.
170 padding:
171 Any additional padding to add to the final shape.
172 axes:
173 The axes that are being transformed.
174 use_max:
175 Whether or not to use the maximum of the two shapes,
176 or the sum of the two shapes.
177
178 Returns
179 -------
180 shape:
181 Tuple of the shape to use when the two images are transformed
182 into k-space.
183 """
184 if isinstance(im_or_shape1, np.ndarray):
185 shape1 = np.asarray(im_or_shape1.shape)
186 else:
187 shape1 = np.asarray(im_or_shape1)
188 if isinstance(im_or_shape2, np.ndarray):
189 shape2 = np.asarray(im_or_shape2.shape)
190 else:
191 shape2 = np.asarray(im_or_shape2)
192 # Make sure the shapes are the same size
193 if len(shape1) != len(shape2):
194 msg = (
195 "img1 and img2 must have the same number of dimensions, "
196 f"but got {len(shape1)} and {len(shape2)}"
197 )
198 raise ValueError(msg)
199 # Set the combined shape based on the total dimensions
200 if axes is None:
201 if use_max:
202 shape = np.max([shape1, shape2], axis=0)
203 else:
204 shape = shape1 + shape2
205 else:
206 if isinstance(axes, int):
207 axes = [axes]
208 shape = np.zeros(len(axes), dtype="int")
209 for n, ax in enumerate(axes):
210 shape[n] = shape1[ax] + shape2[ax]
211 if use_max:
212 shape[n] = np.max([shape1[ax], shape2[ax]])
213
214 shape += padding
215 # Use the next fastest shape in each dimension
216 shape = [fftpack.next_fast_len(s) for s in shape]
217 return tuple(shape)
218
219

◆ match_kernel()

Fourier | np.ndarray lsst.scarlet.lite.fft.match_kernel ( np.ndarray | Fourier kernel1,
np.ndarray | Fourier kernel2,
int padding = 3,
int | Sequence[int] axes = (-2, -1),
bool return_fourier = True,
bool normalize = False ) 
Calculate the difference kernel to match kernel1 to kernel2

Parameters
----------
kernel1:
    The first kernel, either as array or as `Fourier` object
kernel2:
    The second kernel, either as array or as `Fourier` object
padding:
    Additional padding to use when generating the FFT
    to supress artifacts.
axes:
    Axes that contain the spatial information for the kernels.
return_fourier:
    Whether to return `Fourier` or array
normalize:
    Whether or not to normalize the input kernels.

Returns
-------
result:
    The difference kernel to go from `kernel1` to `kernel2`.

Definition at line 433 of file fft.py.

440) -> Fourier | np.ndarray:
441 """Calculate the difference kernel to match kernel1 to kernel2
442
443 Parameters
444 ----------
445 kernel1:
446 The first kernel, either as array or as `Fourier` object
447 kernel2:
448 The second kernel, either as array or as `Fourier` object
449 padding:
450 Additional padding to use when generating the FFT
451 to supress artifacts.
452 axes:
453 Axes that contain the spatial information for the kernels.
454 return_fourier:
455 Whether to return `Fourier` or array
456 normalize:
457 Whether or not to normalize the input kernels.
458
459 Returns
460 -------
461 result:
462 The difference kernel to go from `kernel1` to `kernel2`.
463 """
464 if not isinstance(kernel1, Fourier):
465 kernel1 = Fourier(kernel1)
466 if not isinstance(kernel2, Fourier):
467 kernel2 = Fourier(kernel2)
468
469 if kernel1.shape[0] < kernel2.shape[0]:
470 shape = kernel2.shape
471 else:
472 shape = kernel1.shape
473
474 diff = _kspace_operation(kernel1, kernel2, padding, operator.truediv, shape, axes=axes)
475 if return_fourier:
476 return diff
477 else:
478 return np.real(diff.image)
479
480
Generated on Wed Nov 20 2024 09:02:54 for LSST Applications by doxygen 1.10.0