# On the impact of Performance Antipatterns in multi-objective software model
refactoring optimization ††thanks: This research was supported by the AIDOaRt
project (ECSEL-JU program - grant agreement n.101007350).

Vittorio Cortellessa1, Daniele Di Pompeo2, Vincenzo Stoico3, Michele Tucci4
Department of Information Engineering, Computer Science and Mathematics
University of L’Aquila
{1vittorio.cortellessa, 2daniele.dipompeo<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Software quality estimation is a challenging and time-consuming activity, and
models are crucial to face the complexity of such activity on modern software
applications.

One main challenge is that the improvement of distinctive quality attributes
may require contrasting refactoring actions on an application, as for trade-
off between performance and reliability. In such cases, multi-objective
optimization can provide the designer with a wider view on these trade-offs
and, consequently, can lead to identify suitable actions that take into
account independent or even competing objectives.

In this paper, we present an approach that exploits the $\mathtt{NSGA-II}$
multi-objective evolutionary algorithm to search optimal Pareto solution
frontiers for software refactoring while considering as objectives: i)
performance variation, ii) reliability, iii) amount of performance
antipatterns, and iv) architectural distance. The algorithm combines randomly
generated refactoring actions into solutions (i.e., sequences of actions) and
compares them according to the objectives.

We have applied our approach on a train ticket booking service case study, and
we have focused the analysis on the impact of performance antipatterns on the
quality of solutions. Indeed, we observe that the approach finds better
solutions when antipatterns enter the multi-objective optimization. In
particular, performance antipatterns objective leads to solutions improving
the performance by up to 15% with respect to the case where antipatterns are
not considered, without affecting the solution quality on other objectives.

###### Index Terms:

multi-objective optimization, software performance, software reliability,
antipatterns, refactoring.

##  I Introduction

In the last decades, multi-objective optimization techniques have been
successfully applied to many model-driven software development problems [1, 2,
3, 4, 5, 6, 7]. These techniques are evidently more effective on problems
whose objectives can be expressed through quantifiable metrics. Problems
related to non-functional aspects undoubtedly fit into this category, as
witnessed by the vast literature in this domain [8, 9, 10]. Most approaches
are based on evolutionary algorithms [11] that allow exploring the solution
space by (re-)combining solutions.

Software refactoring is a task that can be triggered by different causes, such
as the introduction of additional requirements, the adaptation to new
execution contexts, or the degradation of non-functional properties. The
identification of optimal refactoring actions is a non-trivial task, mostly
due to the large space of solutions, and there is still lack of automated
support in this context. A common aspect of multi-objective optimization
approaches applied to model-driven software refactoring problems is that they
search among design alternatives (_e.g.,_ through architectural tactics [12,
13]).

In this paper, we present an approach based on a multi-objective evolutionary
algorithm (_i.e.,_ $\mathtt{NSGA-II}$) that searches sequences of refactoring
actions to be applied on UML software models, leading to the optimization of
several objectives, namely: i) performance variation (analyzed through Layered
Queueing Networks), ii) reliability (analyzed through a closed-form model),
iii) number of performance antipatterns (automatically detected) and iv)
architectural distance. In particular, we are interested in studying the
impact of performance antipatterns on the quality of refactoring solutions.
Since it has been shown that removing performance antipatterns leads to
systems that show better performance than the ones affected by them [14, 15,
16], we aim at studying if this phenomenon also holds in the context of muti-
objective optimization, where performance improvement is not the only
objective.

Since UML represents a commonly recognized standard in the software modeling
domain, our approach applies to UML models augmented by MARTE [17] and DAM
[18] profiles, which allow to embed performance and reliability properties
into UML models. However, UML does not provide native support for performance
analysis, thus we introduce a novel model-to-model transformation that
generates Layered Queueing Networks (LQN) from annotated UML models. The
solution of LQN models feeds the performance variation objective.

Here we consider refactoring actions that are designed to improve performance
in most cases. Since such actions may also have an impact on other non-
functional properties, we introduce the reliability among the optimization
objectives to study whether their combination could keep satisfactory levels
of performance and reliability at the same time. To quantify the reliability
objective, we adopt an existing model for component-based software systems
[19] which can be generated from UML models.

In this context, a multi-objective optimization should also consider the
distance between the initial model and the models resulting from applying
refactoring actions. Indeed, without an objective that minimizes such
distance, the proposed solutions could be impractical because they would
require to completely disassemble and re-assemble a model architecture. This
is the rationale behind the inclusion of the architectural distance as an
objective to minimize in our optimization.

We recently investigated the effect of performance in a multi-objective
optimization problem by presenting _EASIER_ [20], where we also observed the
contribution of performance antipattern and architectural distance in the
generation of Pareto frontiers. However, _EASIER_ relies on a performance-
oriented architectural description language equipped with its own performance
solver, thus limiting its application to more complex cases.

Summarizing, the main contributions of this paper are the following:

  * •

we adopt UML as modeling notation, instead of DSLs that have been adopted in
existing multi-objective approaches [9, 21, 10, 22], thus we move toward a
widely adopted standard;

  * •

we introduce a novel model-to-model transformation to generate Layered
Queueing Networks from annotated UML models;

  * •

we introduce an automated refactoring engine equipped with four refactoring
actions that can be applied to UML models;

  * •

we adopt a probabilistic concept of antipattern occurrence [23], which allows
us to analyze the sensitivity of solutions to the antipattern detection
capability;

  * •

we consider measures of performance, reliability, distance, and the number of
performance antipatterns as competing objectives of our evolutionary process.

We extensively applied our approach to the Train Ticket Booking Service [24,
25] in order to answer the following three research questions:

  * •

_RQ1_ : Does antipattern detection contribute to find better solutions
compared to the case where antipatterns are not considered at all?

  * •

_RQ2_ : Does the probabilistic nature of fuzzy antipatterns detection help to
include higher quality solutions in Pareto frontiers with respect to
deterministic one?

  * •

_RQ3_ : Is the approach able to keep the system reliability satisfactorily
high?

Indeed, the experimentation lasted approximately 80 hours and generated more
than 31,000 solutions. Our results show that, by considering the reduction of
performance antipatterns as an objective, we are able to obtain refactoring
solutions that improve the performance by up to 15% with respect to the case
where antipatterns are not considered.

The structure of the paper is the following: Section II introduces basic
concepts, Section III describes the approach, in Section IV we evaluate our
approach on the case study and discuss the results and threats to validity.
Section V reports related work, and Section VI concludes the paper.

##  II Background

We identify the four competing objectives of our evolutionary approach as
follows: perfQ is a performance quality indicator that quantifies the
performance improvement/detriment between an initial model and one obtained by
applying the refactoring actions of a solution (Section II-A); reliability is
a measure of the reliability of the software system (Section II-B); #pas is a
metric that quantifies the amount of performance antipattern occurrences while
considering the intrinsic uncertainty rising from thresholds used by the
detection mechanism (Section II-C); archDist represents the distance between
an initial model and one obtained by applying the refactoring actions of a
solution (Section II-D).

We employ the Non-dominated Sorting Algorithm II ($\mathtt{NSGA-II}$) as our
genetic algorithm [26], since it is extensively used in the software
engineering community, _e.g.,_ [12, 27]. $\mathtt{NSGA-II}$ randomly creates
an initial population from that the offspring population is created by
applying the _Crossover_ operator with probability $P_{crossover}$, the
_Mutation_ operator with probability $P_{Mutation}$. The union of the initial
and the offspring populations is sorted by the _Non-dominated sorting_
operator, which identifies different Pareto frontiers with respect to
considered objectives. Finally, the _Crowding distance_ operator cuts the
worse half the sorted union off, and the remaining population becomes the
initial population for the next step.

###  II-A Performance Quality Indicator (perfQ)

Table I: Detectable performance antipatterns in our approach. Performance antipattern |  Description
---|---
Pipe and Filter |  Occurs when the slowest filter in a “pipe and filter” causes the system to have unacceptable throughput.
Blob |  Occurs when a single component either i) performs the greatest part of the work of a software system or ii) holds the greatest part of the data of the software system. Either manifestation results in excessive message traffic that may degrade performance.
Concurrent Processing System |  Occurs when processing cannot make use of available processors.
Extensive Processing |  Occurs when extensive processing in general impedes overall response time.
Empty Semi-Truck |  Occurs when an excessive number of requests is required to perform a task. It may be due to inefficient use of available bandwidth, an inefficient interface, or both.
Tower of Babel |  Occurs when processes use different data formats and they spend too much time in convert them to an internal format.

perfQ quantifies the performance improvement/detriment between two models, and
it is defined as follows:

$\textsc{perfQ}(M)=\frac{1}{c}\sum\limits_{j=1}^{c}p_{j}\cdot\frac{F_{j}-I_{j}}{F_{j}+I_{j}}$

where $M$ is a model obtained by applying a refactoring solution to the
initial model, $F_{j}$ is the value of a performance index in $M$, and $I_{j}$
is the value of the same index on the initial model. $p\in\\{-1,1\\}$ is a
multiplying factor that holds: i) $1$ if the $j$–th index has to be maximized
(i.e., the higher the value, the better the performance), like the throughput;
ii) $-1$ if the $j$–th index has to be minimized (i.e., the smaller the value,
the better the performance), like the response time.

Notice that, for performance measures representing utilization, $p$ also holds
$1$ but we deﬁne a _utilization correction factor_ $\Delta_{j}$ to be added to
each j–th term above, as defined in [20]. The utilization correction factor
penalizes refactoring actions that push the utilization too close to 1, i.e.,
its maximum value.

Finally, the global perfQ is computed as the average across the number $c$ of
performance indices considered in the performance analysis.

As mentioned in the introduction, in order to obtain performance indices of a
UML model, the analysis has to be conducted on a performance modeling
notation, and we have adopted the Layered Queueing Networks
(LQNs)111http://www.sce.carleton.ca/rads/lqns/LQNSUserMan-jan13.pdf for this
goal. LQNs have been introduced to give a layered structure to Queuing
Networks [28]. Each layer contains software resources providing different
classes of service. These elements are named _tasks_ , while classes of
service are identified by _entries_. Each task implements a queue that is
served by a _processor_. The behavior of tasks is modeled through requests
among entries, or it is specified using _activities_. An activity is the
atomic unit of computation in an LQN. Activities can be related using ordering
operators. An ordering operator indicates precedence between activities.
Figure 2 shows a graphical example of an LQN.

An LQN is solved by decomposing it in several conventional Queuing Network
submodels related to each other, one for each layer of the whole model. Each
submodel is in turn made of two layers: the top layer includes tasks
interpreting the clients, while tasks in the bottom layer are intended to be
servers of the top ones.

###  II-B Reliability model

The reliability model that we adopt here to quantify the reliability objective
is based on the model introduced in [19]. The mean failure probability
$\theta_{S}$ of a software system $S$ is defined by the following equation:

$\theta_{S}=1-\sum\limits_{j=1}^{K}p_{j}\left(\prod\limits_{i=1}^{N}(1-\theta_{i})^{InvNr_{ij}}\cdot\prod\limits_{l=1}^{L}(1-\psi_{l})^{MsgSize(l,j)}\right)$

This model takes into account failure probabilities of components
($\theta_{i}$) and communication links ($\psi_{l}$), as well as the
probability of a scenario to be executed ($p_{j}$). Such probabilities are
combined to obtain the overall reliability on demand of the system
($\theta_{S}$), which represents how often the system is not expected to fail
when its scenarios are invoked.

The model is considered to be composed of $N$ components and $L$ communication
links, whereas its behavior is made of $K$ scenarios. The probability
($p_{j}$) of a scenario $j$ to be executed is multiplied by an expression that
describes the probability that no component or link fails during the execution
of the scenario. This expression is composed of two terms:
$\prod_{i=1}^{N}(1-\theta_{i})^{InvNr_{ij}}$, which is the probability of the
involved components not to fail raised to the power of their number of
invocations in the scenario (denoted by $InvNr_{ij}$), and
$\prod_{l=1}^{L}(1-\psi_{l})^{MsgSize(l,j)}$, which is the probability of the
involved links not to fail raised to the power of the size of messages
traversing them in the scenario (denoted by $MsgSize(l,j)$).

###  II-C Performance Antipatterns

A performance antipattern describes bad design practices that might lead to
performance degradation in a system. Smith and Williams have introduced the
concepts of performance antipatterns in [14]. These textual descriptions were
later translated into first-order logic (FOL) equations [29].

A performance antipattern FOL is a combination of multiple literals, where
each one represents a system aspect (_e.g.,_ the number of connections among
components). These literals must be compared to thresholds in order to reveal
the occurrence of a performance antipattern. The identification of such
thresholds is a non-trivial task, and using deterministic values may result in
an excessively strict detection where the smallest change in the value of a
literal determines the occurrence of the antipattern. For these reasons, we
use the fuzzy threshold concept that has been introduced in [23]. An example
of a performance antipattern fuzzy threshold is the following:

$1-\frac{UB(literal)-literal}{UB(literal)-LB(literal)}$

The upper (UB) and the lower (LB) bounds, in the above equation, are the
maximum and minimum values of the $literal$ computed on the entire system.
Instead of detecting a performance antipattern in a deterministic way, such
thresholds lead to assign probabilities to antipattern occurrences.

