2.1.

Forward Model

The forward model of CUP formulates the process of recording a three-dimensional [3D; i.e., $(x,y,t)$] scene to one or a few 2D snapshots. In general, this forward model can be expressed mathematically using either element-wise or matrix-vector notations.

2.1.1.

Element-wise notation

Many mathematical and scientific libraries are designed to efficiently handle element-wise operations. The dynamic scene and the binary-valued encoding mask are denoted by $\mathbf{F}\in {\mathbb{R}}^{M\times N\times L}$ and $\mathbf{R}\in {\mathbb{R}}^{M\times N}$, respectively. $M$ and $N$ represent the data lengths in the two spatial dimensions, and $L$ signifies the data length in time. The discrete output from the sensor for the time-sheared view (hence the subscript “ts”) can be modeled as

Eq. (1)

$${(G_{\mathrm{ts}})}_{i,j}=\sum _{l=0}^{L-1}{\overline{F}}_{i,j,l}{\overline{R}}_{i,j}+(E_{\mathrm{ts}}{)}_{i,j},$$

where

$${\overline{F}}_{i,j,l}=\{\begin{array}{ll}{F}_{i,j-l,l}& \text{if}l\le j\le [N+(l-1)]\\ 0& \text{otherwise},\end{array}$$

and

$${\overline{R}}_{i,j}=\{\begin{array}{ll}{R}_{i,j-l}& \text{if}l\le j\le [N+(l-1)]\\ 0& \text{otherwise},\end{array}$$

where ${\overline{F}}_{i,j,l}$ is the intensity of the $(i,j,l)$’th element of a right-zero-padded version with a frame-dependent right circular shifting of the dynamic scene’s datacube with $\overline{\mathbf{F}}\in {\mathbb{R}}^{M\times [N+(L-1)]\times L}$. ${\overline{R}}_{i,j}$ stands for the intensity of the $(i,j)$’th element a right-zero-padded version with a frame-dependent right circular shifting of the spatial encoder with $\overline{\mathbf{R}}\in {\mathbb{R}}^{M\times [N+(L-1)]}$. $(G_{\mathrm{ts}}{)}_{i,j}$ is the intensity measured on the $(i,j)$’th element of the sensor with ${\mathbf{G}}_{\mathrm{ts}}\in {\mathbb{R}}^{M\times [N+(L-1)]}$. ${(E_{\mathrm{ts}})}_{i,j}$ stands for the noise of the $(i,j)$’th element in ${\mathbf{G}}_{\mathrm{ts}}$ with ${\mathbf{E}}_{\mathrm{ts}}\in {\mathbb{R}}^{M\times [N+(L-1)]}$.

The discrete output for the time-unsheared view (hence “tu” as the subscript) can be modeled as

Eq. (2)

$${(G_{\mathrm{tu}})}_{i,j}=\sum _{l=0}^{L-1}{F}_{i,j,l}+{(E_{\mathrm{tu}})}_{i,j},$$

where ${(G_{\mathrm{tu}})}_{i,j}$ is the intensity measured at the $(i,j)$’th element of the time-unsheared view with ${\mathbf{G}}_{\mathrm{tu}}\in {\mathbb{R}}^{M\times N}$, and ${(E_{\mathrm{tu}})}_{i,j}$ represents the noise in ${\mathbf{G}}_{\mathrm{tu}}$ with ${\mathbf{E}}_{\mathrm{tu}}\in {\mathbb{R}}^{M\times N}$.

As an example, a Matlab script that simulates dual-view CUP’s forward model, with a linear shearing operator and a pseudorandom binary mask, is shown in Algorithm 1. Moreover, a step-by-step guide with illustrations of the “cell-division” dynamic scene is shown in Fig. 2. The ground truth video was taken from the public “Mouse Embryo Tracking” database⁶⁷ and can be downloaded using the link in Ref. 68.

Algorithm 1

Simulating dual-view CUP’s forward model with the element-wise notation using Matlab.

%% Example of dual-view CUP’s sensing process
% Encoding step (generating R)
load(‘Cell.mat’)	% Loading the example video – F
[M,N,L] = size(F);	% Calculating the video dimensions % (y,x,t)->(M,N,L)
R = 1*(rand(M,N)>0.5);	% Mask initialization with a % transmittance of ∼50%
R = repmat(R,1,1,L);
Gts = F.* R;	% Spatial encoding – Hadamard product
Gts = padarray(Gts,[0,L-1,0],0,‘post’);	% Right column zero padding
for l=0:L-1	% Shearing operation
Gts(:,:,l+1) =circshift(Gts(:,:,l+1),[0,l]);
end
Gts =sum(Gts,3);	% Integration of time-sheared view
Gtu =sum(F,3);	% Integration of time-unsheared view

Fig. 2

Illustrations of simulating dual-view CUP’s forward model with the element-wise notation using Matlab.

2.1.2.

Matrix-vector notation

The element-wise notation of CUP’s forward model, despite owning simple expression and easy comprehension, is inherently limited by its sequential execution. This characteristic engenders surplus calculations within specific functions when the modeling is subjected to extensive datasets or algorithms mandating intricate computations—such as matrix inversion, matrix factorization, eigenvalue decomposition, and low-rank approximation. Thus, most practices use matrix-vector operations by converting $\mathbf{F}$ into vector $\mathbf{f}\in {\mathbb{R}}^{n\times 1}$, where $n=M·N·L$. In this way, CUP’s forward model can be expressed by matrix multiplication, which can be computed by powerful linear algebra methods to concisely formulate solutions.

In particular, dual-view CUP’s forward model, by following a matrix-vector representation, can be expressed as

Eq. (3)

$$\begin{array}{c}\mathbf{g}=\left[\begin{array}{c}{\mathbf{g}}_{\mathrm{ts}}\\ {\mathbf{g}}_{\mathrm{tu}}\end{array}\right]=\mathbf{\Phi }\mathbf{f}=\left[\begin{array}{c}{\mathbf{T}}_{\mathrm{ts}}\mathbf{SC}\\ {\mathbf{T}}_{\mathrm{tu}}\end{array}\right]\mathbf{f},\end{array}$$

where, ${\mathbf{g}}_{\mathrm{ts}}\in {\mathbb{R}}^{m_{\mathrm{ts}}\times 1}$ and ${\mathbf{g}}_{\mathrm{tu}}\in {\mathbb{R}}^{m_{\mathrm{tu}}\times 1}$ are the vectorial version of the time-sheared view and the time-unsheared view with sizes $m_{\mathrm{ts}}=M·[N+(L-1)]$ and $m_{\mathrm{tu}}=M·N$, respectively. $\mathbf{g}\in {\mathbb{R}}^{m\times 1}$ is the vectorial version of the concatenated two views with a size $m=m_{\mathrm{ts}}+m_{\mathrm{tu}}$. $\mathbf{\Phi }\in {\mathbb{R}}^{m\times n}$ is the dual-view CUP’s sensing matrix. $\mathbf{C}\in {\mathbb{R}}^{n\times n}$ is the spatial encoding matrix. $\mathbf{S}\in {\mathbb{R}}^{(m_{\mathrm{ts}}·L)\times n}$ is the temporal shearing matrix. ${\mathbf{T}}_{\mathrm{ts}}\in {\mathbb{R}}^{m_{\mathrm{ts}}\times (m_{\mathrm{ts}}·L)}$ and ${\mathbf{T}}_{\mathrm{tu}}\in {\mathbb{R}}^{m_{\mathrm{tu}}\times n}$ are the spatiotemporal integration matrices of the time-sheared view and the time-unsheared view, respectively. ${\mathbf{\Phi }}_{\mathrm{ts}}={\mathbf{T}}_{\mathrm{ts}}\mathbf{SC}\in {\mathbb{R}}^{m_{\mathrm{ts}}\times n}$ is also defined as the time-sheared sensing matrix.

The entries of $\mathbf{C}$, $\mathbf{S}$, ${\mathbf{T}}_{\mathrm{ts}}$, and ${\mathbf{T}}_{\mathrm{tu}}$ are given as

Eq. (4)

$$C_{i,j}=\{\begin{array}{ll}{r}_v& \text{if}i=j\\ 0& \text{otherwise}\end{array},$$

Eq. (5)

$$S_{i,j}=\{\begin{array}{ll}1& \text{if}i=j+M·L·\lfloor \frac{j}{M·N}\rfloor \\ 0& \text{otherwise},\end{array}$$

Eq. (6)

$${\mathbf{T}}_{\mathrm{ts}}={\mathbf{1}}_L^T\otimes {\mathbf{I}}_{m_{\mathrm{ts}}\times m_{\mathrm{ts}},}$$

Eq. (7)

$${\mathbf{T}}_{\mathrm{tu}}={\mathbf{1}}_L^T\otimes {\mathbf{I}}_{m_{\mathrm{tu}}\times m_{\mathrm{tu}}}.$$

In Eq. (4), $v=\mathrm{mod}(j,M·N)$ with $v\in \mathbb{W}$, and $r_v\in \{\mathrm{0,1}\}$ is the value at the $v$’th position of $\mathbf{r}\in {\mathbb{R}}^{M·N\times 1}$, which is the vectorial version of the encoding mask $\mathbf{R}$. In Eqs. (6) and (7), ${\mathbf{1}}_L\in {\mathbb{R}}^{L\times 1}$ is an all-one vector. $\mathbf{I}$ is the identity matrix. The matrix ${\mathbf{T}}_{\mathrm{tu}}$ has a structure similar to ${\mathbf{T}}_{\mathrm{ts}}$ (i.e., a horizontal concatenation of identity matrices) but with a shorter diagonal dimension. Thus, the sensing matrix $\mathbf{\Phi }=\left[\begin{array}{c}{\mathbf{\Phi }}_{\mathrm{ts}}\\ {\mathbf{T}}_{\mathrm{tu}}\end{array}\right]$ can be directly defined as

Eq. (8)

$${\mathrm{\Phi }}_{i,j}=\{\begin{array}{ll}{r}_v& \mathrm{if}i=v+M·\lfloor \frac{j}{M·N}\rfloor \\ 1& \text{if}i=v+m_{\mathrm{ts}}\\ 0& \text{otherwise}\end{array}.$$

A Matlab script for constructing the matrices $\mathbf{C}$, $\mathbf{S}$, ${\mathbf{T}}_{\mathrm{ts}}$, ${\mathbf{T}}_{\mathrm{tu}}$, and $\mathbf{\Phi }$ is presented in Algorithm 2. In this example, $M\times N\times L=256\times 256\times 25\text{pixels}$. The sensing matrix of the time-sheared view ${\mathbf{\Phi }}_{\mathrm{ts}}$ is created by assembling a series of diagonal patterns that cyclically repeat along the horizontal direction, shifting downward by $M$ rows following each iteration. The sensing matrix $\mathbf{\Phi }$, illustrated schematically in Fig. 3, is created by vertically concatenating ${\mathbf{\Phi }}_{\mathrm{ts}}$ with ${\mathbf{T}}_{\mathrm{tu}}$.

Algorithm 2

Programming matrices of dual-view CUP’s operation using Matlab.

clear all; close all; clc

M = 256; N = 256; L = 25;

R = 1*(rand(M,N)>0.5);

R = repmat(R,1,1,L);

%%

n = M*N*L;

mts = M*(N+L-1);

mtu = M*N;

%% Encoding matrix (C)

i = 0:n-1;

j = 0:n-1;

C = sparse(i+1,j+1,R(:),n,n);

%% Shearing matrix (S)

j = 0:n-1;

i = j+(M*L)*floor(j/(M*N));

S = sparse(i+1,j+1,1,mts*L,n);

%% Integration matrices (Ts and Tu)

Tts = kron(ones(1,L),speye(mts,mts));

Ttu = kron(ones(1,L),speye(mtu,mtu));

%% CUP sensing matrix (\Phi)

Phi = [((Tts*S)*C)’,Ttu’]’;

Fig. 3

Construction of the matrices of dual-view CUP’s operation using Matlab. (a) Spatial encoding matrix. (b) Temporal shearing matrix. (c) Spatiotemporal integration matrix for the time-sheared view. (d) Spatiotemporal integration matrix for the time-unsheared view. (e) Sensing matrix of dual-view CUP. Insets: Zoomed-in views of local regions (indicated by the red boxes with different line types).

2.2.

Image Reconstruction

After data acquisition, the captured 2D snapshots (i.e., ${\mathbf{G}}_{\mathrm{ts}}$ and ${\mathbf{G}}_{\mathrm{tu}}$) are input to an algorithm to reconstruct the dynamic scene. To date, analytical-modeling-based algorithms are dominantly used in CUP’s reconstruction because they can incorporate prior knowledge about the imaging system and the underlying physics of light propagation,^69–74 leading to accurate reconstructions. Before getting into sophisticated analytical-modeling-based reconstruction algorithms for CUP, let us analyze the structure of a basic optimization problem:

Various reconstruction algorithms^75–78 are developed based on Eq. (9). A popular choice, especially in the early stage of CUP’s development, is the two-step iterative shrinkage/thresholding (TwIST) algorithm.⁷⁸ The regularizer, $\phi (\mathbf{f})$, can be set to various forms, including $\Vert {\mathbf{\Psi }\mathbf{f}\Vert }_1$ and ${\Vert \mathbf{f}\Vert }_{\mathrm{TV}}$, where ${\Vert ·\Vert }_1$ represents the ${\ell }_1$-norm, $\mathbf{\Psi }\in {\mathbb{R}}^{n\times n}$ is an arbitrary representation basis matrix, and ${\Vert ·\Vert }_{\mathrm{TV}}$ represents the total-variation (TV) regularization.⁷⁹ The TwIST algorithm combines the shrinkage operation used in iterative soft-thresholding algorithms with a correction step that enforces fidelity to the measurements. It exploits the sparsity naturally embedded in the transient scene via the regularizer. In particular, the ${\ell }_1$-norm requires prior knowledge about the scene to select an adequate representation basis. The TV-norm uses spatiotemporal correlation by removing low variations between neighboring pixels. These characteristics enable the TwIST algorithm to efficiently recover a transient scene from an underdetermined measurement.

Later, more advanced image reconstruction algorithms are developed based on the paradigm of alternating direction method of multipliers (ADMM),^69,80–83 which sets $\phi (\mathbf{f})={\Vert \mathbf{f}-{\mathbf{f}}^{(k)}\Vert }_2^2$, where ${\mathbf{f}}^{(k)}$ stands for the $k$’th reconstruction with $k=\{0,\dots ,K-1\}$ and $K\in \mathbb{N}$ as the total number of the algorithm’s iterations. The ADMM has gained increasing popularity due to its flexibility to customize the optimization steps by incorporating additional constraints (e.g., noise reduction algorithms or neural-network approaches), which is hence selected for this tutorial.

2.2.2.

Deep-learning approaches

Deep-learning approaches have been increasingly featured owing to their faster reconstruction compared with their analytical-modeling-based counterparts. Recent advances have allowed embedding mathematical properties offered by the CS theory by designing custom layers that emulate the forward sensing model or exploiting spatiotemporal sparsity via image-denoising nets.⁹³ Given access to rich available training datasets, many novel methods based on convolutional neural networks (CNNs) have been developed for CUP’s reconstruction as well as for the encoding mask design, including the end-to-end CNN with residual learning,⁹⁴ the U-Net-based DeepCUP,⁹⁵ the hybrid algorithm that combines the AL method with deep learning,⁹⁶ and the snapshot-to-video autoencoder based on a generative adversarial network.^97–103

Here, we review a representative CNN—the deep high-dimensional adaptive net (D-HAN)¹⁰⁴ that offers multifaceted supervision to CUP by optimizing the encoding mask, sensing the shearing operation, and reconstructing the 3D datacubes. The main goal of the D-HAN is to leverage the merits of both the ADMM and the network-based CS methods by mapping one iteration of the ADMM steps to a deep network architecture. For these reasons, the D-HAN will be used as a benchmark to explain how to link the CUP’s forward model with a CNN approach.

Originally designed to use only the time-sheared view, the D-HAN is composed of two cascaded neural networks: a deep-unfolding-based network to embody the sensing model of the time-sheared view in CUP and a U-Net architecture¹⁰⁵ to further improve image reconstruction (Fig. 4) by exploiting the spatiotemporal correlation of the transient scene. The deep-unfolding net and the U-Net manifest the “divide-and-conquer” approach embedded in the ADMM. Then, the time-unsheared view was incorporated to boost the reconstruction performance by using it as an initialization for the deep-unfolding network and a prior restriction in the loss function. This configuration leverages the original D-HAN’s mathematical advantages and the reduction of unknowns via prior information. This design is memory efficient and thus essential for learning to reconstruct high-dimensional datacubes.

Fig. 4

Schematic of the D-HAN for dual-view CUP’s image reconstruction. BN, batch normalization; ReLU, rectified linear activation unit. Adapted with permission from Ref. 104.

In this regard, the ADMM-based inverse problem can be formulated using the ADMM’s scaled form [i.e., Eqs. (19)–(21)]. Note that in Eq. (19), the analytical inverse model of $\mathbf{f}$ refers to a quadratic problem with the closed-form that involves the inversion of a $n\times n$ size matrix [see Eq. (22)]. Toward this goal, the Sherman–Woodbury–Morrison (SWM) matrix inversion lemma¹⁰⁶—a mathematical theorem allowing calculating a matrix’s inverse by converting it into a full rank matrix—and the full-column rank properties are exploited to simplify the process to a smaller-scale matrix inversion and obtain the closed-form solution of the first inverse model in Eq. (22)

Eq. (25)

$$\mathbf{f}={\tilde{\rho }}^{-1}[\mathbf{I}-{\mathbf{\Phi }}_{\mathrm{ts}}^T{[\tilde{\rho }\mathbf{I}+{\mathbf{\Phi }}_{\mathrm{ts}}{\mathbf{\Phi }}_{\mathrm{ts}}^T]}^{-1}{\mathbf{\Phi }}_{\mathrm{ts}}][{\mathbf{\Phi }}_{\mathrm{ts}}^T{\mathbf{g}}_{\mathrm{ts}}+\tilde{\rho }(\mathbf{z}-\mathbf{w})],$$

where ${\mathbf{\Phi }}_{\mathrm{ts}}{\mathbf{\Phi }}_{\mathrm{ts}}^T\in {\mathbb{R}}^{m_{\mathrm{ts}}\times m_{\mathrm{ts}}}$ represents a matrix product resulting in a diagonal matrix, and $\tilde{\rho }=\rho /2$.

To implement the D-HAN, the first step is to define the operators of dual-view CUP’s data acquisition. First, the direct sensing operators of the time-sheared view and the time-unsheared view, denoted by ${\mathcal{G}}_{\mathrm{ts}}$ and ${\mathcal{G}}_{\mathrm{tu}}$, are expressed as

Eq. (26)

$${\mathcal{G}}_{\mathrm{ts}}(\mathbf{F},\mathbf{R})=\sum _{l=0}^{L-1}\mathcal{R}(\mathrm{\Gamma }(l),{\mathbf{F}}_{:,:,l}\circ \mathbf{R}),$$

Eq. (27)

$${\mathcal{G}}_{\mathrm{tu}}(\mathbf{F})=\sum _{l=0}^{L-1}{\mathbf{F}}_{:,:,l}.$$

Here, ${\mathcal{G}}_{\mathrm{ts}}(·):{\mathbb{R}}^{M\times N\times L}\to {\mathbb{R}}^{M\times [N+(L-1)]}$ and ${\mathcal{G}}_{\mathrm{tu}}(·):{\mathbb{R}}^{M\times N\times L}\to {\mathbb{R}}^{M\times N}$. They are shown as the magenta and orange layers in Fig. 4. $\circ$ represents the Hadamard product. The operator $\mathcal{R}(·):{\mathbb{R}}^{M\times N}\to {\mathbb{R}}^{M\times [N+(L-1)]}$ introduces a right-zero-padding (i.e., $[{\mathbf{F}}_{:,:,l},\mathbf{0}]$ with $\mathbf{0}\in {\mathbb{R}}^{M\times (L-1)}$) followed by a right-horizontal circular shifting of $\mathrm{\Gamma }(l)$ pixels. A script of Python to construct the sensing operators ${\mathcal{G}}_{\mathrm{ts}}$ and ${\mathcal{G}}_{\mathrm{tu}}$ are presented in Algorithm 4.

Algorithm 4

Programming the direct sensing operators of the time-sheared view and the time-unsheared view (i.e., Gts and Gtu) using TensorFlow.

## Direct sensing operator time-sheared view

class DirectSensing_ts(tf.keras.layers.Layer):

def __init__(self, L, M, N, **kwargs):

self.L = L

self.M = M

self.N = N

super(DirectSensing_ts, self).__init__(**kwargs)

def get_config(self):

config = super().get_config().copy()

config.update({

‘bands’: self.L})

return config

def call(self, F, R, **kwargs):

F = tf.multiply(R, F)

F = tf.pad(F, [[0, 0], [0, 0], [0, self.L - 1], [0, 0]],name = “padsensing”)

Gts = None

for i in range(0,self.L):

if Gts is not None:

Gts = Gts  + tf.roll(F[:,:,:,i], shift=i, axis=2)

else:

Gts = F[:,:,:,i]

Gts = tf.expand_dims(Gts ,axis=-1)

Gts = tf.math.divide(Gts ,self.L)

return Gts

## → Output: Compressed measurement

## Direct sensing operator time-unsheared view

Gtu = tf.math.reduce_mean(F,axis=-1)

Then, the transpose sensing operators of the time-sheared view and the time-unsheared view, shown as the red and purple layers respectively in Fig. 4, are defined as

Eq. (28)

$${\mathcal{F}}_{\mathrm{ts}}({\mathbf{G}}_{\mathrm{ts}},\mathbf{R})=\mathcal{S}(\mathrm{\Gamma }(l),{\mathbf{G}}_{\mathrm{ts}}\circ \mathcal{R}(\mathrm{\Gamma }(l),\mathbf{R})),$$

Eq. (29)

$${\mathcal{F}}_{\mathrm{tu}}({\mathbf{G}}_{\mathrm{tu}})={({\mathbf{G}}_{\mathrm{tu}})}_{:,:,l}.$$

Here, ${\mathcal{F}}_{\mathrm{ts}}(·):{\mathbb{R}}^{M\times [N+(L-1)]}\to {\mathbb{R}}^{M\times N\times L}$ and ${\mathcal{F}}_{\mathrm{tu}}(·):{\mathbb{R}}^{M\times N}\to {\mathbb{R}}^{M\times N\times L}$. They return a datacube from a 2D compressed measurement. $\mathcal{S}(·)$ is an operator that performs a left-horizontal circular shifting of $\mathrm{\Gamma }(l)$ pixels, followed by the removal of the last $(L-1)$ columns in the resulting shifted matrix to preserve the spatial dimension of the datacube. Algorithm 5 presents a Python script to construct ${\mathcal{F}}_{\mathrm{ts}}$ and ${\mathcal{F}}_{\mathrm{tu}}$.

Algorithm 5

Programming the transpose sensing operators of the time-sheared view and the time-unsheared view (i.e., Fts and Ftu) using TensorFlow.

## Transpose sensing operator of the time-sheared view

class Transposesensing_ts(tf.keras.layers.Layer):

def __init__(self, L, M, N, **kwargs):

self.L = L

self.M = M

self.N = N

super(Transposesensing_ts, self).__init__(**kwargs)

def get_config(self):

config = super().get_config().copy()

config.update({

‘bands’: self.L})

return config

def call(self, Gts, R, **kwargs):

F = None

R = R[0,:,:,0]

Gts = Gts[:,:,:,0]

for i in range(0,self.L):

if F is not None:

Ab = tf.roll(Gts, shift=-i, axis=2)

Ax = tf.expand_dims(tf.multiply(R, Ab[:,:,0:self.N]), -1)

F = tf.concat([F, Ax], axis=-1)

else:

Ab = tf.roll(Gts, shift=0, axis=2)

F = tf.expand_dims(tf.multiply(R,Ab[:,:,0:self.N]), -1)

F = self.L*F

return F

## Transpose sensing operator of the time-unsheared view

Gtu = tf.expand_dims(Gtu,axis=-1)

F_tu = tf.broadcast_to(Gtu, [Gtu.shape[0], M, N, L])

Finally, the inverse operator of the time-sheared view, shown as the brown layer in Fig. 4, is defined as

Eq. (30)

$${\mathcal{I}}_{\mathrm{ts}}({\mathbf{G}}_{\mathrm{ts}},\mathbf{R})={\mathbf{G}}_{\mathrm{ts}}{\circ (\sum _{l=0}^{L-1}\mathcal{R}(\mathrm{\Gamma }(l),{\mathbf{R}}^{°2})+\tilde{\rho }\mathbf{I})}^{°-1},$$

where ${\mathcal{I}}_{\mathrm{ts}}(·):{\mathbb{R}}^{M\times [N+(L-1)]}\to {\mathbb{R}}^{M\times [N+(L-1)]}$. $(·{)}^{°2}$ and $(·{)}^{°-1}$ represent the Hadamard quadratic power and the Hadamard inverse operation, respectively. An example script of Python to construct the inverse operator of the time-sheared view ${\mathcal{I}}_{\mathrm{ts}}$ is summarized in Algorithm 6.

Algorithm 6

Programming the inverse operator of the time-sheared view (i.e., Its) using TensorFlow.

class InverseOperator_ts(tf.keras.layers.Layer):

def __init__(self, L, M, N, **kwargs):

self.L = L

self.M = M

self.N = N

super(InverseOperator_ts, self).__init__(**kwargs)

def get_config(self):

config = super().get_config().copy()

config.update({

‘bands’: self.L})

return config

def build(self, input_shape):

Lambda = tf.constant_initializer(1)

Tau = tf.constant_initializer(1)

Psi = np.zeros([self.M,self.N+self.L-1])

Psi = tf.constant_initializer(Psi)

self.Lambda = self.add_weight(name=“Lbd,” initializer=Lambda, shape=(1),trainable=True)

self.Tau = self.add_weight(name=“Tau,” initializer=Tau, shape=(1),trainable=True,constraint=tf.keras.constraints. MaxNorm(max_value=1, axis=0))

self.Psi = self.add_weight(name=“Psi,” initializer=Psi, shape=(self.M,self.N+self.L-1),trainable=True)

super(InverseOperator_ts, self).build(input_shape)

def call(self, Gts, R, **kwargs):

Gts = Gts[:,:,:,0]

R1 = tf.broadcast_to(R,[1,self.M,self.N,self.L])

Gp = DirectSensing_ts(L=self.L, M=self.M, N=self.N, name=‘DirectPr_InitInv’)(R, R1)

Gp = Gp[:,:,:,0]

Gp = Gp/(self.Lambda) + tf.ones(Gp.shape)

Inv = tf.math.reciprocal(Gp, name=None)

Gts = tf.multiply((self.Tau**2)*Inv+(1-self.Tau**2)*self.Psi,Gts)

Gts = tf.expand_dims(Gts,axis=-1)

F = TransposeSensing_ts(L=self.L, M=self.M,N=self.N, name=‘TransPr_InitInv’)(Gts,R)

#

F = F/(self.Lambda**2)

return F

Following the definition of these five operators, the next step is to model the SWM matrix approach. Toward this goal, Eq. (25) is split into two main equations ${\mathbf{\Phi }}_{\mathrm{ts}}^T{\mathbf{g}}_{\mathrm{ts}}+\tilde{\rho }(\mathbf{z}-\mathbf{w})$ and ${\tilde{\rho }}^{-1}\mathbf{I}-{\tilde{\rho }}^{-1}{\mathbf{\Phi }}_{\mathrm{ts}}^T{[\tilde{\rho }\mathbf{I}+{\mathbf{\Phi }}_{\mathrm{ts}}{\mathbf{\Phi }}_{\mathrm{ts}}^T]}^{-1}{\mathbf{\Phi }}_{\mathrm{ts}}$. In the D-HAN, the first equation is reflected as ${\mathcal{F}}_{\mathrm{ts}}$ coupled to two 2D convolutional layers, each of which with a batch normalization operation (referred to hereafter as a 2D convolutional + batch normalization (BN) layer and shown in cyan in Fig. 4). The output from the 2D convolutional + BN layer is added with an estimate from the time-unsheared view generated by ${\mathcal{G}}_{\mathrm{tu}}$ and ${\mathcal{F}}_{\mathrm{tu}}$. Subsequently, the second equation is represented by two parallel arms. The upper arm, corresponding to ${\tilde{\rho }}^{-1}{\mathbf{\Phi }}_{\mathrm{ts}}^T{[\tilde{\rho }\mathbf{I}+{\mathbf{\Phi }}_{\mathrm{ts}}{\mathbf{\Phi }}_{\mathrm{ts}}^T]}^{-1}{\mathbf{\Phi }}_{\mathrm{ts}}$, is composed of ${\mathcal{G}}_{\mathrm{tu}}$ as a first layer followed by ${\mathcal{I}}_{\mathrm{ts}}$ and ${\mathcal{F}}_{\mathrm{ts}}$ along with four 2D convolutional + BN layers. The bottom arm, which corresponds to ${\tilde{\rho }}^{-1}\mathbf{I}$, has three 2D convolutional + BN layers. The outputs of both arms are subtracted and given as the input to the U-Net in the D-HAN that reflects Eq. (24). In the U-Net, the datacube passes through an encoding pathway comprised of max-pooling layers that simultaneously reduce the spatial dimension and increase the channels. This down step returns a smaller-size datacube with the more meaningful high-level details of the image (e.g., edges, textures, or shapes) linked to the scene’s sparsity. Then, in the decoding step, comprised of upsampling layers, the U-Net reconstructs the full-size datacube (denoted by $\tilde{\mathbf{F}}$) using these learned high-level details.

The loss function $\mathcal{L}(·)$, used to learn the D-HAN’s weights, is established as

Eq. (31)

$$\mathcal{L}(\mathbf{F})=l_1(\mathbf{F},\tilde{\mathbf{F}})+l_1({\mathbf{G}}_{\mathrm{ts}},{\tilde{\mathbf{G}}}_{\mathrm{ts}})+l_1({\mathbf{G}}_{\mathrm{tu}},{\tilde{\mathbf{G}}}_{\mathrm{tu}})+l_{\mathrm{SSIM}}(\mathbf{F},\tilde{\mathbf{F}}),$$

where $\tilde{\mathbf{F}}$ is the D-HAN’s output, ${\tilde{\mathbf{G}}}_{\mathrm{ts}}$ and ${\tilde{\mathbf{G}}}_{\mathrm{tu}}$ are estimations of the compressed measurement from $\tilde{\mathbf{F}}$ using Eqs. (26) and (27), respectively. $l_1(·)$ is the $l_1$-norm operator, and $l_{\mathrm{SSIM}}(·)$ represents the structural similarity (SSIM) index.¹⁰⁷

2.2.3.

Simulation

CUP’s image reconstruction of the “cell-division” scene is simulated using both the analytical-modeling-based algorithm (in Matlab) and deep-learning algorithm (in Python). The dimensions of the datacube were set as $M\times N\times L=256\times 256\times 25\text{pixels}$, and the binary mask holds the structure proposed in Ref. 104. Four popular databases—“SumMe,”¹⁰⁸ “Need for Speed,”¹⁰⁹ “Sports Videos in the Wild,”¹¹⁰ and “Mouse Embryo Tracking”⁶⁷—were used to train the D-HAN. The PnP-ADMM algorithm and a pretrained version of the D-HAN can be downloaded from Ref. 92 (Matlab 2022b) and Ref. 111 (Python, TensorFlow). In addition, a more beginner-friendly Python version is available in Ref. 112, which was trained on the Google Colaboratory (CoLab) application—a free Jupyter Notebook interactive development environment for Python hosted in Google’s cloud.

Six exemplary frames of the scene (as the ground truth) and their corresponding frames reconstructed by single-view and dual-view CUP using the PnP-ADMM and the D-HAN are shown in Fig. 5(a). The movie is shown in Video 1. As shown in Figs. 5(b) and 5(c), results show that implementing the dual-view approach exceeds the reconstruction performance of a single-view CUP in terms of the average peak signal-to-noise ratio (PSNR) defined as $\overline{\mathrm{PSNR}}=\frac{1}{L}\sum _{l=0}^{L-1}\left[10{\log }_{10}\right(\frac{[\max ({\mathbf{F}}_{:,:,l}){]}^2}{m_{\mathrm{tu}}^{-1}{\Vert \mathrm{vec}({\mathbf{F}}_{:,:,l})-\mathrm{vec}({\hat{\mathbf{F}}}_{:,:,l})\Vert }_2^2}\left)\right]$ and the average SSIM index¹¹³ defined as $\overline{\mathrm{SSIM}}=\frac{1}{L}\sum _{l=0}^{L-1}[[\mathrm{Lum}({\mathbf{F}}_{:,:,l},{\hat{\mathbf{F}}}_{:,:,l}){]}^{\alpha }[\mathrm{Cont}({\mathbf{F}}_{:,:,l},{\hat{\mathbf{F}}}_{:,:,l}){]}^{\beta }{[\mathrm{Struc}({\mathbf{F}}_{:,:,l},{\hat{\mathbf{F}}}_{:,:,l})]}^{\gamma }]$. Here, $\mathrm{vec}(·)$ is a vectorization operator, and $\hat{\mathbf{F}}$ is the reconstructed result. The operators $\mathrm{Lum}(·)=\frac{2{\mu }_x·{\mu }_y+C_1}{{\mu }_x^2+{\mu }_y^2+C_1}$, $\mathrm{Cont}(·)=\frac{2{\sigma }_x·{\sigma }_y+C_2}{{\sigma }_x^2+{\sigma }_y^2+C_2}$, and $\mathrm{Struc}(·)=\frac{{\sigma }_{xy}+C_3}{{\sigma }_x{\sigma }_y+C_3}$ measure the similarities in luminance, contrast, and structure, respectively, where ${\mu }_x,{\mu }_y,{\sigma }_x,{\sigma }_y$, and ${\sigma }_{xy}$ are the local means, standard deviation, and cross-covariance for the images. $\{\alpha ,\beta ,\gamma \}>0$ are parameters used to adjust the relative importance of the three components. $C_1,C_2$, and $C_3$ are constants to stabilize the division with weak denominator. For the results shown in Fig. 5, SSIM’s parameters were set as $\alpha =\beta =\gamma =1$, $C_1={0.01}^2$, $C_2={0.03}^2$, and $C_3=C_2/2$. The D-HAN obtains a better average PSNR and a comparable average SSIM to the PnP-ADMM approach in both single-view and dual-view CUP.

Fig. 5

Simulation of CUP’s image reconstruction of the “cell-division” dynamic scene using the PnP-ADMM algorithm and the D-HAN. (a) Six selected frames of the ground truth and the single-view and dual-view reconstructions using the PnP-ADMM and D-HAN algorithms. (b) PSNR of each reconstructed frame. (c) As (b), but showing the SSIM index (Video 1, MP4, 408 KB [URL: https://doi.org/10.1117/1.JBO.29.S1.S11524.s1]).

CUP’s performance decreases with higher noise and stronger compression. Table 1 illustrates this general trend from an ablation analysis of the “cell-division” dynamic scene using the PnP-ADMM algorithm. Higher noise levels reduce spatial resolution. Higher compression ratios result in stronger blurring in the temporal shearing direction, which further decreases the spatial resolution in that direction.^59,91,114 Both factors hamper the reconstruction algorithm’s ability to accurately place the correct amount of intensity from the compressed snapshot to the appropriate spatiotemporal position in the reconstructed datacube.

Table 1

Average PSNRs in reconstructed datacubes with different compression ratios and signal-to-noise ratios (SNRs).

Compression ratio	SNR (dB)
	15	20	25	30	Infinity
10.1×	26.3 ± 0.8	27.9 ± 0.5	29.5 ± 0.7	30.0 ± 0.8	30.4 ± 0.9
11.9×	26.2 ± 0.6	27.7 ± 0.5	29.2 ± 0.6	29.9 ± 0.7	30.2 ± 0.8
16.4×	26.1 ± 0.5	27.5 ± 0.5	28.7 ± 0.5	29.7 ± 0.7	29.9 ± 0.7
22.9×	25.8 ± 0.8	27.4 ± 0.4	28.2 ± 0.6	28.8 ± 0.7	28.9 ± 0.8
45.6×	25.6 ± 0.9	27.1 ± 0.8	28.1 ± 0.9	28.6 ± 1.0	28.8 ± 1.1

3.1.

Spatial Encoder

The selection of a suitable spatial encoder in a CUP system includes encoding pattern design and the encoder’s implementation. Because CUP relies on CS principles, its sensing matrix can be designed based on the restricted isometry property (RIP) to ensure its incoherence to the representation matrix of the scene. Notably, the sensing matrix based on a random pattern has been verified to meet the RIP criterion for a wide range of representation bases.¹¹⁵ Therefore, pseudorandom masks [Fig. 6(a)] are dominantly implemented as spatial encoders in reported CUP systems.

Fig. 6

Representative encoding patterns for CUP. (a) Pseudorandom pattern. (b) Deep-learning-optimized pattern. Insets: Zoomed-in views of local regions.

The RIP also provides a valuable metric for evaluating the encoder’s quality. The general strategy is to reduce the coherence between the sensing matrix and the representation matrices to ensure that the projection of high-dimensional data onto a lower-dimensional space preserves the essential data features.^116,117 It guarantees that the compressed measurements retain sufficient information to accurately reconstruct the original signal. To date, several works have improved the mask via deep learning.^104,118,119 As an example, an encoding mask designed via the D-HAN is shown in Fig. 6(b), where the shearing operation in CUP’s forward model and the training data produce horizontal stripe-like structures.¹⁰⁴

As summarized in Table 2, four approaches have been used to implement CUP’s spatial encoders: digital micromirror devices (DMDs),^{17,66,102,120,121} liquid-crystal spatial light modulators (LC-SLMs),¹³⁴ high-definition printing,^97,135 and photolithography.⁹¹ Among them, the DMD, as a reflective binary-amplitude spatial light modulator,¹³⁶ can provide reconfigurable, stable, and broadband encoding [Fig. 7(a)]. However, due to the micromirror’s tilt angle, the DMD is often required to be placed in the Littrow configuration in CUP systems^17,37,59 to retro-reflect the incident light. Since the DMD is not parallel to the object plane, this design limits the field of view (FOV). Moreover, the DMD’s structure limits light efficiency in three main aspects. First, since the DMD has a $\sim 94\%$ fill factor, a part of incident light is lost in the gaps between neighboring micromirrors. Second, as a 2D diffraction grating,¹³⁶ the DMD has an overall diffraction efficiency of $\sim 86\%$,¹³⁷ which indicates energy loss in high-diffraction orders. Finally, the aluminum coating of the micromirrors has a reflectivity of 89% in the visible spectrum with a dip at around 800 nm corresponding to the absorption of inter-band transitions in aluminum. Consequently, the constructed CUP system may not have an optimal spectral response for the dynamic scenes under investigation.

Table 2

Representative approaches for CUP’s spatial encoders.

Fig. 7

Pseudorandom binary masks displayed on representative spatial encoders. (a) DMD. (b) LC-SLM. (c) Plastic mask fabricated by high-definition printing. (d) Chromium mask made by photolithography. Inset: Zoom-in view of a local region.

Another choice of reconfigurable spatial encoding is LC-SLMs. They have been widely implemented in coded optical imaging.^18,138–141 LC-SLMs can simultaneously modulate amplitude and phase in grayscale.¹³⁴ In the context of CUP, they can provide both reflective and transmissive spatial encoding [Fig. 7(b)]. Nonetheless, LC-SLMs could also bring in some limitations in spatial encoding. Its modulation is sensitive to both the wavelength and polarization. Moreover, a relatively low fill factor of transmissive LC-SLMs (e.g., 58%)¹⁴² and the flicker noise could limit pattern quality and encoding stability.¹⁴³

Besides using the programmable devices, an encoded mask can be directly fabricated on a substrate. As a representative approach, high-definition printing can manufacture encoding masks at up to 50,800 dots per in. resolutions, with up to $\sim 30\mathrm{in.}\times 30\mathrm{in.}$ in size, at $16.7 per ${\mathrm{in.}}^2$ (Ref. 144) [Fig. 7(c)]. In one printing task, users can pack multiple masks with different encoding pixel sizes down to $7\mu \mathrm{m}$ and different pattern types as well as calibration patterns, such as single pinholes, pinhole arrays, and slits. As another approach, photolithography can produce spatial encoders with nanometer-level encoding pixel sizes over inches [Fig. 7(d)]. As an example, a 3-in. × 3-in. mask with 125-nm resolution can be fabricated at ∼$6,000. As a well-established fabrication technique, photolithography can be used with various materials to target different spectral bands.¹⁴⁵ These fabricated coded masks can be directly inserted in CUP systems, which conserves space for a more compact system design. Although capable of providing broadband and transmissive encoding, these two approaches can only prepare fixed spatial encoders. In addition, the almost unavoidable defect pixels in the fabricated encoder request careful calibration to build an accurate sensing matrix.

3.2.

Temporal Shearing Unit

Depending on the necessity of external power, temporal shearing units can be classified into passive units and active units (Table 3). The former deflects the temporal information transferred to certain photon tags (e.g., wavelengths) by exploiting the properties in these tags (e.g., color dispersion). Being jitter-free, these compact units bring in stable operation without increasing the control complexity.⁵⁹ The active units are driven by time-varying electric signals to trigger deflection. Usually integrated into the detection side of the imaging systems, they enable receive-only detection, which is specifically suited for capturing self-luminescent and color-selective events.^16,18,66

Table 3

Representative methods for CUP’s temporal shearing units.

Figure 8 shows two examples of passive temporal shearing units. Both need to team up with a chirped ultrashort probe pulse, which maps the temporal information of the event to its spectral band. As shown in Fig. 8(a), the modulated chirped pulse is spatially dispersed by a grating.⁵⁹ In recent years, the development of metamaterials has made metalenses a potential passive temporal shearing unit. They consist of an array of waveguide structures with a subwavelength size, with resonant metamaterial elements etched into the surface [Fig. 8(b)].^152,153 Metalenses can strongly disperse light while manipulating its phase, amplitude, and polarization.¹⁵⁴ This property has been exploited in hyperspectral imaging.¹⁵⁵ Grafting this sensing paradigm in CUP, a metalens integrates imaging and temporal shearing, which greatly reduces the system’s size and complexity.¹⁴⁶ Besides the aforementioned two units, other dispersive optical elements such as kinoforms,¹⁵⁶ zone plates,¹⁵⁷ and diffractive optical elements (DOEs)¹⁵⁸ could also be used for passive temporal shearing of chirped pulses.

Fig. 8

Representative passive temporal shearing units for CUP. The temporal information is mapped to the spectrum and deflected to different spatial positions by (a) a diffraction grating and (b) a metalens. $t_1$ to $t_n$, temporal information.

Active temporal shearing units have also been featured in many CUP systems. As an example, the image-converter streak tube is shown in Fig. 9(a). Such a device works by directing the dynamic scene onto a photocathode, where the incident photons are converted to photoelectrons. After being accelerated by a pulling voltage added on a metal mesh, these photoelectrons are temporally sheared by a varying electric field produced by applying a voltage to a pair of sweep electrodes. Then, the photoelectric signal is amplified by a microchannel plate. Finally, the photoelectrons bombard a phosphor screen and are converted back to photons.^59,91 The configuration of the image-converter streak tube takes advantage of the movement of electrons under high-voltage electric fields, enabling ultrafast shearing for the CUP system to provide up to femtosecond-level temporal resolution.^58,59 However, this operation is inevitably affected by the electronic jitter. Moreover, due to the space-charge effect in electronic imaging,¹⁶⁰ a trade-off needs to be made between the incident light intensity and the signal gain, which limits the imaging quality of the streak tube-based CUP systems.^37,151 The efficiency of image-converter streak tubes is also inherently limited by the photon–electron–photon conversion. The quantum yield of the photocathode is moderate for the visible light and decreases dramatically for near-infrared light.^161,162 The phosphor screen also has a relatively low conversion efficiency, especially for the fast-responding types.¹⁶³ The limited overall efficiency makes the image-converter streak tubes less suitable for imaging faint transient events.

Fig. 9

Representative active temporal shearing units for CUP. (a) Image-converter streak tube. (b) Rotating mirror. (c) TDI mode of a CCD camera. (d) Electro-optical deflector. (e) Ultrashort-pulse-induced ${\mathrm{CO}}_2$ molecule deflector. $\alpha$, Deflection angle. (b) Reprinted with permission from Ref. 135. (d) Adapted with permission from Ref. 159. (e) Adapted with permission from Ref. 151.

Rotating mirrors are another popular choice of active temporal shearing units for CUP. The mirror rotation continuously alters the angle of incidence, hence shearing the reflected light [Fig. 9(b)]. Rotating mirrors are preferred to be placed at the Fourier plane of a $4f$-system so that after the second lens, the chief rays of all temporal frames can propagate and enter the sensor perpendicularly, which avoids aberrations introduced by the field curvature.¹⁶⁴ Producing tunable temporal resolutions typically from hundreds of nanoseconds to microseconds, they are much slower than the image-converter streak tube. However, the all-optical operation avoids the space-charge effect, which enables optics-limited spatial resolution and high dynamic ranges.³⁷ Moreover, by circumventing the photon-to-photoelectron conversion in a photocathode, rotating-mirror-based CUP systems can employ sensors in matching responsive bands to sense photons with relatively low energy (e.g., in the infrared range). Leveraging high reflectivity coatings (e.g., $>95\%$ at 0.4 to $20.0\mu \mathrm{m}$),¹⁶⁵ these CUP systems are attractive candidates for high-sensitivity transient imaging at broad spectral bands.

Besides these two popular approaches, other specialized optical and/or electronic devices have been implemented as CUP’s temporal shearing units. As an example, Fig. 9(c) shows the operation principle of the time-delay-integration (TDI) mode of a CCD camera. Initially developed to visualize moving objects under extremely low light levels, the TDI configuration employs a long exposure during which the generated photoelectrons shift down row by row before eventually reading out.¹⁶⁶ In this way, the read-out data are the integration of information from different rows at different time points. Such a mechanism enables TDI cameras to combine the operations of temporal shearing and spatiotemporal integration, which considerably reduces the system’s complexity.^134,150 Recently, electron-transfer-based temporal shearing has also been implemented in a streak-camera sensor.^167,168 Hundreds of sampling and storage cells are placed underneath a line of photodiodes. During the sensor’s exposure, the 1D signal is sampled and sequentially stored at a temporal resolution of 500 ps. Although its 1D FOV excludes its implementation with CUP, this highly integrated device marks its potential to be further developed for future CUP systems.

Electro-optic crystals can also be used as the temporal shearing unit of CUP systems. As shown in Fig. 9(d), a time-varying electric field is applied to modulate the gradient of the refractive index of an electro-optic crystal. In this way, this electro-optic deflector (EOD) can direct the incident light to different propagation directions according to its time of arrival.^159,169 The EOD is currently the only all-optical shearing unit capable of achieving $50\times {10}^9$ frames per second (fps) in a CUP system.¹²² However, the shortcomings of small numerical aperture, high operating voltage, and limited deflection angle still hinder EODs for further applications in CUP.

Finally, transient materials’ behaviors have been proposed as CUP’s temporal deflectors. Figure 9(e) depicts how the transient alignment of ${\mathrm{CO}}_2$ molecules excited by an ultrashort laser pulse can induce a time-varying refractive-index gradient, resulting in different deflection angles to temporally shear the dynamic scene.¹⁵¹ Although having not been experimentally demonstrated, this mechanism could open a new avenue of transient-event-assisted ultrafast imaging. The fast responses of properly selected materials could push CUP’s imaging speed to the quadrillion fps level.¹⁷⁰

3.3.

2D Sensor

After being spatially encoded and temporally sheared, the dynamic scene is spatiotemporally integrated over each pixel by a 2D sensor. Most of the current commercial cameras (e.g., CCD, CMOS, scientific CMOS, and electron-multiplying CCD cameras) have been implemented to construct a CUP system. Nonetheless, used as the last component in a CUP system, the 2D sensor needs to be carefully selected to accommodate the characteristics of the dynamic scenes, the spatial encoders, and the temporal shearing units. In terms of spectral responsiveness, the 2D sensors are desired to have the highest sensitivity at the corresponding spectra of the dynamic scene. However, it might be restricted by the device. For example, for the image-converter streak tube, the quantum yield of the deployed camera should be peaked at the wavelength of the phosphor screen (e.g., 540 nm of a P43 phosphor screen). The pixel size of these sensors is required to sufficiently sample each encoding pixel for the given system’s magnification.

The shutter type is another important factor in the sensor selection. Overall, the global shutter is much preferred for CUP operation compared with the rolling shutter. Figure 10(a) shows a simulated dynamic scene of a rotating spinner with constant intensity. For a rolling-shutter sensor, the exposure of each row starts sequentially from the top to the bottom for the same period and ends at a different time point. The induced rolling-shutter effect distorts the image of fast-moving objects.¹⁷¹ CUP can overcome this distortion by putting the information back in the correct spatiotemporal position. However, due to the different starting times of exposure, only a part of FOV can be reconstructed for images at the beginning and the end of the movie, as shown in Fig. 10(b). This issue can be bypassed by limiting the occurrence of the dynamic event when all rows are under exposure. In contrast, the global shutter, which can be implemented in both CCD and CMOS sensors, allows capturing the dynamic scene over the full FOV [Fig. 10(c)] and thus avoid the time-windowing effect of the rolling shutter.

Fig. 10

Comparison between the rolling shutter and the global shutter in CUP’s data acquisition. (a) 10 representative frames of a simulated rotating-spinner scene. (b) Rolling shutter’s operating principle (top-left panel), the produced 2D snapshots (top-middle and top-right panels), and the illustration of CUP’s reconstructed frames (bottom panel). (c) As (b), but showing the results produced by the global shutter. F1 to F10, frame indices.

3.4.

Calibrations

In this section, we outline a few important calibration steps in CUP’s operation. They are necessary for both physical data acquisition and computational reconstruction of the dynamic scene.

3.4.1.

Co-registration of multiple views

Due to the difference among individual imaging arms, the acquired snapshots may have different aberrations. It is thus indispensable to accurately co-register all the views for accurate image reconstruction. Toward this goal, a static image of the time-sheared view is acquired by turning off the shearing unit (Fig. 11). In Matlab, the co-registration can be carried out using “control point registration” in the “Image Processing Toolbox.”¹⁷² The function “cpselect” opens a window for the user to select at least four pairs of control points in both views. Then, the function “fitgeotform2d” estimates the transformation matrix that best aligns the control points. Finally, the function “imwarp” applies the transformation matrix to complete the co-registration. The co-registered time-unsheared view and the time-sheared view are then fed into the reconstruction algorithms.

Fig. 11

Co-registration for dual-view CUP.

3.4.2.

Acquisition of the encoding mask

The experimentally captured encoding mask image produces better reconstruction than the design file of the used pattern because it takes into consideration various practical imperfections. For example, the DMD’s micromirrors may have a different orientation than those of the sensors. The fixed encoders fabricated by high-definition printing or photolithography also have defective pixels or membrane curving. These cases cannot be eliminated even if the imaging system is tuned with a proper magnification that matches the size of the encoding pixels to that of the sensor’s pixels. Thus, a mask image is captured by tuning off the shearing unit and then binarized for CUP’s image reconstruction (Fig. 12). Besides background subtraction and white-field correction, threshold selection and edge detection are combined to optimize binarization.⁵⁷ This calibration can also reduce aberrations and field curvature.¹²³ In practice, the FOV and the maximum shearing distance are also limited to ensure high quality in captured images.

Fig. 12

Binarization of the captured encoding mask image. (a) Section cropped from the acquired mask image. (b) Cropped section after background subtraction and white-field correction. (c) Image binarization by applying a threshold to (b). (d) Image binarization by detecting edges in panel (b). (e) Combining (c) and (d) using OR operation. Adapted with permission from Ref. 57.

3.4.3.

Linearity test of shearing operation

Linear temporal shearing is used in CUP’s forward model (see Sec. 2.1). However, various experimental factors could deviate the linear temporal shearing operation, including misalignment, jitter, and imperfect instrument responses.¹⁷³ Therefore, a linearity test is required to assess the system’s performance and to compensate for these factors. An example of a rotating-mirror-based CUP system is shown in Fig. 13(a). About 100 frames containing number indices and short lines were displayed on a DMD at 20 kHz. From the recorded snapshot, the displacements between the centroids of the consecutive short lines were calculated to determine the rotating mirror’s shearing operation. In this example, the shearing deviates from a linear function by 2 pixels over 100 frames.¹⁰⁴

Fig. 13

Linearity test of the temporal shearing operation. (a) Test of a rotating-mirror-based temporal shearing. Left panel: composite of the 100 frames with consecutive short lines and frame indices. Right panel: analysis of displacement in the time-sheared view. (b) Test of a diffraction-grating-based temporal shearing. Left panel: setup (top) of generating pulses with selected wavelengths (bottom). Right panel: result of the linearity test by illuminating a square pattern with the narrow-band pulses. (c) Test of a streak-tube-based temporal shearing. Left panel: setup of the test. Right panel: cross-section in the streak measurement of the pulse train generated by the etalon. (a) Adapted with permission from Ref. 104. (b) Adapted with permission from Ref. 59. (c) Adapted with permission from Ref. 174.

As another example, Fig. 13(b) shows the linearity test of a diffraction-grating-based CUP system.⁵⁹ A tunable bandpass filter was built based on a rotating grating [top-left panel in Fig. 13(b)] to produce pulses with a selected wavelength [bottom-left panel in Fig. 13(b)]. The generated narrowband pulses illuminated a small square pattern, whose positions in the streak images were measured to obtain the relationship with wavelengths [right panel in Fig. 13(b)]. Finally, an example of the linearity test of an image-converter streak tube is shown in Fig. 13(c). Following a calibration protocol similar to that of the diffraction grating-based CUP system, a pulse train with a known interval was generated by an etalon. The linearity was computed by measuring the deflected pulses’ positions.¹⁷⁴

4.1.

Neuroimaging

Monitoring the spatiotemporal dynamics of neuron signaling is essential to the understanding of the brain’s structure and function. Direct visualization can aid researchers and clinicians in studying neurological disorders, cognitive processes, and brain development. Frame rates at the level of one billion fps are demanded to image the propagation of action potentials (APs) in myelinated axons ($\sim 100\mathrm{m}/\mathrm{s}$) with high spatial resolution and in real time. Unreachable by conventional electronic sensors, this requirement poses a considerable technical challenge to neuroimaging research.

Overcoming this challenge, CUP imaged phase and lifetime dynamics evoked by neuronal activities. As an example, by combining Mach–Zehnder interferometry and utilizing its ultrafast imaging speed and large sequence depth, differentially enhanced CUP (Diff-CUP) imaged internodal current flow in myelinated axons from the sciatic nerves of Xenopus laevis frog at $200\times {10}^9\mathrm{fps}$⁶⁵ [Fig. 14(a)]. The high phase sensitivity of Diff-CUP enables simultaneously capturing the substantial cellular deformations and consequent phase alterations induced by passive current flows (i.e., without the amplification of the electrical current)^175,176 resulting from a 10-V and $1\text{-}\mu \mathrm{s}$ pulse injected into the axon [Fig. 14(b)]. The reconstructed correlation curves of each segment of the FOV (labeled with numbers 1 to 8) reveal the microsecond-level phase changes induced by the propagating internodal current flow [Fig. 14(c)], whose conduction speeds in myelinated axons were calculated to be $100\pm 26\mathrm{m}/\mathrm{s}$. To date, Diff-CUP is the fastest imaging-based approach for assessing AP-related conduction.

Fig. 14

CUP of neuronal activities. (a) Schematic of Diff-CUP. Inset: Adhesion microscope slide for Diff-CUP imaging. BS, beam splitter; DG, delay generator; $E(t)$, transient field stimulation; HWP, half-wave plate; LN, lithium niobate; OB, objective lens; PBS, polarizing beam splitter; PG, pulse generator; SC, streak camera. (b) Spatiotemporal interferogram of a propagating internodal current flow in a myelinated axon captured by Diff-CUP. (c) Reconstruction of the current flow signals based on the stimulus interferogram. Black dashed lines indicate the signal region of the internodal current flow. T, the propagation time of the internodal current flow within the FOV. (d) Schematic of compressed FLIM. (e) Six representative frames from the reconstruction of a cultured hippocampal neuron upon potassium stimulation at 100 fps. (f) Time-lapsed lifetime and intensity curves of a cultured hippocampal neuron. (g) Intensity (top panel) and lifetime (bottom panel) waveforms of neural spikes (black lines) and their means (green lines) for a cultured hippocampal neuron under stimulation. (a)–(c) Adapted with permission from Ref. 65. (d)–(g) Adapted with permission from Ref. 60.

CUP is also implemented as a CS-based fluorescence lifetime imaging microscopy (FLIM)⁶⁰ to record high-resolution 2D lifetime images of immunofluorescently stained neurons [Fig. 14(d)]. With an imaging speed of $\sim 10\times {10}^9\mathrm{fps}$, this CUP-based FLIM system captured the fluorescence intensity decay in real time, which produced a 2D lifetime map. Leveraging the intrinsic frame rate of the internal CMOS camera, lifetime maps were generated at 100 fps. This technique visualized neural spike dynamics via the fluorescence intensity and donor lifetime decrease during Förster resonance energy transfer.¹⁷⁷ Figure 14(e) illustrates six representative lifetime images of a cultured hippocampal neuron at 100 fps. The time courses of the averaged fluorescence intensity variation and lifetime of this sample over 1 s are plotted in Fig. 14(f). Finally, the hippocampal neuron’s fluorescence intensity and lifetime waveforms of single AP and their means [black lines and green line in Fig. 14(g), respectively] were acquired experimentally. This analysis revealed that a single spiking event led to an average relative fluorescence intensity change ($\mathrm{\Delta }F/F$) of $-2.9\%$ and a lifetime change of $-0.7\mathrm{ns}$.

4.2.

Temperature Sensing

Temperature, as an important biomarker, is linked to many biological processes (e.g., metabolism¹⁷⁸) and medical procedures (e.g., photothermal therapy¹⁷⁹). Accurate and real-time temperature sensing is important to pathology diagnostics, physiology monitoring, and therapeutical efficiency. Photoluminescence thermometry presents an emerging method by utilizing the temperature-sensitive optical emissions of photoluminescent materials as well as optical detections at high spatial resolution. Its merits include noncontactness, high adaptability to a broad temperature range, high accuracy, flexibility in sample selection, and suitability for diverse environments.¹⁸⁰ Thus, photoluminescence thermometry is increasingly featured in recent advances in optical temperature measurements.

The success of photoluminescence thermometry depends on two essential constituents: temperature indicators and optical imaging instruments. Recent advances in biochemistry, materials science, and molecular biology have unveiled numerous labeling indicators for photoluminescence thermometry.^162,181,182 From semiconductor quantum dots¹⁸³ and organic fluorophores¹⁸⁴ to rare-earth-doped phosphors,¹⁸⁵ the diversity of these agents allows for tailored temperature sensing across different thermal sensitivities, optical properties, and response times for biomedical applications.^186–189 For example, lanthanide-doped upconverting nanoparticles (UCNPs), which can sequentially absorb two (or more) low-energy near-infrared photons and convert them to one higher-energy photon, enable biocompatible temperature sensing with low excitation power densities and high sensitivity.^190,191

CUP has enabled wide-field temperature mapping using photoluminescence lifetimes of UCNPs.¹³⁵ In the schematic shown in Fig. 15(a), near-infrared pulses, generated by a 980-nm continuous-wave (CW) laser and an optical chopper, are focused on the back focal plane of an objective lens to form wide-field illumination. The excited UCNPs on the sample emit visible upconversion luminescence. After passing the filter, the dynamic photoluminescence of a selected emission band is imaged by a rotating mirror-based dual-view CUP system at 33,000 fps. The reconstructed lifetime images in the UCNPs’ two upconversion emission bands at different temperatures are shown in Fig. 15(b). The averaged intensity decays [Fig. 15(c)] enable the establishment of the temperature-lifetime relationship [Fig. 15(d)]. Furthermore, the system tracked the 2D temperature of a moving onion epidermis sample labeled by UCNPs at a rate of 20 lifetime-maps per second [Fig. 15(e)]. The intensity decays of four selected areas [labeled in the top-left panel in Fig. 15(e)] are shown in Fig. 15(f). It is worth noting that the fluences of the four selected areas are different but the measured photoluminescence lifetimes remain stable, showing the lifetime-based approach contributed by CUP is more reliable in accurate temperature sensing.

Fig. 15

CUP of temperature sensing. (a) Schematic of wide-field photoluminescence lifetime thermometry based on a dual-view rotating-mirror CUP system. L1 to L5, lenses. (b) Lifetime maps of the two emission bands (i.e., ${{}^4\mathrm{S}}_{3/2}\to {{}^4\mathrm{I}}_{15/2}$ and ${{}^4\mathrm{F}}_{9/2}\to {{}^4\mathrm{I}}_{15/2}$) of the used UCNPs under different temperatures. (c) Normalized photoluminescence decays of the two emission bands after averaging over the FOV. (d) Temperature-lifetime relationship of both emission bands. (e) Selected time-unsheared views (top row) and reconstructed lifetime maps (bottom row) of a moving onion epidermis cell sample labeled by UCNPs. (f) Intensity decays at four selected areas with different intensities marked in (e). Adapted with permission from Ref. 135.

4.3.

Microfluidics

A rotating-mirror CUP system has been applied to the video recording of complex fluid dynamics and interactions at the microscale.¹⁴⁸ A schematic of rotating-mirror-based CUP is shown in Fig. 16(a). This system observed flow droplet samples within a microfluidic chip.¹⁴⁹ Two immiscible liquids of transparent oil and chemical dye, injected through a motorized dispenser, flowed in the chip channels at $0.9\mathrm{m}/\mathrm{s}$. Three separate measurements are recorded at 3000 fps, 50,000 fps, and 120,000 fps [Fig. 16(b)]. These experimental results show a high reconstruction quality with well-preserved edge features in the frames, showing clear and distinguishable droplets flowing in the microfluidic chip. These results show CUP’s potential to visualize cell-shape changes in response to rapid external stimuli or internal dynamics in microfluidics,^192–194 which will provide new insights into cellular biomechanical properties that are closely linked to cellular function and disease development.

Fig. 16

CUP of microfluidics. (a) Schematic of a rotating-mirror-based CUP system. (b) Snapshot of flowing immiscible liquids (i) and representative frames from the reconstructed videos at (ii) 3000, (iii) 50,000, and (iv) 120,000 fps. Adapted with permission from Ref. 148 and Ref. 149.

4.4.

Photoacoustic Imaging

CUP can also contribute to photoacoustic (PA) imaging. Figure 17 shows a simulation study on implementing CUP with optical interferometric detection of PA waves.^195,196 In the proposed system schematically shown in Fig. 17(a), a pulsed laser illuminates a biological sample. The induced PA effect generates thermoelastic initial pressure, which is detected at the surface of the sample with a Fabry–Pérot etalon (FPE). This interaction of the ultrasonic waves with the surface of the FPE results in the modulation of the reflected CW laser beam on the opposite side of the FPE.¹⁹⁷ The modulated CW laser beam is then imaged by a CUP system based on a DMD and a galvanometer scanner. Figure 17(b) shows a simulation of this method to image the initial pressure distribution of 12 vessel-like structures.

Fig. 17

CUP of photoacoustic imaging. (a) Proposed system schematic. CAM, camera; COL, collimator; CW, continuous-wave; DMD, digital micromirror device; FPE, Fabry–Pérot etalon; L, lens; $\lambda /4$, quarter wave plate; LP, linear polarizer; OI, optical isolator; PBS, polarizing beam splitter; SMF, single-mode fiber. (b) Simulation of image reconstruction of initial pressure distribution. Adapted with permission from Ref. 195.

5.1.

Faster

As currently the world’s fastest optical imaging technology, CUP naturally carries a mission to explore even higher speeds in optical imaging. Since CUP’s invention, innovation in temporal shearing units has been a focus of its technical improvement. From $100\times {10}^9\mathrm{fps}$ of the original CUP system,¹⁷ various image-converter streak tubes have been deployed to increase its imaging speed to $10\times {10}^{12}\mathrm{fps}$,¹¹⁴ which currently holds the world record for single-shot receive-only ultrafast imaging. However, at this speed, the image quality is considerably affected by the space-charge effect,¹⁶⁰ posing challenges in further improvement of frame rates. In the future, transient perturbation in reflective index induced by a temporally modulated ultrashort laser pulse or molecular orientation could bring higher imaging speeds and circumvent image degradation.^151,170

Leveraging the advance in chirped pulse illumination, CUP systems using passive temporal shearing units have boosted the imaging speed to $3.85\times {10}^{12}\mathrm{fps}$,¹²¹ $70\times {10}^{12}\mathrm{fps}$,⁵⁷ $219\times {10}^{12}\mathrm{fps}$,¹⁹⁸ and $256\times {10}^{12}\mathrm{fps}$.⁶³ The last value marks the fastest speed in single-shot optical imaging. In the future, by synergizing the ultra-broadband ultrashort pulses¹⁹⁹ and photonic streaking in gas,^170,200 CUP’s imaging speed could top quadrillion fps, entering the attosecond-level imaging regime.

5.2.

Clearer

A higher spatial resolution allows CUP to visualize fine details. In biomedicine, this ability can transfer to informative depiction of cellular and tissue morphology, accurate diagnostics, and precise treatment. Nonetheless, in CUP’s operation, both spatial encoding in data acquisition and denoising in image reconstruction could reduce the effective system bandwidth. To visualize the targeted spatial details, a common practice is to magnify the scene, which unavoidably reduces the FOV. Thus, how to regain the lost bandwidth to achieve a diffraction-limited spatial resolution is an important research direction of CUP. One potential approach is subpixel shifting.²⁰¹ In particular, a DOE could be used to duplicate the dynamic scene to multiple bands, each would be encoded with the same encoding mask but with a different subpixel shift. A joint image reconstruction using all the captured snapshots could recover the original optical bandwidth.

Another interesting research direction is super-resolution CUP. As many ultrafast phenomena also occur at the nanoscale, overcoming the diffraction limit in the CUP system will likely open avenues for new studies not possible before, including temperature dynamics in mitochondria,²⁰² conformational transitions of protein,²⁰³ evolutions of membrane fragments produced by cellular lysis.²⁰⁴ Toward this goal, CUP can be incorporated into existing super-resolution microscopy techniques (e.g., structured illumination microscopy²⁰⁵) or bypass the optical diffraction limit by electron imaging (e.g., transmission electron microscopy¹¹⁷).

5.3.

Broader Spectrum

Although CUP has been experimentally demonstrated in the ultraviolet, visible, and near-infrared spectral ranges, extending its imaging capability to a broader spectrum will likely continue in future application-driven development. Toward this goal, the spatial encoder should have high contrast in the desired spectrum. Besides the popular broadband metallic masks made from aluminum, silver, or chromium,^17,91,147 photonic crystals with broad tunable bandgaps can selectively block specific wavelengths,²⁰⁶ giving them the potential to be used as a spatial encoder in certain spectra. For temporal shearing units, the photocathode in the image-converter streak tube excludes high-sensitivity imaging for wavelengths of $>950\mathrm{nm}$. Contrarily, leveraging its all-optical functionality, rotating mirrors can fill out this gap, which will likely lead to the development of CUP for deep ultraviolet, mid-infrared, and far-infrared spectra. Moreover, advanced design and fabrication of metasurfaces and metalenses could potentially extend CUP to a spectrum from extreme ultraviolet to terahertz.^207,208

5.4.

Smarter

Many deep learning-based approaches have been used in CUP’s image reconstruction.^94–104 Harnessing the power of artificial intelligence, they have unlocked new capabilities for analyzing ultrafast events, such as real-time data processing, on-device analysis, and on-time feedback. It is expected that these deep-learning algorithms will provide multifacet supervision to CUP systems in the future. For example, the next-generation systems could autonomously adjust the patterns loaded on the spatial encoder according to the initial classification of the dynamic scene. These systems could also monitor the nonlinear shearing operation and adaptively compensate for it in image reconstruction or system alignment.¹⁰⁴

5.5.

Higher Dimensions

Recent developments in CUP have explored high-dimensional ultrafast imaging. To date, several advanced systems—such as multispectral CUP,¹²⁰ stereo-polarimetric CUP,⁵⁷ and spectral-volumetric CUP⁶²—have pushed the overall sensing capability to four dimensions and even five dimensions. In the future, by extending the configuration used in stereo-polarimetric CUP⁵⁷ to generate multiple perspectives of the dynamic scene, light-field imaging could be incorporated into CUP. Ultimately, single-shot imaging of seven-dimensional plenoptic function would be within reach. Using CUP to sense other photon tags that are not included in the conventionally defined plenoptic function is also a future direction. CUP has already enabled amplitude and phase imaging of a femtosecond laser pulse.⁶³ CUP could also be combined with other existing technologies, such as the transport-of-intensity equation ²⁰⁹ and coherent modulation,²¹⁰ for ultrafast quantitative phase imaging. Finally, recent advancements in on-chip polarization imaging and metasurface-based angular moment separation could incorporate these parameters into CUP’s measurement scope.^211,212

5.6.

Smaller

CUP systems with a compact size are important to studies that require restricted weight and space. For biomedicine, it will offer the ability to mount the system the same way as conventional cameras on microscopes and hand-held systems as well as in operating rooms. An innovative optical design that folds the optical path could reduce the system size.²¹³ Selecting a multifunctional component (e.g., a metalens and a TDI sensor) provides another approach to reducing the number of optical elements in CUP systems. Advances in sensor design and nanofabrication could provide the streak imaging sensor^167,168 with a 2D FOV. All efforts will contribute to engineering compact and even miniature CUP systems in the future.

5.7.

Cheaper

Besides reducing the size of a CUP system, making an economical CUP system carries considerable value from both research and commercialization perspectives. To manufacture a fixed spatial encoder, high-resolution printing offers the lowest cost (i.e., <$20 per ${\mathrm{in.}}^2$). For reconfigurable spatial encoders, a 0.47″ DMD chip ($1920\times 1080$ micromirrors; $5.4\text{-}\mu \mathrm{m}$ pitch) costs ∼$120.²¹⁴ Future development could develop a DMD controller specifically tailored for CUP with much-reduced functionality compared with existing ones to decrease the cost. For the temporal shearing unit, the first approach to reduce the cost is to find a replacement for expensive image-converter streak tubes. Electro-optic modulators have made their debut in this direction,¹²² which produced an imaging speed of $50\times {10}^9\mathrm{fps}$. For a rotating-mirror CUP system, a viable strategy is to add the spatial encoder and an affordable rotating mirror (e.g., galvanometer scanners and polygonal mirrors) in front of the CCD/CMOS cameras existing in the system. It is envisaged that a minor addition in cost could endow ultrahigh-speed imaging to existing cameras while retaining their inherent advantages (e.g., in sensitivity and sensing spectrum).