In this study we detect the performance antipatterns listed in Table I.

###  II-D Architectural distance

The architectural distance archDist quantifies the distance of the model
obtained by applying refactoring actions to the initial one. On one side, a
_baseline refactoring factor (BRF)_ is associated with each refactoring action
in our portfolio, and expresses the refactoring effort to be spent when
applying the action. On the other side, an _architectural weight (AW)_ is
associated with each model element on the basis of the number of connections
to other elements in the model. Hence, the effort needed to perform a
refactoring is quantified as the product between the _baseline refactoring
factor_ of an action and the _architectural weight_ of the model element on
which that action is applied. archDist is obtained by summing the efforts of
all refactoring actions contained in a solution.

As an example, let us assume that a refactoring sequence is made up of two
refactoring actions: A1 with $BRF(A1)=1.23$, and A2 with $BRF(A2)=2.3$. For
each refactoring action the algorithm randomly selects a target element in the
model. For instance, let those target elements be: E1 with $AW(E1)=1.43$, and
E2 with $AW(E2)=1.32$. The resulting archDist of A1 and A2 would be:

$\textsc{archDist}(A1,A2)=1.23\cdot 1.43+2.3\cdot 1.32$

Details about the _baseline refactoring factor_ for each considered
refactoring action are in Section IV-B.

##  III Approach

Our process starts from an initial UML model and, through refactoring actions,
generates alternative models that are functionally equivalent to the initial
one, but that show different values for the non-functional attributes
considered as objectives. Once an initial population is available, the
evolutionary process starts seeking the optimal Pareto frontier. Figure 1 is a
graphical representation of the approach we present in this paper.

Figure 1: Our multi-objective evolutionary approach

###  III-A Assumptions on UML models

In our approach, we consider UML models including three views, namely _static,
dynamic_ and _deployment_ views. The static view is modeled by a UML Component
diagram in which static connections among components are represented by
interface realizations and their usages. The dynamic view is described by UML
Use Case and Sequence diagrams. A Use Case diagram defines user scenarios,
while a Sequence diagram describes the behavior inside a single scenario
through component operations (as defined in their interfaces) and interactions
among them. A Deployment diagram is used to model platform information and map
Components to Deployment Nodes. As mentioned before, we use an augmented UML
notation by embedding two existing profiles, namely MARTE [17] that expresses
performance concepts, and DAM [18] that expresses reliability concepts.

###  III-B The Refactoring Engine

The automated refactoring of software models is a key point when evolutionary
algorithms are employed in order to optimize some model attributes. For the
sake of full automation of our approach, we have implemented a refactoring
engine that applies refactoring actions on UML models.

Each solution that our evolutionary algorithm produces is a sequence of
refactoring actions that, once applied to an initial model, leads to a model
alternative that shows different non-functional properties. Since our
refactoring actions are combined during the evolutionary approach, we exploit
the feasibility engine that verifies in advance whether a sequence of
refactoring actions is feasible or not [30].

The engine is based on the Cinnéide and Nixon [31] framework. Each refactoring
action is equipped with a pre- and post-condition, which describe the state of
the subject model before and after the application of the action. The engine
reduces a sequence of refactoring actions to a single refactoring action,
which includes all the changes (see Equation (1a)). For example, considering
two refactoring actions ($M_{i}$, and $M_{j}$), then the global pre-condition
is obtained by logical ANDing the first action pre-condition
(${{}^{P_{r}}}M_{i}$) and all the parts of $M_{j}$ pre-condition that are not
yet verified by $M_{i}$ post-conditions ($M_{j}^{P_{r}}\ /\
M_{i}{{}^{P_{o}}}$) (see Equations (1b)). Since the status of the model after
a refactoring is synthesized by its post-condition, we can discard the parts
of a subsequent refactoring pre-condition that, by construction, are already
verified by its post-condition. The global post-condition is obtained by
logical ANDing all post-conditions within the sequence ($M_{i}^{P_{o}}\wedge
M_{j}^{P_{o}}$) (see Equation (1c)).

$\displaystyle{{}^{P_{r}}}M_{i}^{P_{o}}\wedge{{}^{P_{r}}}M_{j}^{P_{o}}\longmapsto{{}^{P_{r}}}M^{P_{o}}$
(1a) $\displaystyle{{}^{P_{r}}}M_{i}\wedge M_{j}^{P_{r}}\ /\
M_{i}{{}^{P_{o}}}\longmapsto{{}^{P_{r}}}M$ (1b) $\displaystyle
M_{i}^{P_{o}}\wedge M_{j}^{P_{o}}\longmapsto M^{P_{o}}$ (1c)

Our feasibility engine also allows to reduce the number of invalid refactoring
sequences, thus reducing the computational time. In the following, we report
the four refactoring actions considered within the optimization process to
generate model alternatives.

#### Clone a Node

This action is aimed at introducing a replica of a Node. Adding a replica
means that every deployed artifact and every connection of the original Node
has to be in turn cloned. Stereotypes and their tagged values are cloned as
well. The rationale of this action is to introduce a replica of a platform
device with the aim of reducing its utilization.

#### Move an Operation to a new Component deployed on a new Node

This action is in charge of randomly selecting an operation and moving it to a
new Component. All the elements related to the moving operation (_e.g.,_
links) will move as well. Since we adopt a multi-view model, and coherence
among views has to be preserved, this action has to synchronize dynamic and
deployment views. A lifeline for the newly created Component is added in the
dynamic view, and messages related to the moved operation are forwarded to it.
In the deployment view, instead, a new Node, a new artifact, and related links
are created. The rationale of this action is to lighten the load of the
original Component and Node.

#### Move an Operation to a Component

This action is in charge of randomly selecting and transferring an Operation
to an arbitrary existing target Component. The action consequently modifies
each UML Use Case in which the Operation is involved. Sequence Diagrams are
also updated to include a new lifeline representing the Component owning the
Operation, but also to re-assign the messages invoking the operation to the
newly created lifeline. The rationale of this action is quite similar to the
previous refactoring action, but without adding a new UML Node to the model.

#### Deploy a Component on a new Node

This action simply modifies the deployment view by redeploying a Component to
a newly created Node. In order to be consistent with the initial model, the
new Node is connected with all other ones directly connected to the Node on
which the target Component was originally deployed. The rationale of this
action is to lighten the load of the original UML Node by transferring the
load of the moving Component to a new UML Node.

###  III-C UML to LQN Transformation

As mentioned before, for a full automation of our approach we have implemented
a model transformation (UML2LQN) that generates a Layered Queuing Network
(LQN) starting from a UML model compliant with our assumptions. The
performance metrics that are required from the evolutionary algorithm come out
from the LQN analysis and are reported back in the UML model. The analysis is
conducted using the LQN solver 222https://github.com/layeredqueuing/V5
provided by Franks _et al._ [32].

The transformation engine is based on the Extensible Platform for Integrated
Languages for model management (Epsilon) [33].

UML2LQN mappings are listed in Table II, and schematically exemplified in
Figure 2. The first column of Table II lists the UML Elements that are
consumed by the transformation, while the third column lists the corresponding
LQN Elements generated by the transformation. The second column, instead,
lists the stereotypes that filter elements within the UML model.

The transformation maps each UML Node to an LQN Processor. A node is included
in the transformation if stereotyped with MARTE _GaExecHost_. The
transformation considers all nodes deploying at least an interacting
component. A component is considered to be interacting if there is at least a
UML Lifeline describing its behavior in at least one of the considered
scenarios (_i.e.,_ UML Use Cases).

Each Sequence Diagram should contain at least a lifeline representing the
element (e.g., a user or another system) that triggers the scenario. This
element is annotated as MARTE _GaWorkloadEvent_. In our context, each scenario
is activated by a UML Actor. Components and Actors are the elements that send
and receive requests in the Sequence Diagram that describes the behavior
originated by triggering a Use Case. For this reason, Components and Actors
are mapped to LQN Tasks. Each lifeline has one or more units of behavior
represented by a UML Behavior Execution Specification (BES). A BES represents
a set of events. Each event stands for the reception or the emission of a
message and is modeled by a UML Message Occurrence Specification (MOS). MOSs
are associated with the requests made in an LQN, and consequently mapped to
LQN _synch-call_. It follows that BESs are mapped to LQN Entries.

A UML Message is transformed into an LQN Activity if stereotyped by MARTE
_GaStep_. Each LQN Entry holds at least an LQN Activity that models a
synchronous call to a different LQN Entry. The _GaStep.rep_ tagged value
denotes the number of requests of an LQN Activity. Instead, the
_GaStep.execTime_ tagged value indicates the service time of a UML Message.
This information is needed to calculate the demand for a call.

In an LQN, a request may be executed by multiple processors. UML2LQN exploits
this feature by calculating the number of UML Nodes deploying the same
component. These nodes are represented, in the obtained LQN, as a single
processor having multiplicity larger than one.

Our approach is similar to the one from Altamini  _et al._ [34], but we could
not adopt the latter because it starts from UML Activity Diagrams, whereas we
introduce Sequence Diagrams to more explicitly express the interactions among
components. Both approaches report the results produced by the LQN solver back
in the original UML model. In Altamini  _et al._ , this is accomplished using
a trace model generated during the transformation. UML2LQN exploits an Epsilon
Object Language (EOL) script that reads the solver output and directly
introduces the results in the original UML model. The generation of a trace
model would have added unnecessary complexity to the transformation. Hence,
given our mapping, it is sufficient to use an EOL script that reports
performance results back into the original model.

Table II: The stereotype application on UML Elements and their mapping to LQN
ones.

UML Element Stereotype LQN Element Model - lqnmodel Node MARTE::GaExecHost
processor CommunicationLink DAM::DaConnector - Component DAM::DaComponent task
Actor MARTE::GaWorkloadEvent task Behavior Execution Specification - entry
Message MARTE::GaStep activity Message Occurrence Specification - synch-call
Use Case MARTE::GaScenario -

Figure 2: Graphical representation of the UML to LQN mapping.

###  III-D Computing reliability on UML models

The reliability parameters of the model introduced in Section II-B are
annotated on UML models by means of the MARTE-DAM profile. The probability of
executing a scenario ($p_{j}$) is specified by annotating UML Use Cases with
the _GaScenario_ stereotype. This stereotype has a tag named _root_ that is a
reference to the first _GaStep_ in a sequence. We use the _GaScenario.root_
tag to point to the triggering UML Message of a Sequence Diagram and the
_GaStep.prob_ to set the execution probability. Failure probabilities of
components ($\theta_{i}$) are defined by applying the _DaComponent_ stereotype
on each UML Component and by setting, in the _failure_ tag, a _DaFailure_
element with the failure probability specified in the _occurrenceProb_ tag.
Analogously, failure probabilities of links ($\psi_{l}$) are defined in the
_failure.occurrenceProb_ tag of the _DaConnector_ stereotype that we apply on
UML CommunicationPath elements. Such elements represent the connection links
between UML Nodes in a Deployment Diagram. Sequence Diagrams are traversed to
obtain the number of invocations of a component $i$ in a scenario $j$ (denoted
by $InvNr_{ij}$ in our reliability model), but also to compute the total size
of messages passing over a link $l$ in a scenario $j$ (denoted by
$MsgSize(l,j)$). The size of a single UML Message is annotated using the
_GaStep.msgSize_ tag. The Java implementation of the reliability model is
available online333https://github.com/SEALABQualityGroup/uml-reliability.

##  IV Evaluation

In this section, we apply our approach to the Train Ticket Booking Service
(TTBS) case study [24, 25]. We report results of Pareto frontiers shapes, with
a special focus on the role of performance antipatterns in the solution
quality444Due to lack of space, we evaluated our approach on a single case
study. However, its characteristics allowed a non-trivial validation
exercise..

###  IV-A Case Study

TTBS is a web-based booking application, and its architecture is based on the
microservice paradigm. The system is made up of 40 microservices, and it
provides different scenarios through users that can perform realistic
operations, _e.g.,_ book a ticket or watch trip information like intermediate
stops. The application employs a docker container for each microservice, and
connections among them are managed by a central pivot container.

Our UML model of TTBS is available
online.555https://github.com/SEALABQualityGroup/SEAA-replication-
package/tree/main/case-study The static view is made of 11 UML Components,
where each component represents a microservice. In the deployment view,
instead, we consider 11 UML Nodes, each one representing a docker container.

Among all available TTBS scenarios shown in [24], in this paper we have
considered 3 UML Use Cases, namely _login_ , _update user details_ and
_rebook_. We selected these three scenarios because they commonly represent
performance-critical ones in a ticketing booking service. In particular, each
scenario is described by a UML Sequence Diagram. Furthermore, the model
comprises two user categories: simple and admin users. The simple user
category can perform the login and the rebook scenarios, while the admin
category can perform the login and the update user details scenarios.

###  IV-B Experimental setting

A configuration is defined by the combination of parameters related to the
genetic algorithm, the performance antipattern detection, and the refactoring
engine. In order to investigate which configuration produces better Pareto
frontiers, we have executed multiple tuning runs to find a set of optimal
configurations.

Table III: The baseline refactoring factor for actions belonging to our portfolio. Refactoring action | BRF
---|---
Clone UML Node | 1.23
Move Operation to a New Component to a new Node | 1.80
Move Operation to a Component | 1.64
Deploy a Component to a new Node | 1.45

To set parameters related to the genetic algorithm, we have performed a tuning
phase with the intent of increasing the quality of the Pareto frontiers. In
particular, we have set the length of refactoring sequences to four actions,
which represents a good compromise on the number of refactoring actions
usually applied by a designer in a single session. We have set the
$P_{crossover}$ and $P_{mutation}$ probabilities to 0.8 and 0.2, respectively,
following common configurations [35]. The higher the values of these two
probabilities, the greater the chance of generating an unfeasible sequence of
refactoring actions, which in turn causes a longer simulation time due to a
higher number of discarded sequences. For example, increasing $P_{crossover}$
could cause a lot of permutation among sequences, and it might lead to wrong
or unfeasible sequences of refactoring actions. The initial population size
has been set to 16 elements (i.e., 16 different UML alternative models).We
have chosen this size for the initial population because it represents, for
this case study, a good compromise between the diversity of solutions
considered and the time required to create them.

Also, we have considered {0.55, 0.80, 0.95} as fuzziness thresholds to study
the impact of performance antipatterns on Pareto frontiers. Regarding
parameters related to refactoring actions, we have set the BRF of each
refactoring action as reported in Table III. Since a tuning phase of the BRF
is out of the scope of this paper, we did not investigate the impact of other
BRF values on the overall quality. Thus, for the sake of the evaluation, we
have chosen the values reported in the table, which represent the effort for
applying those actions onto a model.

Furthermore, in order to reduce the randomness of the genetic algorithm, a
variable number of multiple runs have been executed for each configuration. In
the following, a Reference Pareto Frontier (RPF) refers to those solutions
outcoming from the same problem configuration that are not dominated by other
ones.

Our experimental settings on the TTBS case study has generated 31,404 model
alternatives and it has taken 82.4 hours of computation. We performed our
experiments on a server equipped with two Intel Xeon E5-2650 v3 CPUs at
2.30GHz, 40 cores and 80GB of RAM. The replication package including the Java
implementation of the approach, the case study and the results of our
experimentation are available
online666https://github.com/SEALABQualityGroup/SEAA-replication-package.

###  IV-C Results and discussion

Results presented in this section are aimed at answering the aforementioned
three research questions.

RQ1: Does antipattern detection contribute to find better solutions compared
to the case where antipatterns are not considered at all?

In order to answer this research question, we have conducted additional
experimentation on the same problem configurations, where we have removed
performance antipattern occurrences from the fitness function, thus reducing
the optimization to the remaining three objectives.

Figure 3: Reference Pareto Fronts obtained when considering performance
antipatterns (with PAs) and without them (without PAs).

Figure 3 reports the solutions in Reference Pareto Fronts across different
problem configurations, where two objectives are represented on the axes and
the architectural distance is represented by the symbol size of each solution.
Hence, better solutions are the small sized ones located in the upper right
corner of the figure. Two RPFs are reported in the figure: gray diamonds are
solutions without PAs, whereas pink triangles with PAs, with a 0.95 fuzziness
value. The initial model is identified by a black x.

Figure 3 provides clear evidence of higher quality solutions in the case with
PAs, where antipattern detection helps to drive the evolutionary algorithm
towards solutions that achieve perfQ improvements even twice higher than the
ones without PAs (i.e., 0.25 vs 0.12).

On the basis of our experimentation, we can state that considering the
performance antipattern occurrences in the optimization leads to better
solutions than ignoring them.

RQ2: Does the probabilistic nature of fuzzy antipatterns detection help to
include higher quality solutions in Pareto frontiers with respect to
deterministic one?

In order to answer this research question, we run different problems while
varying the values of the fuzziness threshold of the performance antipattern
detection within {0.50, 0.80, 0.95}. The obtained RPFs are reported in Figure
4.

Figure 4: Reference Pareto Fronts obtained when varying the fuzziness
threshold of the performance antipatterns detection.

The figure shows that moving from blue squared to green round solutions (i.e.,
in an increasing fuzziness value direction) already leads better solutions, at
least in terms of reliability. By further moving from green round to pink
triangle solutions it is possible to identify even better solutions,
especially in terms of perfQ performance objective.

Hence, it seems that having a more deterministic antipattern detection (i.e.,
higher values of fuzziness) is better than a probabilistic one. However, the
deterministic detection has the drawback of relying on fixed thresholds that
must be computed in advance for each model alternative. For this reason, the
trade-off between better quality solution and the effort to bind thresholds is
likely domain-dependent and worth to be more investigated.

On the basis of our experimentation, we can state that performance antipattern
fuzzy detection does not help to improve the quality of Pareto frontiers.

Table IV: Details of the values of perfQ, reliability and archDist obtained when varying the fuzziness threshold of the performance antipatterns detection. | perfQ | Reliability | archDist
---|---|---|---
Fuzziness | # sol | min | max | median | mean | min | max | median | mean | min | max | median | mean
.55 | 16 | -0.038 | 0.116 | 0.052 | 0.041 | 0.302335 | 0.847111 | 0.787196 | 0.713601 | 1.85 | 10.41 | 6.06 | 5.82
.80 | 25 | -0.058 | 0.116 | 0.053 | 0.040 | 0.309063 | 0.847808 | 0.788695 | 0.740425 | 1.97 | 9.18 | 4.83 | 4.96
.95 | 13 | 0.000 | 0.256 | 0.138 | 0.134 | 0.549760 | 0.847111 | 0.785697 | 0.717795 | 3.03 | 14.67 | 5.49 | 6.56
no PAs | 29 | -0.058 | 0.116 | 0.066 | 0.056 | 0.591957 | 0.872976 | 0.788695 | 0.778037 | 2.46 | 9.75 | 5.40 | 5.88

RQ3: Is the approach able to keep the system reliability satisfactorily high?

Figure 4 shows that, regardless of the fuzziness value, most solutions improve
reliability with respect to the initial model (i.e., the black x marker in the
figure). When looking at the mean values of reliability in Table IV, it is
clear that, even if in few cases the reliability of solutions decreases, on
average we obtain an improvement from 8% to 18% (in absolute terms) with
respect to the value of the initial solution, _i.e.,_ 0.657925. Moreover, a
median higher than the mean also confirms that solutions with improved
reliability are more likely to be part of a Pareto frontier. In the best
cases, the proposed refactoring actions are able to produce alternative models
that achieve an improvement in reliability by over 28% (in absolute terms).

Our experiments show that, in the majority of cases, the approach is able to
derive solutions with satisfactorily high reliability.

Table IV reports a summary of the values of perfQ, reliability and archDist of
solutions within different Pareto frontiers, both when performance
antipatterns have been excluded and when a fuzzy detection is employed. In
particular, in order to achieve the highest performance and reliability
together, the algorithm has combined refactoring actions that, in some cases,
have produced quite poor solutions in terms of architectural distance (_i.e.,_
14.67), as also shown by large triangles in the upper right corner of the
figure. If, instead, we move towards lower values of reliability objective
than the initial architecture one (_i.e.,_ towards the lower right side of the
figure), then the algorithm has generated closer models in terms of
architectural distance in the range of {1.85, 4.96}, as also shown by small
triangles in the $0.5$ to $0.7$ range on y-axis of Figure 4). Therefore,
solutions that are able to achieve high performance and reliability usually
show larger values of the distance metric, whereas the ones where either
performance or reliability is pursued exhibit lower distance values.

###  IV-D Threats to validity

We follow the Wholin et al. [36] classification for the following threats.

#### Conclusion validity

Our results might be affected by _Conclusion validity_ threats, since our
considerations might change with better-tuned parameters for the
$\mathtt{NSGA-II}$. Also, problem parameter configurations might threat our
conclusion. We did not perform an extensive tuning phase for the latter due to
the long duration of each run, while we used common parameters for the
$\mathtt{NSGA-II}$, which should mitigate these threats.

#### Construct validity

The way we have designed our problem and our experimentation might be affected
by _Construct validity_ threats. In particular, the role played by the
architectural distance objective on the combination of refactoring actions
might affect the refactoring actions selection.

#### External validity

Our approach might be affected by _external validity_ threats, because we have
used a single modelling notation. We cannot generalize our results to other
modelling notations, which could imply using a different portfolio of
refactoring actions. In fact, the syntax and semantics of the modelling
notation determine the amount and nature of refactoring actions that can be
performed. We could mitigate these threats, for example, by using another
modelling notation.

##  V Related Work

In the last decade, software architecture multi-objective optimization studies
have been introduced to optimize various quality attributes (e.g.,
reliability, and energy [9, 21, 37, 38]) with different degrees of freedom to
modify the architecture (e.g., service selection [39, 40]). A systematic
literature review on architecture optimization can be found in [8]. We
consider here as related work those approaches that directly involve multi-
objective evolutionary algorithms and, on the other hand, approaches that
exploit LQN as performance-oriented modelling notation [41, 42, 43, 12].

###  V-A Software Architecture optimization

Menasce _et al._ have presented a whole framework for architectural design and
quality optimization [44], where architectural patterns are used to support
the searching process (e.g., load balancing, fault tolerance). Two limitations
affects the approach: the architecture has to be designed in a tool-related
notation and not in a standard modelling language (as we do in this paper),
and it uses equation-based analytical models for performance indices that
could be too simple to capture architectural details and resource contention.

Aleti _et al._ [10] have presented an approach for modeling and analyzing ADDL
architectures [45]. They have also introduced a tool aimed at optimizing
different quality attributes while varying the architecture deployment and the
component redundancy. Instead, our work relies on UML models and offers more
complex refactoring actions as well as different target attributes for the
fitness function. Besides, we investigate the role of performance antipatterns
in the context of multi-objective software architecture refactoring
optimization.

In the context of software architecture optimization, Cortellessa and Di
Pompeo studied the sensitivity of multi-objective software architecture
refactoring to configuration characteristics [22]. They compared two genetic
algorithms in terms of Pareto frontiers quality dealing with architectures
defined in AEmilia, which is a performance-oriented ADL. In this paper, we
change the modelling notation from AEmilia to UML, and we add the reliability
as a new objective. Both approaches provide a refactoring engine, however, in
this paper, the refactoring engine offers more complex refactoring actions
since UML is more expressive than AEmilia.

###  V-B Layered Queueing Network approaches

Koziolek _et al._ have presented PerOpteryx [9, 12], a performance-oriented
multi-objective optimization problem. In PerOpteryx the optimization process
is guided by tactics referring to component reallocation, faster hardware, and
more hardware, which do not represent structured refactoring actions, as we
intend to do in this paper. Moreover, PerOpteryx supports architectures
specified in Palladio Component Model [46] and produces, through model
transformation, a LQN for of performance analysis.

Rango _et al._ have presented SQuAT [47], an extensible platform aimed at
including flexibility in the definition of an architecture optimization
problem. SQuAT supports models conforming to Palladio Component Model
language, exploits LQN for performance evaluation, and PerOpteryx tactics for
architecture.

We differ to PerOpteryx and SQuAt because we use UML as modelling notation.
Therefore, our approach could be applied from a large audience of software
designers since we rely on a standard notation instead of a specific language.

Model-to-model (M2M) transformations from UML to LQN notations have been
presented in [41, 42, 43, 34]. In contrast with these approaches, we present a
novel M2M transformation mapping that employs UML Sequence Diagrams as the
behavioral view of software architectures instead of UML Activity Diagrams.
Considering UML Sequence Diagrams has two benefits: they are adopted more
frequently than UML Activity Diagrams during the software design [48], and
they explicitly define method calls, while UML Activity Diagrams usually focus
on workflows and processes. Therefore, our approach supports a more detailed
behavioral representation in terms of time intervals between method calls.

##  VI Conclusions and Future Work

In this work, we used $\mathtt{NSGA-II}$ to optimize UML models with respect
to performance and reliability properties, as well as the number of detected
performance antipatterns and the architectural distance. We focused our study
on the impact that performance antipatterns may have on the quality of
refactoring solutions resulting from such optimization setting.

From our experimentation, we gathered interesting insights about the quality
of the generated solutions and the role of perfomance antipatterns as an
objective of the algorithm. In this regard, we showed that, by including the
detection of performance antipatterns in the optimization process, we are able
to obtain better solutions in terms of performance and reliability. Moreover,
we also showed that, the more we increase the probability of detecting a
performance antipattern using the fuzziness threshold, the better the quality
of the refactoring solutions. Another important aspect of our study was to
ensure that the reliability of solutions remained within satisfactory levels.
In this respect, our experiments showed that we were in fact able to increase
the reliability of refactoring solutions, with respect to the initial model,
in the majority of cases.

As future work, we intend to tackle the threats to validity discussed before.
In particular, we intend to investigate the influence of the experimental
settings on the quality of Pareto frontiers. Furthermore, we are interested in
the role played by the architectural distance, and specifically in studying
the effect of the baseline refactoring factor on the combination of
refactoring actions. In the same context, alternative models and additional
parameters to evaluate the architectural distance can also be considered.
Similarly, for the reliability objective, we intend to include fault tolerance
refactoring actions [49] in the refactoring portfolio, and to extend the
reliability model to also take into account error propagation [50].

## References

  * [1] A. Ramírez, J. R. Romero, and S. Ventura, “A survey of many-objective optimisation in search-based software engineering,” _JSS_ , vol. 149, pp. 382–395, 2019.
  * [2] T. Mariani and S. R. Vergilio, “A systematic review on search-based refactoring,” _JIST_ , vol. 83, pp. 14–34, Mar. 2017.
  * [3] A. Ouni, M. Kessentini, and K. Inoue, “Search-based Web Service Antipatterns Detection,” _TEVC_ , vol. 10, no. 4, pp. 603–617, 2017.
  * [4] A. Ouni, M. Kessentini, K. Inoue, and M. Ó Cinnéide, “Search-Based Web Service Antipatterns Detection,” _TSC_ , vol. 10, no. 4, pp. 603–617, 2017.
  * [5] M. Ray and D. P. Mohapatra, “Multi-objective test prioritization via a genetic algorithm,” _Innov. Syst. Softw. Eng._ , vol. 10, no. 4, pp. 261–270, 2014.
  * [6] G. Bavota, M. Di Penta, and R. Oliveto, “Search Based Software Maintenance: Methods and Tools,” in _Evol. Soft. Sys._ , 2014, pp. 103–137.
  * [7] M. Kessentini, H. Sahraoui, M. Boukadoum, and O. Ben Omar, “Search-based model transformation by example,” _SOSYM_ , vol. 11, no. 2, pp. 209–226, 2012.
  * [8] A. Aleti, B. Buhnova, L. Grunske, and I. Meedeniya, “Software Architecture Optimization Methods: A Systematic Literature Review,” _TSE_ , vol. 39, no. 5, pp. 658–683, 2013.
  * [9] A. Martens, H. Koziolek, S. Becker, and R. H. Reussner, “Automatically improve software architecture models for performance, reliability, and cost using evolutionary algorithms,” in _ICPE_ , 2010, pp. 105–116.
  * [10] A. Aleti, S. Björnander, L. Grunske, and I. Meedeniya, “ArcheOpterix: An extendable tool for architecture optimization of AADL models,” in _ICSE MOMPES Workshop_ , 2009, pp. 61–71.
  * [11] C. Blum and A. Roli, “Metaheuristics in combinatorial optimization: Overview and conceptual comparison,” _ACM Comput. Surv._ , vol. 35, no. 3, pp. 268–308, 2003.
  * [12] A. Koziolek, H. Koziolek, and R. H. Reussner, “PerOpteryx: automated application of tactics in multi-objective software architecture optimization,” in _QoSA_ , 2011, pp. 33–42.
  * [13] Y. Ni, X. Du, P. Ye, L. L. Minku, X. Yao, M. Harman, and R. Xiao, “Multi-objective software performance optimisation at the architecture level using randomised search rules,” _JIST_ , vol. 135, p. 106565, 2021.
  * [14] C. U. Smith and L. G. Williams, “Software performance antipatterns; common performance problems and their solutions,” in _CMG Conf._ , 2001, pp. 797–806.
  * [15] ——, “More New Software Performance Antipatterns: Even More Ways to Shoot Yourself in the Foot,” in _CMG Conf._ , 2003, pp. 717–725.
  * [16] D. Arcelli, V. Cortellessa, and D. Di Pompeo, “Performance-driven software model refactoring,” _JIST_ , vol. 95, pp. 366–397, 2018.
  * [17] O. M. Group, “A UML profile for MARTE: modeling and analysis of real-time embedded systems,” Object Management Group, 2008.
  * [18] S. Bernardi, J. Merseguer, and D. C. Petriu, “A dependability profile within MARTE,” _Softw. Syst. Model._ , vol. 10, no. 3, pp. 313–336, 2011.
  * [19] V. Cortellessa, H. Singh, and B. Cukic, “Early reliability assessment of UML based software models,” in _WOSP@ISSTA_ , 2002, pp. 302–309.
  * [20] D. Arcelli, V. Cortellessa, M. D’Emidio, and D. Di Pompeo, “EASIER: an Evolutionary Approach for multi-objective Software archItecturE Refactoring,” in _ICSA_ , 2018, pp. 1–10.
  * [21] R. Li, R. Etemaadi, M. T. M. Emmerich, and M. R. V. Chaudron, “An evolutionary multiobjective optimization approach to component-based software architecture design,” in _CEC_. IEEE, 2011, pp. 432–439.
  * [22] V. Cortellessa and D. Di Pompeo, “Analyzing the sensitivity of multi-objective software architecture refactoring to configuration characteristics,” _JIST_ , vol. 135, p. 106568, 2021.
  * [23] D. Arcelli, V. Cortellessa, and C. Trubiani, “Performance-based software model refactoring in fuzzy contexts,” in _FASE_ , 2015, pp. 149–164.
  * [24] D. Di Pompeo, M. Tucci, A. Celi, and R. Eramo, “A microservice reference case study for design-runtime interaction in MDE,” in _STAF MDE@DeRun Workshop_ , vol. 2405, 2019, pp. 23–32.
  * [25] X. Zhou, X. Peng, T. Xie, J. Sun, C. Ji, W. Li, and D. Ding, “Fault analysis and debugging of microservice systems: Industrial survey, benchmark system, and empirical study,” _TSE_ , vol. 47, no. 2, pp. 243–260, 2021.
  * [26] K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan, “A fast and elitist multiobjective genetic algorithm: NSGA-II,” _TEVC_ , vol. 6, no. 2, pp. 182–197, Apr. 2002.
  * [27] U. Mansoor, M. Kessentini, M. Wimmer, and K. Deb, “Multi-view refactoring of class and activity diagrams using a multi-objective evolutionary algorithm,” _Software Quality Journal_ , pp. 1–29, 2015.
  * [28] G. Franks, T. Al-Omari, M. Woodside, O. Das, and S. Derisavi, “Enhanced Modeling and Solution of Layered Queueing Networks,” _TSE_ , vol. 35, no. 2, pp. 148–161, 2009.
  * [29] V. Cortellessa, A. Di Marco, and C. Trubiani, “An approach for modeling and detecting software performance antipatterns based on first-order logics,” _Softw. Syst. Model._ , vol. 13, no. 1, pp. 391–432, 2014.
  * [30] D. Arcelli, V. Cortellessa, and D. Di Pompeo, “A metamodel for the specification and verification of model refactoring actions,” in _Proceedings of the 2nd International Workshop on Refactoring, IWoR@ASE 2018, Montpellier, France, September 4, 2018_ , A. Ouni, M. Kessentini, and M. Ó. Cinnéide, Eds. IWoR@ACM, 2018, pp. 14–21.
  * [31] M. Ó Cinnéide and P. Nixon, “Composite refactorings for Java programs,” in _Formal Techniques for Java Programs_ , 2000, pp. 129–135.
  * [32] G. Franks, P. Maly, M. Woodside, D. Petriu, and A. Hubbard, _Layered Queueing Network Solver and Simulator User Manual_ , 2013.
  * [33] D. Kolovos, L. Rose, R. Paige, and A. Garcıa-Domınguez, “The epsilon book,” _Structure_ , vol. 178, pp. 1–10, 2010.
  * [34] T. Altamimi, M. H. Zargari, and D. C. Petriu, “Performance analysis roundtrip: automatic generation of performance models and results feedback using cross-model trace links,” in _CASCON_ , 2016, pp. 208–217.
  * [35] A. Arcuri and G. Fraser, “Parameter tuning or default values? an empirical investigation in search-based software engineering,” _Empirical Software Engineering_ , vol. 18, no. 3, pp. 594–623, 2013.
  * [36] C. Wohlin, P. Runeson, M. Höst, M. C. Ohlsson, B. Regnell, and A. Wesslén, _Experimentation in software engineering_. Springer Science & Business Media, 2012.
  * [37] I. Meedeniya, B. Buhnova, A. Aleti, and L. Grunske, “Architecture-Driven Reliability and Energy Optimization for Complex Embedded Systems,” in _QoSA_. Springer, 2010, pp. 52–67.
  * [38] A. Martens, D. Ardagna, H. Koziolek, R. Mirandola, and R. H. Reussner, “A Hybrid Approach for Multi-attribute QoS Optimisation in Component Based Software Systems,” in _Research into Practice – Reality and Gaps_ , 2010, pp. 84–101.
  * [39] F. Rosenberg, M. B. Müller, P. Leitner, A. Michlmayr, A. Bouguettaya, and S. Dustdar, “Metaheuristic Optimization of Large-Scale QoS-aware Service Compositions,” in _SCC_ , 2010, pp. 97–104.
  * [40] V. Cardellini, E. Casalicchio, V. Grassi, F. Lo Presti, and R. Mirandola, “QoS-driven Runtime Adaptation of Service Oriented Architectures,” in _ESEC/FSE_ , 2009, pp. 131–140.
  * [41] C. M. Woodside, D. C. Petriu, D. B. Petriu, H. Shen, T. Israr, and J. Merseguer, “Performance by unified model analysis (PUMA),” in _WOSP_ , 2005, pp. 1–12.
  * [42] C. Li, T. Altamimi, M. H. Zargari, G. Casale, and D. C. Petriu, “Tulsa: A tool for transforming UML to layered queueing networks for performance analysis of data intensive applications,” in _QEST_ , vol. 10503, 2017, pp. 295–299.
  * [43] T. Altamimi and D. C. Petriu, “Incremental change propagation from UML software models to LQN performance models,” in _CASCON_ , 2017, pp. 120–131.
  * [44] D. A. Menascé, J. M. Ewing, H. Gomaa, S. Malek, and J. P. Sousa, “A framework for utility-based service oriented design in SASSY,” in _WOSP/SIPEW_ , 2010, pp. 27–36.
  * [45] P. H. Feiler and D. P. Gluch, _Model-Based Engineering with AADL - An Introduction to the SAE Architecture Analysis and Design Language_ , ser. SEI series in software engineering. Addison-Wesley, 2012.
  * [46] S. Becker, H. Koziolek, and R. H. Reussner, “The Palladio component model for model-driven performance prediction,” _Systems and Software_ , vol. 82, no. 1, pp. 3–22, Jan. 2009.
  * [47] A. Rago, S. Vidal, J. A. Diaz-Pace, S. Frank, and A. van Hoorn, “Distributed quality-attribute optimization of software architectures,” in _SBCARS_ , 2017\.
  * [48] J. Erickson and K. Siau, “Can UML be simplified? Practitioner use of UML in separate domains,” in _EMMSAD_ , vol. 7, 2007, pp. 87–96.
  * [49] V. Cortellessa, R. Eramo, and M. Tucci, “From software architecture to analysis models and back: Model-driven refactoring aimed at availability improvement,” _Inf. Softw. Technol._ , vol. 127, p. 106362, 2020.
  * [50] V. Cortellessa and V. Grassi, “A modeling approach to analyze the impact of error propagation on reliability of component-based systems,” in _CBSE_ , vol. 4608, 2007, pp. 140–156.


# Sparse Sampling in Fractional Fourier Domain:
Recovery Guarantees and Cramér–Rao Bounds

Václav Pavlíček and Ayush Bhandari A. Bhandari’s work is supported by the UK
Research and Innovation council’s _Future Leaders Fellowship_ program “Sensing
Beyond Barriers” (MRC Fellowship award no. MR/S034897/1).The authors are with
the Dept. of Electrical and Electronic Engineering, Imperial College London,
South Kensington, London SW7 2AZ, UK. (Emails:
<EMAIL_ADDRESS>and ayush@alum.mit.edu).Manuscript submitted
on January 2, 2024; accepted April 20, 2024. Date of the current version April
29, 2024.

###### Abstract

Sampling theory in fractional Fourier Transform (FrFT) domain has been studied
extensively in the last decades. This interest stems from the ability of the
FrFT to generalize the traditional Fourier Transform, broadening the
traditional concept of bandwidth and accommodating a wider range of functions
that may not be bandlimited in the Fourier sense. Beyond bandlimited
functions, sampling and recovery of sparse signals has also been studied in
the FrFT domain. Existing methods for sparse recovery typically operate in the
transform domain, capitalizing on the spectral features of spikes in the FrFT
domain. Our paper contributes two new theoretical advancements in this area.
First, we introduce a novel time-domain sparse recovery method that avoids the
typical bottlenecks of transform domain methods, such as spectral leakage.
This method is backed by a sparse sampling theorem applicable to arbitrary
FrFT-bandlimited kernels and is validated through a hardware experiment.
Second, we present Cramér–Rao Bounds for the sparse sampling problem,
addressing a gap in existing literature.

###### Index Terms:

Annihilation, Cramér–Rao Bounds, Fractional Fourier Transform, Sparse
Sampling.

###### Contents

  1. I Introduction
  2. II Mathematical Preliminaries
  3. III Recovery Algorithm and Guarantees
  4. IV Cramér–Rao Bounds for Sparse Sampling
  5. V Conclusion

##  I Introduction

The early 20th century saw the emergence of the fractional Fourier transform
(FrFT) as a significant development in harmonic analysis, pioneered by the
works of Wiener [1] and Condon [2]. The FrFT generalizes the traditional
Fourier transform (FT) by introducing the concept of _fractionalization_
explained next. Consider a time-domain function $f(t)$ and its FT,
$\mathscr{F}_{f}(\omega)=\widehat{f}(\omega)$. Let
$\mathscr{F}^{\left(n\right)}_{f}=\mathscr{F}^{\left(n-1\right)}_{f}\circ\mathscr{F}_{f},n\geqslant
1$ define the recursive application of the FT operator, where $\circ$ denotes
operator composition and $\mathscr{F}^{\left(0\right)}_{f}=\mathrm{I}$
(identity). It can be verified via basic Fourier analysis that,

$\begin{array}[]{*{20}{c}}\underbrace{\mathscr{F}^{\left(1\right)}_{f}=\widehat{f}}_{n=1}&\underbrace{\mathscr{F}^{\left(2\right)}_{f}=\widetilde{f}}_{n=2}&\underbrace{\mathscr{F}^{\left(3\right)}_{f}=\widetilde{\widehat{f}}}_{n=3}&\underbrace{\mathscr{F}^{\left(4\right)}_{f}=\mathscr{F}^{\left(0\right)}_{f}}_{n=4}\end{array},$

where $\widetilde{f}(x)=f(-x)$. Iterating the Fourier operator shows that
$\\{\mathscr{F}^{\left(n\right)}_{f}\\}_{n\geqslant 1}$ is 4-periodic
automorphism or the FT is _cyclic on a group of four_ [2]. This periodicity of
the Fourier operator is visually depicted in Fig. 1.

$\theta$$\underbrace{\mathscr{F}_{f}^{(\theta)}=\widehat{f}_{\theta}(\omega)}_{\mathsf{FrFT\
Domain}}$$\mathscr{F}_{f}^{(0)}=f(t)$$\mathscr{F}_{f}^{(2)}=f(-t)$$\mathscr{F}_{f}^{(3)}=\widehat{f}(-\omega)$$\mathscr{F}_{f}^{(1)}=\widehat{f}(\omega)$$=\mathscr{F}_{f}^{(4)}$TimeFrequency

Figure 1: Visual illustration of “fractionalizing” the Fourier transform (FT)
which introduces the FrFT for non-integer orders, $\theta\in\mathbb{R}$. While
the traditional FT is a 4-periodic automorphism [2] (applying FT 4 times
returns the original function), this property holds only for integer orders,
$n\in\mathbb{Z}$. The FrFT generalizes this concept to arbitrary real orders
$\theta$. When $\theta=\pi/2$, the FrFT simplifies back to the FT.

Given this cyclical nature of the FT operator, an intriguing question then is
how to define a Fourier-like mapping on a continuous circle. Condon’s work [2]
provides the answer, introducing an integral transformation that applies to
any arbitrary point $\theta\in\mathbb{R}$ on the circle. This leads to the
mathematical definition of the FrFT which is defined as a mapping,

$\mathscr{F}^{\left(\theta\right)}_{f}:f\to\widehat{f}_{\theta},\qquad\widehat{f}_{\theta}\left(\omega\right)=\int{f\left(t\right)\kappa_{\theta}^{\ast}\left({t,\omega}\right)dt},$
(1)

where the FrFT kernel parameterized by $\theta\in\mathbb{R}$ is written as,

${\kappa_{\theta}}\left({t,\omega}\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\begin{cases}K_{\theta}{e^{-\jmath\left({\frac{{{t^{2}}+{\omega^{2}}}}{2}\cot\theta-\frac{{\omega
t}}{{\sin\theta}}}\right)}}&\theta\not=n\pi\\\
\delta\left({t-\omega}\right)&\theta=2n\pi\\\
\delta\left({t+\omega}\right)&\theta+\pi=2n\pi,\\\ \end{cases}$ (2)

and $K_{\theta}$ is a known constant. The inverse-FrFT is the FrFT with
$-\theta$, which can be interpreted as rotation in $t$–$\omega$ plane (see
Fig. 1). Clearly, at $\theta=\pi/2$, FrFT simplifies to the FT.

Since its inception, the FrFT has been widely applied in various fields. For
comprehensive coverage, the reference book [3] by Ozaktas & co-workers is an
essential resource. Notably, the FrFT offers advantages in digital
communications [4] and radar systems [5], where chirp-like basis functions
(see (2)) are more effective than traditional sinusoidal basis (Fourier).
Other areas of application include quantum mechanics [6], harmonic analysis
[7, 8], time-frequency representations [9, 10], and optics and imaging [11].

In the signal processing and harmonic analysis fields, a fundamental interest
lies in understanding the connection between continuous and discrete
representations through sampling theory. In the context of the FrFT, key focus
has been on the following two fronts.

Shannon’s Sampling. This involves understanding how to represent bandlimited
functions in the FrFT domain as point-wise samples, considering that these
functions are typically non-bandlimited in the FT domain. Several proofs of
the FrFT sampling theorem have been presented [12, 13, 14, 15, 16, 17, 18, 19,
20]. Additionally, sampling theory has expanded beyond bandlimited spaces to
include sparse signals [21, 22, 23] and shift-invariant spaces [19, 20], and
more recently to the Unlimited Sensing Framework [24, 25, 26]. The main result
in the area is as follows.

######  Theorem 1 (Sampling Theorem for FrFT Domain [12, 13, 14, 15, 16, 17,
18, 19, 20]).

If the function $f(t)$ contains no frequencies higher than $\Omega$ in the
fractional Fourier domain, then $f(t)$ is completely determined by its
equidistant samples spaced $T\leqslant\left(\pi\sin\theta\right)/\Omega$.

Note that, at $\theta=\pi/2$, Theorem 1 reduces to the Shannon’s sampling
theorem for the Fourier domain. Recovery of the continuous-time function from
samples $\\{f\left(nT\right)\\}_{n\in\mathbb{Z}}$ is facilitated by the FrFT
interpolation formula,

$f\left(t\right)={e^{-\jmath\frac{{{t^{2}}}}{2}\cot\theta}}\sum\limits_{{n}\in\mathbb{Z}}{f\left({{n}T}\right){e^{\jmath\tfrac{{{{\left({{n}T}\right)}^{2}}}}{2}\cot\theta}}}\operatorname{sinc}\left({\frac{t}{T}-{n}}\right).$

Discretization of Continuous Transforms via Sampling. Here, the focus has been
to achieve discrete representations for the continuous transform. Several
papers (see [27, 28, 29], and references therein) have worked towards a
definition of the discrete fractional Fourier transform (DFrFT), extending the
discrete Fourier transform (DFT) in a manner analogous to how the continuous
FrFT expands upon the continuous FT.

Motivation and Contributions. Sampling theory in the FrFT domain is a
prominent research area, as evidenced by the density of papers on this theme
[12, 13, 14, 15, 16, 17, 18, 19, 20]. Within this framework, the sampling and
reconstruction of _sparse signals_ [21, 22, 23] has been studied in a variety
of contexts due to the pervasive presence of such signals [30, 31] in diverse
application areas. The generalization of sparse recovery methods beyond
Fourier domain [23, 29] is an interesting study in its own right. Current
methods [21, 22, 23] suffer from two research gaps.

  1. $\bullet$

Sparse recovery methods rely on transform domain representations, where an
unknown sparse signal corresponds to a sum-of-complex exponentials in DFrFT
representation [21, 22, 23], with recovery depending on spectral estimation
methods [30, 31]. However, this approach is limited due to spectral leakage in
discrete transforms like DFrFT [27, 28, 29] and DFT, especially in sub-Nyquist
sampling scenarios. A time-domain method would not suffer with such
limitations.

  2. $\bullet$

Current works [21, 22, 23] are focused on sampling criterion with ideal
measurements. Study of performance gurantees, a key tool for benchmarking
algorithmic performance in the presence of noise, has not been discussed in
the literature.

This paper covers the theoretical aspects of the above research gaps. Our main
contributions are as follows:

  1. $\bullet$

We present a novel, time-domain sparse sampling approach that does not require
any DFrFT operations and hence, its side-effects _e.g._ spectral leakage, are
eliminated. Our time-domain recovery approach is backed by a sampling theorem.

  2. $\bullet$

We derive Cramér–Rao Bounds (CRB) for the sparse parameter estimation problem.
This serves as a performance guarantee for the recovery problem in the
presence of noise.

##  II Mathematical Preliminaries

Leading up to our main results, we will first revisit the basic mathematical
operations associated with the FrFT domain.

######  Definition 1 (Chirp Modulation).

For a given FrFT order $\theta\in\mathbb{R}$, we define the quadratic phase,
chirp modulation function by
${\xi_{\theta}}\left(t\right)={\exp\left(\jmath\frac{{\cot\theta}}{2}{t^{2}}\right)}$.
For any function $f\left(t\right)$ we define up-chirp and down-chirp
operations, respectively, by

${\overset{\lower
1.99997pt\hbox{$\smash{\scriptstyle\hpup}$}}{{f}}}\left({t}\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\xi_{\theta}\left(t\right){f}\left(t\right)\quad\mbox{
and }\quad{\overset{\lower
1.99997pt\hbox{$\smash{\scriptstyle\hpdn}$}}{{f}}}\left({t}\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\xi^{*}_{\theta}\left(t\right){f}\left(t\right).$
(3)

When working with low-pass filters, we will be using the fractional
convolution operator introduced by Zayed [32].

######  Definition 2 (Fractional Convolution Operator [32]).

Given functions $f$ and $g$, let $*$ denote the conventional convolution
operator, that is
$\left({f*g}\right)\left(t\right)=\int{f\left(x\right)g\left({t-x}\right)dx}$.
Then, the FrFT convolution operator denoted by $*_{\theta}$, is defined as

$h\left(t\right)=\left(f*_{\theta}g\right)\left(t\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}C_{\theta}\xi^{*}_{\theta}\left(t\right)\left({\overset{\lower
5.0pt\hbox{$\smash{\scriptstyle\hpup}$}}{{f}}}*{\overset{\lower
5.0pt\hbox{$\smash{\scriptstyle\hpup}$}}{{g}}}\right)\left(t\right),$ (4)

where
${C_{\theta}}=\sqrt{\left({1-\jmath\cot\theta}\right)/\left({2\pi}\right)}$.

The $*_{\theta}$ operator (4) admits convolution-multiplication duality. Given
$h=f*_{\theta}g$ in the time domain yields
$\widehat{h}_{\theta}\left(\omega\right)=\xi^{*}_{\theta}\left(\omega\right)\widehat{f}_{\theta}\left(\omega\right)\widehat{g}_{\theta}\left(\omega\right)$
in FrFT domain [32]. In our work, for simplicity of exposition, we will omit
$C_{\theta}$ in our computations.

Signal Model and Measurements. The unknown, continuous-time $K$-sparse signal
is defined as,

$s_{K}(t)=\sum\limits_{k=0}^{K-1}c_{k}\delta\left(t-t_{k}\right),\quad\sum\nolimits_{k}|c_{k}|<\infty.$
(5)

We will follow the measurement setup in [21, 22, 23] where one observes low-
pass filtered version of the spikes in (5) in the FrFT sense. This is
consistent with the conventional sparse sampling and super-resolution problems
[30, 31] in the Fourier domain. Concretely, we will assume that the sampling
kernel is $\phi_{M}\left(t\right)={\overset{\lower
5.0pt\hbox{$\smash{\scriptstyle\hpdn}$}}{{\psi}}}_{M}(t)$ is bandlimited in
the FrFT domain. We will be working with _arbitrary kernels_ of the form,
$\psi_{M}(t)=\sum\nolimits_{m=0}^{M-1}{{p_{m}}\operatorname{sinc}\left({\frac{t}{T}-m}\right)}$
where the weights $\\{p_{m}\\}_{m=0}^{M-1}$ are assumed to be known. The low-
pass filtered measurements simplify to,

$\displaystyle y\left(t\right)$
$\displaystyle=\left(s_{K}*_{\theta}\phi_{M}\right)\left(t\right)$ (6)
$\displaystyle\stackrel{{\scriptstyle(\ref{eq:frc})}}{{=}}\xi^{*}_{\theta}\left(t\right)\sum\limits_{k=0}^{K-1}c_{k}\int_{-\infty}^{\infty}\xi_{\theta}\left(\tau\right){\delta(\tau-
t_{k})\phi_{M}(t-\tau)}d\tau$
$\displaystyle=\xi^{*}_{\theta}\left(t\right)\sum\limits_{k=0}^{K-1}c_{k}\xi_{\theta}\left(t_{k}\right)\sum\limits_{m=0}^{M-1}p_{m}\operatorname{sinc}\left(\frac{t-t_{k}}{T}-m\right).$

Problem Statement. Given $N$ samples, $y\left[n\right]=y\left(nT\right),T>0$,
our goal is to estimate the unknown sparse signal $s_{K}\left(t\right)$.

##  III Recovery Algorithm and Guarantees

Our sparse recovery result is as follows.

######  Theorem 2 (Sampling Criterion).

Let $y\left(t\right)=\left(s_{K}*_{\theta}\phi_{M}\right)\left(t\right)$, as
in (6) and $\phi_{M},\theta$ be known. Suppose we are given $N$ samples of
$y\left[n\right]=y\left(nT\right),T>0$, then, $N\geqslant 2KM$ guarantees
recovery of the unknown signal $s_{K}\left(t\right)$ in (5).

###### Proof.

We provide a proof by construction. Let
$\overline{t}_{k}\stackrel{{\scriptstyle\rm{def}}}{{=}}t_{k}/T$. Given samples
$y\left[n\right]=y\left(nT\right),T>0$, by modulating both sides with
$\xi_{\theta}\left(nT\right)$ we obtain,

$\underbrace{y\left[n\right]{e^{\jmath\frac{{\cot\theta}}{2}{\left(nT\right)^{2}}}}}_{{\overset{\lower
1.39998pt\hbox{$\smash{\scriptstyle\hpup}$}}{{y}}}\left[n\right]}\stackrel{{\scriptstyle(\ref{eq:yt})}}{{=}}\sum\limits_{k=0}^{K-1}c_{k}\xi_{\theta}\left(t_{k}\right)\sum\limits_{m=0}^{M-1}p_{m}\operatorname{sinc}\left(n-\overline{t}_{k}-m\right).$
(7)

Since $\forall n\in\mathbb{Z}$,
$\sin\left({\pi\left({n-x}\right)}\right)={\left({-1}\right)^{n+1}}\sin\left({\pi
x}\right)$ and
$\operatorname{sinc}\left(n-\overline{t}_{k}-m\right)=\frac{(-1)^{n-m+1}\sin\left(\pi\overline{t}_{k}\right)}{\pi\left(n-m-\overline{t}_{k}\right)}$,
we obtain via (7),

${\overset{\lower
1.99997pt\hbox{$\smash{\scriptstyle\hpup}$}}{{y}}}\left[n\right]\stackrel{{\scriptstyle(\ref{eq:yt})}}{{=}}(-1)^{n+1}\sum\limits_{k=0}^{K-1}c_{k}\xi_{\theta}\left(t_{k}\right)\sum\limits_{m=0}^{M-1}p_{m}\frac{(-1)^{-m}\sin(\pi\overline{t}_{k})}{\pi(n-m-\overline{t}_{k})}.$

A re-arrangement of the terms leads to,

$y\left[n\right]\stackrel{{\scriptstyle\rm{def}}}{{=}}\left(\frac{\pi
y\left[n\right]{e^{\jmath\frac{{\cot\theta}}{2}{\left(nT\right)^{2}}}}}{(-1)^{n+1}}\right)\stackrel{{\scriptstyle(\ref{eq:smodel})}}{{=}}\sum\limits_{k=0}^{K-1}\breve{c}^{\theta}_{k}\sum\limits_{m=0}^{M-1}\frac{(-1)^{-m}p_{m}}{\left(n-m-\overline{t}_{k}\right)},$
(8)

where
$\breve{c}^{\theta}_{k}\stackrel{{\scriptstyle\rm{def}}}{{=}}c_{k}\sin\left(\pi\overline{t}_{k}\right){\exp\left(\jmath\frac{{\cot\theta}}{2}{t_{k}^{2}}\right)}$.
Let $\mathsf{Q}\left(z\right)$ be a polynomial of order $KM$ with its roots
$\\{\overline{t}_{k}\\}$, of multiplicity $M$,

$\mathsf{Q}\left(z\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\sum\limits_{k=0}^{\left(KM\right)}q\left[k\right]z^{k}=\prod\limits_{k=0}^{K-1}\prod\limits_{m=0}^{M-1}\left(z-m-\overline{t}_{k}\right).$
(9)

With $p_{m}^{\prime}\stackrel{{\scriptstyle\rm{def}}}{{=}}(-1)^{-m}p_{m}$, we
may now write,

$\underbrace{y\left[n\right]\mathsf{Q}\left(n\right)}_{y^{\mathsf{Q}}\left[n\right]}=\sum\limits_{k=0}^{K-1}\breve{c}^{\theta}_{k}\underbrace{\sum\limits_{m=0}^{M-1}p_{m}^{\prime}\frac{\prod\limits_{k=0}^{K-1}\prod\limits_{m=0}^{M-1}\left(n-m-\overline{t}_{k}\right)}{n-m-\overline{t}_{k}}}_{P_{n,k}},$
(10)

which can be written as
$\mathbf{y}^{\mathsf{Q}}=\mathbf{P}\breve{\mathbf{c}}^{\theta}$ in vector-
matrix notation. Note that,
$y^{\mathsf{Q}}\left[n\right]=\sum\nolimits_{k=0}^{K-1}\breve{c}^{\theta}_{k}P_{n,k}$
is a polynomial of order $\left(KM-1\right)$ because
$\left(n-m-\overline{t}_{k}\right)$ in the denominator of (10) peels away a
monomial. Let $\left(\Delta f\right)[n]=f[n+1]-f[n]$ whose recursive
application results in $\Delta^{L}=\Delta^{L-1}\Delta$. The implication is
that
$\Delta^{KM}\left(\sum\nolimits_{k=0}^{K-1}\breve{c}^{\theta}_{k}P_{n,k}\right)\left[n\right]=0$.
Consequently,

$\left(\Delta^{KM}y^{\mathsf{Q}}\right)\left[n\right]=\sum\limits_{k=0}^{KM}{q\left[k\right]}\underbrace{\Delta^{KM}\left(n^{k}y\left[n\right]\right)}_{D_{n,k}}=0\equiv\mathbf{Dq}=\mathbf{0},$

where $\mathbf{D}$ is the matrix purely depending on data samples $y[n]$,

$\left[\mathbf{D}\right]_{n,k}=\Delta^{KM}\left(n^{k}y\left[n\right]\right)\stackrel{{\scriptstyle(\ref{eq:z_n_definition})}}{{=}}\Delta^{KM}\left(n^{k}\frac{\pi{e^{\jmath\frac{{\cot\theta}}{2}{\left(nT\right)^{2}}}}}{(-1)^{n+1}}y\left[n\right]\right),$

and known FrFT order $\theta$. Clearly, given $\mathbf{y}\in\mathbb{C}^{N}$,
$\Delta^{KM}$ leads to a reduction of $KM$ samples and
$\mathbf{D}\in\mathbb{C}^{\left(N-KM\right)\times\left(KM+1\right)}$. Let
$\ker{\left(\mathbf{M}\right)}=\\{\mathbf{z}\in\mathbb{C}^{N}\ |\
\mathbf{Mz}=\mathbf{0}\\}$ and
$\left[\mathbf{q}\right]_{k}\stackrel{{\scriptstyle(\ref{eq:polyKM})}}{{=}}q\left[k\right]$,
$\mathbf{q}\in\mathbb{R}^{KM+1}$. For $t_{k}\not=t_{\ell},\forall k\not=\ell$,
a non-zero $\mathbf{q}\in\ker{\left(\mathbf{D}\right)}$ can be obtained
provided that $\mathbf{D}$ is rank-deficient or whenever $N-KM\geqslant
KM\Rightarrow N\geqslant 2KM$, which gives the sampling criterion. The
solution to $\mathbf{Dq}=\mathbf{0}\mapsto\mathbf{q}$ allows for construction
of $\mathsf{Q}\left(z\right)$ in (9); its roots then lead to the estimates of
$\\{t_{k}\\}_{k=0}^{K-1}$. It remains to estimate $\\{c_{k}\\}_{k=0}^{K-1}$
which can be obtained by solving linear system of equations in in (7), thus
leading to exact recovery of unknown $s_{K}\left(t\right)$ in (6).
$\blacksquare$

\begin{overpic}[width=216.81pt]{hw_data.png}
\put(95.0,335.0){\parbox{21.68121pt}{\centering\footnotesize$c_{0}$\@add@centering}}
\put(95.0,300.0){\parbox{21.68121pt}{\centering\scriptsize
0.748\@add@centering}}
\put(190.0,330.0){\parbox{21.68121pt}{\centering\footnotesize$c_{1}$\@add@centering}}
\put(190.0,300.0){\parbox{21.68121pt}{\centering\scriptsize
0.891\@add@centering}}
\put(285.0,330.0){\parbox{21.68121pt}{\centering\footnotesize$t_{0}\;\text{(s)}$\@add@centering}}
\put(285.0,300.0){\parbox{21.68121pt}{\centering\scriptsize
0.500\@add@centering}}
\put(380.0,330.0){\parbox{21.68121pt}{\centering\footnotesize$t_{1}\;\text{(s)}$\@add@centering}}
\put(380.0,300.0){\parbox{21.68121pt}{\centering\scriptsize
0.830\@add@centering}}
\put(95.0,250.0){\parbox{21.68121pt}{\centering\footnotesize$\tilde{c}_{0}$\@add@centering}}
\put(95.0,220.0){\parbox{21.68121pt}{\centering\scriptsize
0.733\@add@centering}}
\put(190.0,250.0){\parbox{21.68121pt}{\centering\footnotesize$\tilde{c}_{1}$\@add@centering}}
\put(190.0,220.0){\parbox{21.68121pt}{\centering\scriptsize
0.843\@add@centering}}
\put(285.0,250.0){\parbox{21.68121pt}{\centering\footnotesize$\tilde{t}_{0}\;\text{(s)}$\@add@centering}}
\put(285.0,220.0){\parbox{21.68121pt}{\centering\scriptsize
0.501\@add@centering}}
\put(380.0,250.0){\parbox{21.68121pt}{\centering\footnotesize$\tilde{t}_{1}\;\text{(s)}$\@add@centering}}
\put(380.0,220.0){\parbox{21.68121pt}{\centering\scriptsize
0.834\@add@centering}} \end{overpic} Figure 2: Hardware-based experimental
validation of sparse recovery method. Digital samples with $8$-bit resolution
are marked by $\bullet$.

Note that when $\theta=\pi/2$, Theorem 2 reduces to the conventional Fourier
domain case and $\ast_{\theta},\theta=\pi/2$ reduces to the conventional
notion of low-pass filtering or convolution. Furthermore, when $M=1$, our
results coincide with the Fourier domain recovery result in [31].

Hardware Experiment. Although our algorithm is not designed to tackle noise,
we validate its empirical robusteness via hardware experiment. To this end, we
generate $y(t)$ in (6) with FrFT order $\theta=\pi/4$, $\psi_{M}(t)$ with
$M=1$ and $T=0.062$ s and $s_{K}(t)$ with $t_{k}=\\{0.50,0.83\\}$ and
$c_{k}=\\{0.748,0.891\\}$ using PicoScope 2204A. To capture the complex-valued
signal, we utilize in-phase and quadrature channels, digitizing each channel
using an 8-bit ADC. We use $N=16$ samples annotated by “$\bullet$” (blue and
orange) in Fig. 2. In this setting, the measurements are naturally affected by
thermal and quantization noise, among other factors. Despite this, our
algorithm estimates the sparse signal with mean squared error (MSE) for
$c_{k}\propto 10^{-3}$ and $t_{k}\propto 10^{-6}$, respectively.

##  IV Cramér–Rao Bounds for Sparse Sampling

To benchmark the recovery performance in the presence of noise, we derive CRB
for the sparse parameter estimation problem. To this end, we define,

$z\left[n\right]=y\left[n\right]+\epsilon\left[n\right],\quad\epsilon\sim\mathsf{CN}\left(0,\sigma^{2}\right),$
(11)

where $\mathsf{CN}$ denotes the white complex Gaussian noise distribution. The
CRB is obtained by inverting the Fisher Information Matrix (FIM) parametrized
by $\boldsymbol{\Theta}=\left[t_{0},\cdots,t_{K-1}\ |\
c_{0},\cdots,c_{K-1}\right]^{\top}$. The FIM is given by[33]

$\mathsf{J}\left(\boldsymbol{\Theta}\right)={\mathbf{G}}^{H}\mathbf{R}^{-1}\mathbf{G}=\frac{1}{\sigma^{2}}{\mathbf{G}}^{H}\mathbf{G},$
(12)

where $\mathbf{G}=\frac{\partial
z\left[n\right]}{\partial\boldsymbol{\Theta}}$ is a matrix of partial
derivatives of $z\left[n\right]$ and $\mathbf{R}$ is the covariance matrix of
the $\mathsf{CN}$. CRB matrix is obtained by inverting
$\mathsf{J}\left(\boldsymbol{\Theta}\right)$, yielding
$\mathsf{cov}\left(\theta\right)\geqslant\left(\mathsf{J}\left(\boldsymbol{\Theta}\right)\right)^{-1}$
where,

$\displaystyle\mathsf{J}\left(\boldsymbol{\Theta}\right)=\frac{1}{\sigma^{2}}\begin{bmatrix}{\left(\frac{\partial
z\left[n\right]}{\partial t_{k}}\right)}^{H}\frac{\partial
z\left[n\right]}{\partial t_{m}}&{\left(\frac{\partial
z\left[n\right]}{\partial t_{k}}\right)}^{H}\frac{\partial
z\left[n\right]}{\partial c_{m}}\\\ {\left(\frac{\partial
z\left[n\right]}{\partial c_{k}}\right)}^{H}\frac{\partial
z\left[n\right]}{\partial t_{m}}&{\left(\frac{\partial
z\left[n\right]}{\partial c_{k}}\right)}^{H}\frac{\partial
z\left[n\right]}{\partial c_{m}}\\\ \end{bmatrix},$ (13)

with $k,m\in\\{0,1,\ldots,K-1\\}$. With $\phi\left(t\right)={\overset{\lower
5.0pt\hbox{$\smash{\scriptstyle\hpdn}$}}{{\psi}}}(t)$, we can simplify the
matrix elements in the above as follows:

$\displaystyle\frac{{\partial z}\left[n\right]}{{\partial{c_{k}}}}$
$\displaystyle=\psi\left({n-\frac{{{t_{k}}}}{T}}\right){e^{-\jmath\frac{{\cot(\theta)}}{2}\left({{{\left({nT}\right)}^{2}}-t_{k}^{2}}\right)}}$
(14) $\displaystyle\frac{{\partial z}\left[n\right]}{{\partial{t_{k}}}}$
$\displaystyle={c_{k}}\left({-\frac{{{\psi^{\left(1\right)}}\left({n-{t_{k}}/T}\right)}}{{T\psi\left({n-{t_{k}}/T}\right)}}+\jmath{t_{0}}\cot(\theta)}\right)\left(\frac{{\partial
z\left[n\right]}}{{\partial{c_{k}}}}\right),$

where
${{\psi^{\left(1\right)}}}\left(t\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\partial_{t}\psi\left(t\right)$.
For the general case of arbitrary $K$, the analytical results are intractable.
Hence, we derive the closed-form CRB for a single spike case.

Analytical Bounds: Single Spike Case. When $K=1$, the FIM in (13) simplifies
to,
$\mathsf{J}\left(\boldsymbol{\Theta}\right)\stackrel{{\scriptstyle(\ref{eq:fim-
ders})}}{{=}}\frac{1}{\sigma^{2}}\begin{bmatrix}\mathsf{J}_{1}^{\boldsymbol{\Theta}}&\mathsf{J}_{2}^{\boldsymbol{\Theta}}\\\
\mathsf{J}_{3}^{\boldsymbol{\Theta}}&\mathsf{J}_{4}^{\boldsymbol{\Theta}}\end{bmatrix}$.
Next, with $\boldsymbol{\Theta}=\left[t_{0}\ |\ c_{0}\right]^{\top}$, we
investigate each matrix element separately. Using (14), we simplify the
expressions for $\\{\mathsf{J}_{m}^{\boldsymbol{\Theta}}\\}_{m=1}^{m=4}$ read,

$\displaystyle\left\\{{\begin{array}[]{*{20}{l}}\mathsf{J}_{1}^{\boldsymbol{\Theta}}=\frac{{c_{0}^{2}}}{{{T^{2}}}}\left({{T^{2}}{{\cot}^{2}}\left(\theta\right)t_{0}^{2}\Phi_{\psi,\psi}^{N}\left(\overline{t}_{0}\right)+\Phi_{{\psi^{\left(1\right)}},{\psi^{\left(1\right)}}}^{N}\left(\overline{t}_{0}\right)}\right)\\\
\mathsf{J}_{2}^{\boldsymbol{\Theta}}=-\frac{c_{0}}{T}\left(\Phi_{\psi,\psi^{(1)}}^{N}\left(\overline{t}_{0}\right)+\jmath
T\cot\left(\theta\right){t_{0}}\Phi_{\psi,\psi}^{N}\left(\overline{t}_{0}\right)\right)\\\
\mathsf{J}_{3}^{\boldsymbol{\Theta}}=\frac{c_{0}}{T}\left(\jmath
T\cot\left(\theta\right){t_{0}}\Phi_{\psi,\psi}^{N}\left(\overline{t}_{0}\right)-\Phi_{\psi,\psi^{(1)}}^{N}\left(\overline{t}_{0}\right)\right)\\\
\mathsf{J}_{4}^{\boldsymbol{\Theta}}=\Phi_{\psi,\psi}^{N}\left(\overline{t}_{0}\right)\end{array}}\right.\mbox{
where }\
\Phi_{\psi_{1},\psi_{2}}^{N}\left(t\right)=\sum\limits_{n=0}^{N-1}\psi_{1}\left(n-t\right)\psi_{2}\left(n-t\right).$
(15)

This leads to a method for the numerical computation for the CRB with
arbitrary $\psi$. Note that
$\\{\mathsf{J}_{m}^{\boldsymbol{\Theta}}\\}_{m=1}^{m=4}$ implicitly depend on
$\Phi_{\psi_{1},\psi_{2}}^{\infty}$ in (15). With regards to Theorem 2, we can
compute analytical bounds by setting $\psi=\operatorname{sinc}$ in the
asymptotic case, $N\to\infty$. Hence, we compute the different
parametrizations of $\Phi_{\psi_{1},\psi_{2}}^{\infty}$ to obtain the FIM in
(12). To this end, we define,

$\mathcal{S}_{m}\left(x\right)\stackrel{{\scriptstyle\rm{def}}}{{=}}\sum\limits_{n\in\mathbb{Z}}{\frac{1}{{{{\left({x/T-n}\right)}^{m}}}}}=\sum\limits_{n\in\mathbb{Z}}{\frac{1}{{{{\left({\overline{x}-n}\right)}^{m}}}}},$
(16)

which will be a common factor to all of our computations.

Computation of ${\Phi_{\psi,\psi}^{\infty}\left(\overline{t}_{0}\right)}$. Let
$x=\overline{t}_{0}\in\mathbb{R}$ where $\overline{t}_{0}=t_{0}/T$. Since
$\forall
n\in\mathbb{Z},\sin\left({\pi\left({\overline{x}-n}\right)}\right)={\left({-1}\right)^{n}}\sin\left({\pi\overline{x}}\right)$
is the constituent term in
${\Phi_{\psi,\psi}^{\infty}\left(\overline{x}\right)}$, we can write,

${\left.{\Phi_{\psi,\psi}^{\infty}\left(\overline{x}\right)}\right|_{x={t_{0}}}}=\frac{{{{\sin}^{2}}\left({\pi\overline{x}}\right)}}{{{\pi^{2}}}}\underbrace{\frac{1}{{\sum\nolimits_{n\in\mathbb{Z}}{{{\left({\overline{x}-n}\right)}^{2}}}}}}_{\mathcal{S}_{2}\left(x\right)}=1\to\boxed{\mathsf{J}_{4}^{\boldsymbol{\Theta}}}$

because $\mathcal{S}_{2}\left(x\right)={\pi^{2}}{\csc^{2}}\left({\pi
x}\right)$, which is a well-known identity in complex analysis (see [34], page
189).

Computation of ${\Phi_{\zeta,\zeta}^{\infty}\left(\overline{x}\right)}$ where
$\zeta=\psi^{\left(1\right)}$. We begin by noting,
$\zeta\left({\overline{x}}\right)=\pi\left({\frac{{\cos\left({\pi\overline{x}}\right)}}{{\pi\overline{x}}}-\frac{{\sin\left({\pi\overline{x}}\right)}}{{{\pi^{2}}{{\overline{x}}^{2}}}}}\right).$
The constituent term,

$\displaystyle\underbrace{\left({\zeta}\left({\overline{x}-n}\right)\right)^{2}}_{\zeta=\psi^{\left(1\right)}}=$
$\displaystyle\frac{{{{\sin}^{2}}\left({\pi\left({\overline{x}-n}\right)}\right)}}{{{\pi^{2}}{{\left({\overline{x}-n}\right)}^{4}}}}+\frac{{{{\cos}^{2}}\left({\pi\left({\overline{x}-n}\right)}\right)}}{{{{\left({\overline{x}-n}\right)}^{2}}}}\hfill-\frac{{2\sin\left({\pi\left({\overline{x}-n}\right)}\right)\cos\left({\pi\left({\overline{x}-n}\right)}\right)}}{{\pi{{\left({\overline{x}-n}\right)}^{3}}}},$

can be simplified via basic trigonometric identities since,
${\sin^{2}}\left({\pi\left({\overline{x}-n}\right)}\right)={\sin^{2}}\left({\pi\overline{x}}\right)$,
${\cos^{2}}\left({\pi\left({\overline{x}-n}\right)}\right)={\cos^{2}}\left({\pi\overline{x}}\right)$,
and
$2\sin\left({\pi\left({\overline{x}-n}\right)}\right)\cos\left({\pi\left({\overline{x}-n}\right)}\right)=\sin\left({2\pi\overline{x}}\right)$,
yielding,

$\Phi_{\zeta,\zeta}^{\infty}\left({\overline{x}}\right)=\sum\limits_{n\in\mathbb{Z}}{\left({\frac{{{{\sin}^{2}}\left({\pi\overline{x}}\right)}}{{{\pi^{2}}{{\left({\overline{x}-n}\right)}^{4}}}}+\frac{{{{\cos}^{2}}\left({\pi\overline{x}}\right)}}{{{{\left({\overline{x}-n}\right)}^{2}}}}-\frac{{\sin\left({2\pi\overline{x}}\right)}}{{\pi{{\left({\overline{x}-n}\right)}^{3}}}}}\right)}.$
(17)

In the above, note that the terms corresponding to
$\\{\sin^{2}\left(\pi\overline{x}\right),\cos^{2}\left(\pi\overline{x}\right),\sin\left(2\pi\overline{x}\right)\\}$
are linked with the sums
$\\{\mathcal{S}_{4}\left(x\right),\mathcal{S}_{2}\left(x\right),\mathcal{S}_{3}\left(x\right)\\}$
in (16), respectively. Known that,
$\mathcal{S}_{2}\left(x\right)={\pi^{2}}{\csc^{2}}\left({\pi\overline{x}}\right)$
(see [34], page 189). To compute $\mathcal{S}_{3}\left(x\right)$ and
$\mathcal{S}_{4}\left(x\right)$, we proceed by differentiating
$\mathcal{S}_{2}\left(x\right)$,

$\displaystyle\underbrace{{\partial_{x}\mathcal{S}_{2}\left(x\right)}}_{\mathcal{S}_{2}^{\left(1\right)}\left(x\right)}=-\frac{2}{T}\underbrace{\sum\limits_{n\in\mathbb{Z}}{\frac{1}{{{{\left({\overline{x}-n}\right)}^{3}}}}}}_{\mathcal{S}_{3}\left(x\right)}=\underbrace{-\frac{{2{\pi^{3}}}}{T}\cot\left({\pi\overline{x}}\right){{\csc}^{2}}\left({\pi\overline{x}}\right)}_{{\partial_{x}}\left({{\pi^{2}}{{\csc}^{2}}\left({\pi\overline{x}}\right)}\right)}$
$\displaystyle\Rightarrow\mathcal{S}_{3}\left(x\right)={\pi^{3}}\cot\left({\pi\overline{x}}\right){\csc^{2}}\left({\pi\overline{x}}\right).$

Similarly,

$\displaystyle\underbrace{{\partial^{2}_{x}\mathcal{S}_{2}\left(x\right)}}_{\mathcal{S}_{2}^{\left(2\right)}\left(x\right)}=\left(\frac{6}{T^{2}}\right){\mathcal{S}_{4}\left(x\right)}=\underbrace{\frac{{2{\pi^{4}}}}{{{T^{2}}}}\left({\cos\left({2\pi\overline{x}}\right)+2}\right){\csc^{4}}\left({\pi\overline{x}}\right)}_{{\partial^{2}_{x}}\left({{\pi^{2}}{{\csc}^{2}}\left({\pi\overline{x}}\right)}\right)}$
$\displaystyle\Rightarrow\mathcal{S}_{4}\left(x\right)=\tfrac{\pi^{4}}{3}\left({\cos\left({2\pi\overline{x}}\right)+2}\right){\csc^{4}}\left({\pi\overline{x}}\right).$

Given,
$\Phi_{\zeta,\zeta}^{\infty}\left({\overline{x}}\right)=\left({\sin}\left({\pi\overline{x}}\right)/\pi\right)^{2}\mathcal{S}_{4}\left(x\right)+{\cos^{2}}\left({\pi\overline{x}}\right)\mathcal{S}_{2}\left(x\right)-\left(\sin\left({2\pi\overline{x}}\right)/\pi\right)\mathcal{S}_{3}\left(x\right)$
(see (17)), substituting for $\mathcal{S}_{m}\left(x\right),m=\\{2,3,4\\}$ in
(17) with $\overline{x}=x/T$, we obtain,

$\boxed{\Phi_{{\psi^{\left(1\right)}},{\psi^{\left(1\right)}}}^{\infty}\left({\overline{x}}\right)}=\sum\limits_{n\in\mathbb{Z}}\left({{\psi^{\left(1\right)}}\left({\frac{x}{T}-n}\right)}\right)^{2}=\frac{{{\pi^{2}}}}{{3}}.$
(18)

Computation of ${\Phi_{\zeta,\psi}^{\infty}\left(\overline{x}\right)}$ where
$\zeta=\psi^{\left(1\right)}$. This term turns out to be zero. To see this, we
expand the constituent term in
$\Phi_{{\psi^{\left(1\right)}},\psi}^{\infty}\left({\overline{x}}\right)$,
that is,
$\psi\left({\overline{x}-n}\right){\psi^{\left(1\right)}}\left({\overline{x}-n}\right)=\nicefrac{{{\sin(\pi(\overline{x}-n))\cos(\pi(\overline{x}-n))}}}{{{\pi{{(\overline{x}-n)}^{2}}}}}-\nicefrac{{{{{\sin}^{2}}(\pi(\overline{x}-n))}}}{{{{\pi^{2}}{{(\overline{x}-n)}^{3}}}}}$.
Using
$\sin(\pi(\overline{x}-n))\cos(\pi(\overline{x}-n))=\sin(2\pi\overline{x})/2$,
we conclude,

$\displaystyle\boxed{\Phi_{{\psi^{\left(1\right)}},\psi}^{\infty}\left({\overline{x}}\right)}$
$\displaystyle=\sum\limits_{n\in\mathbb{Z}}{\frac{{\sin(2\pi\overline{x})}}{{2\pi{{(n-\overline{x})}^{2}}}}}-\frac{{{{\sin}^{2}}(\pi\overline{x})}}{{{\pi^{2}}{{(n-\overline{x})}^{3}}}}=\underbrace{\frac{{\sin(2\pi\overline{x})}}{{2\pi}}\mathcal{S}_{2}\left(x\right)}_{\pi\cot\left(\pi\overline{x}\right)}-\underbrace{\frac{{{{\sin}^{2}}(\pi\overline{x})}}{{{\pi^{2}}}}\mathcal{S}_{3}\left(x\right)}_{\pi\cot\left(\pi\overline{x}\right)}=0.$

In summary, our computations show,

$\begin{array}[]{*{20}{c}}{\Phi_{\psi,\psi}^{\infty}\left(\overline{t}_{0}\right)}=1&{\Phi_{\psi^{(1)},\psi^{(1)}}^{\infty}\left(\overline{t}_{0}\right)}=\frac{\pi^{2}}{3}&\Phi_{{\psi^{\left(1\right)}},\psi}^{\infty}\left({\overline{t}_{0}}\right)=0.\end{array}$

CRB Computation. Since we have all the values of
$\Phi_{\psi_{1},\psi_{2}}^{\infty}$ in (15) that define
$\\{\mathsf{J}_{m}^{\boldsymbol{\Theta}}\\}_{m=1}^{m=4}$, we can now obtain
the FIM by back substituting the obtained values of
$\Phi_{\psi_{1},\psi_{2}}^{\infty}$ yielding the FIM elements, (i)
$\mathsf{J}_{1}^{\boldsymbol{\Theta}}=c_{0}^{2}\left({{{\cot}^{2}}\left(\theta\right)t_{0}^{2}+\frac{{{\pi^{2}}}}{{3{T^{2}}}}}\right)$,
(ii) $\mathsf{J}_{2}^{\boldsymbol{\Theta}}=-\jmath
c_{0}\cot\left(\theta\right){t_{0}}$, (iii)
$\mathsf{J}_{3}^{\boldsymbol{\Theta}}=\jmath
c_{0}\cot\left(\theta\right){t_{0}}$, and, (iv)
$\mathsf{J}_{4}^{\boldsymbol{\Theta}}=1$. The corresponding FIM takes the form
of,

$\displaystyle\mathsf{J}\left(\boldsymbol{\Theta}\right)=\frac{1}{\sigma^{2}}\begin{bmatrix}c_{0}^{2}\left({{{\cot}^{2}}\left(\theta\right)t_{0}^{2}+\tfrac{{{\pi^{2}}}}{{3{T^{2}}}}}\right)&-\jmath
c_{0}\cot\left(\theta\right){t_{0}}\\\ \jmath
c_{0}\cot\left(\theta\right){t_{0}}&1\\\ \end{bmatrix},$

with
${\operatorname{det}\left(\mathsf{J}\left(\boldsymbol{\Theta}\right)\right)}=\left({{\pi^{2}}c_{0}^{2}}\right)/\left({3{T^{2}}}\right)$.
The CRB is given by,

$\mathsf{cov}\left(\boldsymbol{\Theta}\right)=\frac{\sigma^{2}}{\operatorname{det}\left(\mathsf{J}\left(\boldsymbol{\Theta}\right)\right)}\begin{bmatrix}\mathsf{J}_{4}^{\boldsymbol{\Theta}}&-\mathsf{J}_{2}^{\boldsymbol{\Theta}}\\\
-\mathsf{J}_{3}^{\boldsymbol{\Theta}}&\mathsf{J}_{1}^{\boldsymbol{\Theta}}\end{bmatrix}$
(19)

and using $\mathsf{cov}\left(\boldsymbol{\Theta}\right)$ and defining
$\mathsf{PSNR}=c_{0}^{2}/\sigma^{2}$, we obtain,

$\displaystyle\mathsf{cov}\left(\boldsymbol{\Theta}\right)$
$\displaystyle=\left(\frac{3{T^{2}}}{\pi^{2}\mathsf{PSNR}}\right)\begin{bmatrix}1&\jmath
c_{0}\cot\left(\theta\right){t_{0}}\\\ -\jmath
c_{0}\cot\left(\theta\right){t_{0}}&c_{0}^{2}\left({{{\cot}^{2}}\left(\theta\right)t_{0}^{2}+\tfrac{{{\pi^{2}}}}{{3{T^{2}}}}}\right)\end{bmatrix},$

giving,

$\displaystyle\mathsf{var}\left(t_{0}\right)\geqslant\frac{3T^{2}}{\pi^{2}\mathsf{PSNR}}$
$\displaystyle\mathsf{var}\left(c_{0}\right)\geqslant\frac{3c_{0}^{2}T^{2}}{\pi^{2}\mathsf{PSNR}}\left(\left(t_{0}\cot\left(\theta\right)\right)^{2}+\frac{{{\pi^{2}}}}{{3{T^{2}}}}\right).$

##  V Conclusion

This paper introduces a new time-domain method for recovery of sparse signals
from low-pass filtered measurements in the Fractional Fourier Transform (FrFT)
domain. The advantage of this method is that unlike previous methods that work
in the FrFT domain and are susceptible to spectral leakage, the new method
overcomes restrictions. We also derive Cramér–Rao Bounds (CRB) for the sparse
estimation problem that was missing in previous literature. Future work in
this area includes recovery with different classes of sampling kernels,
particularly non-bandlimited ones and their performance evaluation in terms of
Cramér–Rao Bounds.

## References

  * [1] N. Wiener, “Hermitian polynomials and Fourier analysis,” _MIT J. Math. Phys._ , vol. 18, pp. 70–73, 1929.
  * [2] E. U. Condon, “Immersion of the Fourier transform in a continuous group of functional transformations,” _Proc. Natl. Acad. Sci. U.S.A._ , vol. 23, no. 3, p. 158, 1937.
  * [3] H. M. Ozaktas, Z. Zalevsky, and M. A. Kutay, _The Fractional Fourier Transform_ , ser. Series in Pure and Applied Optics. John Wiley & Sons, 2001.
  * [4] M. Martone, “A multicarrier system based on the fractional Fourier transform for time-frequency-selective channels,” _IEEE Trans. Commun._ , vol. 49, no. 6, pp. 1011–1020, 2001.
  * [5] S. I. Yetik and A. Nehorai, “Beamforming using the fractional Fourier transform,” _IEEE Trans. Signal Process._ , vol. 51, no. 6, pp. 1663–1668, 2003.
  * [6] V. Namias, “The fractional order Fourier transform and its application to quantum mechanics,” _J. Inst. Math. Appl._ , vol. 25, p. 241, 1961.
  * [7] A. I. Zayed, “A class of fractional integral transform: A generalization of the fractional Fourier transform,” _IEEE Trans. Signal Process._ , vol. 50, pp. 619–627, Mar. 2002.
  * [8] ——, “On the relationship between the Fourier and fractional Fourier transforms,” _IEEE Signal Process. Lett._ , vol. 3, no. 12, pp. 310–311, Dec 1996.
  * [9] L. B. Almeida, “The fractional Fourier transform and time-frequency representations,” _IEEE Trans. Signal Process._ , vol. 42, no. 11, pp. 3084–3091, Nov 1994.
  * [10] S.-C. Pei and J.-J. Ding, “Relations between Gabor transforms and fractional Fourier transforms and their applications for signal processing,” _IEEE Trans. Signal Process._ , vol. 55, no. 10, pp. 4839–4850, oct 2007\.
  * [11] A. S. Amein and J. J. Soraghan, “The fractional Fourier transform and its application to high resolution SAR imaging,” in _Proc. IEEE Int. Geoscience and Remote Sensing Symp._ , 2007, pp. 5174–5177.
  * [12] X. G. Xia, “On bandlimited signals with fractional Fourier transform,” _IEEE Signal Process. Lett._ , vol. 3, no. 3, pp. 72–74, Mar 1996.
  * [13] R. Torres, P. F. Pellat, and Y. Torres, “Sampling theorem for fractional bandlimited signals: A self-contained proof. application to digital holography,” _IEEE Signal Process. Lett._ , vol. 13, no. 11, pp. 676–679, Nov. 2006.
  * [14] R. Tao, B. Deng, Z.-Q. Wei, and Y. Wang, “Sampling and sampling rate conversion of band limited signals in the fractional Fourier transform domain,” _IEEE Trans. Signal Process._ , vol. 56, no. 1, pp. 158–171, Jan. 2008.
  * [15] C. Candan and H. M. Ozaktas, “Sampling and series expansion theorems for fractional Fourier and other transforms,” _Signal Process._ , vol. 83, no. 11, pp. 2455–2457, 2003.
  * [16] A. G. Garcia, “Orthogonal sampling formulas: A unified approach,” _SIAM Rev._ , vol. 42, no. 3, pp. 499–512, jan 2000.
  * [17] D. Wei, Q. Ran, and Y. Li, “Generalized sampling expansion for bandlimited signals associated with the fractional Fourier transform,” _IEEE Signal Process. Lett._ , vol. 17, no. 6, pp. 595–598, 2010.
  * [18] T. Erseghe, P. Kraniauskas, and G. Carioraro, “Unified fractional Fourier transform and sampling theorem,” _IEEE Trans. Signal Process._ , vol. 47, no. 12, pp. 3419–3423, Dec 1999.
  * [19] A. Bhandari and A. I. Zayed, “Shift-invariant and sampling spaces associated with the fractional Fourier transform domain,” _IEEE Trans. Signal Process._ , vol. 60, no. 4, pp. 1627–1637, Apr. 2012.
  * [20] ——, “Shift-invariant and sampling spaces associated with the special affine Fourier transform,” _Appl. Comput. Harmon. Anal._ , vol. 47, no. 1, pp. 30–52, Jul. 2019.
  * [21] A. Bhandari and P. Marziliano, “Sampling and reconstruction of sparse signals in fractional Fourier domain,” _IEEE Signal Process. Lett._ , vol. 17, no. 3, pp. 221–224, 2010.
  * [22] A. Bhandari, Y. C. Eldar, and R. Raskar, “Super-resolution in phase space,” in _IEEE Intl. Conf. on Acoustics, Speech and Sig. Proc._ IEEE, Apr. 2015, pp. 4155–4159.
  * [23] A. Bhandari and Y. C. Eldar, “Sampling and super resolution of sparse signals beyond the Fourier domain,” _IEEE Trans. Signal Process._ , vol. 67, no. 6, pp. 1508–1521, Mar. 2019.
  * [24] P. Zhang and A. Bhandari, “Unlimited sampling in phase space,” in _Proc. IEEE Intl. Conf. on Acoustics, Speech and Signal Proc. (ICASSP)_. IEEE, Jun. 2023.
  * [25] A. Bhandari, F. Krahmer, and R. Raskar, “On unlimited sampling and reconstruction,” _IEEE Trans. Signal Process._ , vol. 69, pp. 3827–3839, Dec. 2020.
  * [26] ——, “On unlimited sampling,” in _International Conference on Sampling Theory and Applications (SampTA)_ , Jul. 2017.
  * [27] S.-C. Pei, M.-H. Yeh, and C.-C. Tseng, “Discrete fractional Fourier transform based on orthogonal projections,” _IEEE Trans. Signal Process._ , vol. 47, no. 5, pp. 1335–1348, May 1999.
  * [28] C. Candan, M. Kutay, and H. Ozaktas, “The discrete fractional Fourier transform,” _IEEE Trans. Signal Process._ , vol. 48, no. 5, pp. 1329–1337, May 2000.
  * [29] S. Liu, T. Shan, R. Tao, Y. D. Zhang, G. Zhang, F. Zhang, and Y. Wang, “Sparse discrete fractional Fourier transform and its applications,” _IEEE Trans. Signal Process._ , vol. 62, no. 24, pp. 6582–6595, Dec. 2014.
  * [30] R. J. P. De Figueiredo and C.-L. Hu, “Waveform feature extraction based on Tauberian approximation,” _IEEE Trans. Pattern Anal. Mach. Intell._ , no. 2, pp. 105–116, Mar. 1982.
  * [31] M. Vetterli, P. Marziliano, and T. Blu, “Sampling signals with finite rate of innovation,” _IEEE Trans. Signal Process._ , vol. 50, no. 6, pp. 1417–1428, Jun. 2002.
  * [32] A. I. Zayed, “A convolution and product theorem for the fractional Fourier transform,” _IEEE Trans. Signal Process._ , vol. 5, no. 4, pp. 101–103, Apr 1998.
  * [33] B. Porat and B. Friedlander, “Computation of the exact information matrix of gaussian time series with stationary random components,” _IEEE Transactions on Acoustics, Speech, and Signal Processing_ , vol. 34, no. 1, pp. 118–130, 1986.
  * [34] L. Ahlfors, _Complex Analysis: An Introduction to The Theory of Analytic Functions of One Complex Variable_. McGraw-Hill Education, 1979.