# Decoupling of the structure functions in momentum space based on the Laplace
transformation

G.R.Boroun<EMAIL_ADDRESS>Department of Physics, Razi University,
Kermanshah 67149, Iran  Phuoc Ha<EMAIL_ADDRESS>Department of Physics,
Astronomy and Geosciences, Towson University, Towson, MD 21252

###### Abstract

Using Laplace transform techniques, we describe the determination of the
longitudinal structure function $F_{L}(x,Q^{2})$, at the leading-order
approximation in momentum space, from the structure function $F_{2}(x,Q^{2})$
and its derivative with respect to ${\ln}Q^{2}$ in a kinematical region of low
values of the Bjorken variable $x$. Since the $x$ dependence of
$F_{2}(x,Q^{2})$ and its evolution with $Q^{2}$ are determined much better by
the data than $F_{L}(x,Q^{2})$, this method provides both a direct check on
$F_{L}(x,Q^{2})$ where measured, and a way of extending $F_{L}(x,Q^{2})$ into
regions of $x$ and $Q^{2}$ where there are currently no data. In our
calculations, we ultilize the Block-Durand-Ha parametrization for the
structure function $F_{2}(x,Q^{2})$ [M. M. Block, L. Durand and P. Ha,
Phys.Rev.D 89, 094027 (2014)]. We find that the Laplace transform method in
momentum space provides correct behaviors of the extracted longitudinal
structure function $F_{L}(x,Q^{2})$ and that our obtained results are in line
with data from the H1 Collaboration and other results for $F_{L}(x,Q^{2})$
obtained using Mellin transform method.

****

##  I Introduction

Recently, evolution of the longitudinal and transversal structure functions in
momentum space has been considered in [1]. Structure functions measurable in
deep inelastic scattering (DIS) are formulated in the momentum-space
Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution equations [2, 3,
4]. Scheme independent evolution equations for the structure functions
$F_{i}(x,Q^{2})$ proposed some time ago in [5, 6], as the physical
observables, read

$\displaystyle\frac{{\partial}F_{i}(x,Q^{2})}{{\partial}{\ln}Q^{2}}=\sum_{j}P_{ij}\otimes
F_{j}(x,Q^{2})$ (1)

where anomalous dimensions, $P_{ij}$, are computable in perturbative QCD
(pQCD). Determination of the longitudinal structure function in the nucleon
from the proton structure function, based on a form of the deep inelastic
lepton-hadron scattering structure function proposed by Block-Durand-Ha (BDH)
in [7], is considered in [8, 9, 10, 11]. Parametrization of the proton
structure function proposed in [7] describes the available experimental data
on the reduced cross sections at low $x$ and provides a behavior of the
hadron-hadron cross sections ${\sim}{\ln^{2}}s$ at large $s$ in a full
accordance with the Froissart predictions [12] ($s$ is the Mandelstam variable
denoting the square of the total invariant energy of the process).

Deep inelastic scattering (DIS) is characterized by structure functions
$F_{k}$ that depend on kinematic variables Bjorken $x$ and momentum transfer
$Q$ by the following form

$\displaystyle
F_{k}(x,Q^{2})=<e^{2}>\sum_{a=s,g}\bigg{[}C_{k,a}(x,\alpha_{s}){\otimes}xf_{a}(x,Q^{2})\bigg{]},~{}k=2,L~{}$
(2)

where $C_{k,a}(x,\alpha_{s})$ are the known Wilson coefficient functions in
the order of the perturbation theory, $\alpha_{s}$ is the strong coupling, and
$<e^{2}>$ is the average of the charge $e^{2}$ for the active quark flavors,
$<e^{2}>=n_{f}^{-1}\sum_{i=1}^{n_{f}}e_{i}^{2}$ with $n_{f}$ as the number of
considered flavors. The symbol $\otimes$ denotes convolution according to the
usual prescription and $xf_{a=2,g}(x,Q^{2})$ are the singlet-quark and gluon
densities respectively (the non-singlet quark distributions at small $x$
become negligibly small in comparison with the singlet distributions). The
DGLAP equations, which describe how the parton distribution functions (PDFs)
vary as the energy scale of the scattering process changes, are important for
understanding a wide range of high-energy processes, including DIS, hadron
collisions, and deep-inelastic scattering of heavy ions at future colliders
(Large Hadron electron Collider (LHeC) [13] and Electron-Ion Collider (EIC)
[14, 15]).

Numerical and analytical methods (which extract the PDFs from the experimental
data) to solve the DGLAP evolution equations have been extensively studied in
the literature [16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30,
31, 32]. The solutions to these equations provide a theoretical prediction for
the PDFs used in the interpretation and description of the Hadron-Electron
Ring Accelerator (HERA) data on the total and diffractive cross-sections in
deep inelastic electron-proton scattering. They serve as a mean to test our
understanding of QCD and extract information about the structure of the
proton.

In this paper, we extend the method using a Laplace-transform technique and
obtain an analytical method for the solution of the momentum-space of the
DGLAP equations for $F_{L}(x,Q^{2})$ in terms of $F_{2}(x,Q^{2})$ and known
derivative $dF_{2}(x,Q^{2})/d{\ln}Q^{2}$ in the kinematical region of low
values of the Bjorken variable $x$. The parameterization of the structure
function $F_{2}(x,Q^{2})$ in [7] is obtained from a combined fit to HERA data
in a wide range of the kinematical variables $x$ and $Q^{2}$ by the following
explicit expression as

$\displaystyle F^{\gamma
p}_{2}(x,Q^{2})=D(Q^{2})(1-x)^{n}\bigg{[}C(Q^{2})+A(Q^{2})\ln(\frac{1}{x}\frac{Q^{2}}{Q^{2}+\mu^{2}})+B(Q^{2})\ln^{2}(\frac{1}{x}\frac{Q^{2}}{Q^{2}+\mu^{2}})\bigg{]},$
(3)

where

$\displaystyle A(Q^{2})$ $\displaystyle=$ $\displaystyle
a_{0}+a_{1}{\ln}(1+\frac{Q^{2}}{\mu^{2}})+a_{2}{\ln}^{2}(1+\frac{Q^{2}}{\mu^{2}}),$
$\displaystyle B(Q^{2})$ $\displaystyle=$ $\displaystyle
b_{0}+b_{1}{\ln}(1+\frac{Q^{2}}{\mu^{2}})+b_{2}{\ln}^{2}(1+\frac{Q^{2}}{\mu^{2}}),$
$\displaystyle C(Q^{2})$ $\displaystyle=$ $\displaystyle
c_{0}+c_{1}{\ln}(1+\frac{Q^{2}}{\mu^{2}}),$ $\displaystyle D(Q^{2})$
$\displaystyle=$ $\displaystyle\frac{Q^{2}(Q^{2}+\lambda
M^{2})}{(Q^{2}+M^{2})^{2}},$ (4)

where $M$ is the effective mass and $\mu^{2}$ is a scale factor defined by the
Block-Halzen fit to the real photon-proton cross section [33] in Table I. In
the following, we apply this parametrization function to test the consistency
of the longitudinal structure function owing to the momentum-space of the
DGLAP equations with HERA data on deep inelastic electron-proton scattering.

The paper is organized as follows: in section II, we present the basics of the
momentum-space of the DGLAP equations. Section III summarizes the Laplace
transform method for obtaining an analytical solution for the longitudinal
structure function. In section IV, the numerical results are obtained and
compared with the available H1 Collaboration data [34, 35, 36] and the Large
Hadron electron Collider (LHeC) [13] simulated errors. Conclusions are given
in Sec. V. Some detailed calculations are relegated to Appendix.

##  II Momentum Space

The DIS structure function $F_{2}$ at low $x$ is defined [1] into the singlet
and gluon distribution functions by the following form

$\displaystyle F_{2}(x,Q^{2})$ $\displaystyle=$
$\displaystyle<e^{2}>\bigg{\\{}C^{(0)}_{2,s}+\frac{\alpha_{s}(\mu_{r}^{2})}{2\pi}\bigg{[}C^{(1)}_{2,s}-\ln{\bigg{(}}\frac{\mu_{r}^{2}}{Q^{2}}{\bigg{)}}C^{(0)}_{2,s}{\otimes}P_{qq}\bigg{]}\bigg{\\}}{\otimes}x\Sigma(x,\mu_{r}^{2})$
(5)
$\displaystyle+2\sum_{i=1}^{n_{f}}e_{i}^{2}\frac{\alpha_{s}(\mu_{r}^{2})}{2\pi}\bigg{[}C^{(1)}_{2,g}-\ln{\bigg{(}}\frac{\mu_{r}^{2}}{Q^{2}}{\bigg{)}}C^{(0)}_{2,g}{\otimes}P_{qg}\bigg{]}{\otimes}xg(x,\mu_{r}^{2}),$

where $xg(x,\mu_{r}^{2})$ and
$x\Sigma(x,\mu_{r}^{2}){\equiv}xf_{s}(x,\mu_{r}^{2})=\sum_{q}[xq(x,\mu_{r}^{2})+x\overline{q}(x,\mu_{r}^{2})]$
are the gluon and singlet distribution function at the renormalization scale
$\mu_{r}^{2}$, respectively. The authors in [1] inverted the leading non-zero
order part of Eq.(5). As a result, the singlet and gluon densities in terms of
the structure functions read

$\displaystyle\Sigma(x,\mu_{r}^{2})$ $\displaystyle=$
$\displaystyle\frac{1}{<e^{2}>}\frac{F_{2}(x,Q^{2})}{x},$ $\displaystyle
g(x,\mu_{r}^{2})$ $\displaystyle=$
$\displaystyle\frac{1}{\sum_{i=1}^{n_{f}}e_{i}^{2}}\bigg{(}\frac{C_{F}}{4T_{R}}\delta(1-x){\otimes}x\frac{d}{dx}\frac{F_{2}(x,Q^{2})}{x}-\frac{C_{F}}{2T_{R}}\delta(1-x){\otimes}\frac{F_{2}(x,Q^{2})}{x}+\frac{1}{8T_{R}}\frac{2\pi}{\alpha_{s}(\mu_{r}^{2})}\delta(1-x){\otimes}x^{2}\frac{d^{2}}{dx^{2}}\frac{F_{L}(x,Q^{2})}{x}$
(6)
$\displaystyle-\frac{1}{4T_{R}}\frac{2\pi}{\alpha_{s}(\mu_{r}^{2})}\delta(1-x){\otimes}x\frac{d}{dx}\frac{F_{L}(x,Q^{2})}{x}+\frac{1}{4T_{R}}\frac{2\pi}{\alpha_{s}(\mu_{r}^{2})}\delta(1-x){\otimes}\frac{F_{L}(x,Q^{2})}{x}\bigg{)},$

with the color factors $C_{A}=3$, $T_{R}=1/2$ and $C_{F}=4/3$ associated with
the color group SU(3). Using the leading-order (LO) renormalization group
equation

$\displaystyle\mu_{r}^{2}\frac{d\alpha_{s}(\mu_{r}^{2})}{d\mu_{r}^{2}}=-(\frac{11C_{A}-4T_{R}n_{f}}{6})\alpha^{2}_{s}(\mu_{r}^{2}),$
(7)

the evolution equation of the structure function $F_{2}(x,Q^{2})$ at the first
non-zero order in $\alpha_{s}$ reads

$\displaystyle\frac{dF_{2}(x,Q^{2})}{d{\ln}Q^{2}}$ $\displaystyle=$
$\displaystyle\frac{\alpha_{s}(Q^{2})}{2\pi}x\bigg{\\{}\frac{1}{4}\bigg{(}\frac{2}{x}-\frac{d}{dx}\bigg{)}\frac{2\pi}{\alpha_{s}(Q^{2})}F_{L}(x,Q^{2})+\frac{1}{2}\int_{x}^{1}\frac{dz}{z^{2}}\frac{2\pi}{\alpha_{s}(Q^{2})}F_{L}(z,Q^{2})$
(8)
$\displaystyle+C_{F}\bigg{[}\frac{1}{x}F_{2}(x,Q^{2})-\int_{x}^{1}\frac{dz}{z^{2}}\bigg{(}1-\frac{x}{z}\bigg{)}F_{2}(z,Q^{2})+\frac{1}{x}\int_{x}^{1}dz\frac{1+z^{2}}{(1-z)_{+}}F_{2}(\frac{x}{z},Q^{2})\bigg{]}\bigg{\\}},$

where the renormalization scale is equal to the momentum transfer,
$\mu_{r}^{2}=Q^{2}$, and the plus function is defined as

$\displaystyle\int_{x}^{1}dz\frac{f(z)}{(1-z)_{+}}=\int_{x}^{1}dz\frac{f(z)-f(1)}{1-z}+f(1)\ln(1-x).$
(9)

By writing

$\frac{1+z^{2}}{(1-z)_{+}}=-(1+z)+\frac{2}{(1-z)_{+}},$ (10)

we can rewrite the last integral in Eq. (8) as

$\displaystyle\frac{1}{x}\int_{x}^{1}dz\frac{1+z^{2}}{(1-z)_{+}}F_{2}(\frac{x}{z},Q^{2})$
$\displaystyle=$
$\displaystyle-\frac{1}{x}\int_{x}^{1}dz(1+z)F_{2}(\frac{x}{z},Q^{2})+\frac{1}{x}\int_{x}^{1}dz\frac{2}{(1-z)_{+}}F_{2}(\frac{x}{z},Q^{2}))$
(11) $\displaystyle=$
$\displaystyle-\int_{x}^{1}\frac{dz}{z^{2}}\bigg{(}1+\frac{x}{z}\bigg{)}F_{2}(z,Q^{2})+\frac{1}{x}\int_{x}^{1}dz\frac{2}{(1-z)_{+}}F_{2}(\frac{x}{z},Q^{2})).$

Substituting Eq. (11) into Eq. (8), we get

$\displaystyle\frac{dF_{2}(x,Q^{2})}{d{\ln}Q^{2}}$ $\displaystyle=$
$\displaystyle\frac{\alpha_{s}(Q^{2})}{2\pi}x\bigg{\\{}\frac{1}{4}\bigg{(}\frac{2}{x}-\frac{d}{dx}\bigg{)}\frac{2\pi}{\alpha_{s}(Q^{2})}F_{L}(x,Q^{2})+\frac{1}{2}\int_{x}^{1}\frac{dz}{z^{2}}\frac{2\pi}{\alpha_{s}(Q^{2})}F_{L}(z,Q^{2})$
(12)
$\displaystyle+C_{F}\bigg{[}\frac{1}{x}F_{2}(x,Q^{2})-2\int_{x}^{1}\frac{dz}{z^{2}}F_{2}(z,Q^{2})+\frac{2}{x}\int_{x}^{1}dz\frac{1}{(1-z)_{+}}F_{2}(\frac{x}{z},Q^{2})\bigg{]}\bigg{\\}},$

##  III  Laplace Transformation

In the following, we use the method developed in detail in [37, 38, 39, 40] to
obtain the longitudinal structure function into the proton structure function
and its derivative using a Laplace-transform method. We now rewrite the
momentum-space DGLAP evolution equation for the longitudinal structure
function (i.e., Eq.(12)) in terms of the variables $\upsilon={\ln}(1/x)$ and
$Q^{2}$ instead of $x$ and $Q^{2}$. Using the notation
$\widehat{F_{i}}(\upsilon,Q^{2}){\equiv}F_{i}(e^{-\upsilon},Q^{2})$ for
structure functions, explicitly, from Eq.(12), we find

$\displaystyle\frac{d\widehat{F}_{2}(\upsilon,Q^{2})}{d{\ln}Q^{2}}$
$\displaystyle=$
$\displaystyle\bigg{(}\frac{1}{2}+\frac{1}{4}\frac{d}{d\upsilon}\bigg{)}\widehat{F}_{L}(\upsilon,Q^{2})+\frac{1}{2}\int_{0}^{\upsilon}{dw}~{}e^{-(\upsilon-w)}\widehat{F}_{L}(w,Q^{2})+C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\widehat{F}_{2}(\upsilon,Q^{2})$
(13)
$\displaystyle-2C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\int_{0}^{\upsilon}{dw}~{}e^{-(\upsilon-w)}\widehat{F}_{2}(w,Q^{2})+2C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\int_{0}^{v}dw\ln\left(1-e^{-(v-w)}\right)\frac{\partial\widehat{F}_{2}(w,Q^{2})}{\partial
w}.$

Introducing the notation that the Laplace transform of structure functions
$\widehat{F}_{i}(\upsilon,Q^{2})$ are given by ${f}_{i}(s,Q^{2})$ as
${f}_{i}(s,Q^{2}){\equiv}{\mathcal{L}}[\widehat{F}_{i}(\upsilon,Q^{2});s]$ and
using the fact that the Laplace transform of a convolution factors is simply
the ordinary product of the Laplace transform of the factors by the following
form

$\displaystyle{\mathcal{L}}\bigg{[}\int_{0}^{\upsilon}\widehat{F}_{i}(w,Q^{2})\widehat{H}(\upsilon-w)dw;s\bigg{]}={f}_{i}(s,Q^{2}){\times}h(s),$
(14)

where $h(s){\equiv}{\mathcal{L}}[e^{-\upsilon}\widehat{H}(\upsilon)]$, we find
that the Laplace transform of Eq. (13) is given by

$\displaystyle\frac{d{F}_{2}(s,Q^{2})}{d{\ln}Q^{2}}$ $\displaystyle=$
$\displaystyle\bigg{(}\frac{1}{2}+\frac{1}{4}s+\frac{1}{2}\frac{1}{1+s}\bigg{)}F_{L}(s,Q^{2})-\frac{1}{4}F_{L}(0,Q^{2})$
(15)
$\displaystyle+C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\bigg{(}1-\frac{2}{1+s}-2H_{s}\bigg{)}F_{2}(s,Q^{2})+2C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\frac{H_{s}}{s}F_{2}(0,Q^{2}),$

where $H_{s}=\psi(s+1)-\psi(1)$ with $\psi(s)=\Gamma^{\prime}(s)/\Gamma(s)$.
Here $\psi(1)=-\gamma=-0.5772156...$ is Euler constant. The functions
$F_{i}(0)$ (or exactly $F_{i}(+0)$) are the boundary conditions due to the
Laplace derivatives and defined
$F_{i}(\upsilon=0,Q^{2}){\equiv}F_{i}(x=1,Q^{2})=0$. Then, from Eq. (15), we
find

$\displaystyle F_{L}(s,Q^{2})$ $\displaystyle=$ $\displaystyle
h_{1}(s)\frac{d{F}_{2}(s,Q^{2})}{d{\ln}Q^{2}}+{h_{2}(s)}F_{2}(s,Q^{2}),$ (16)

where $h_{1}(s)=h_{L}^{-1}(s)$ and
$h_{2}(s)=-C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\tilde{h}_{2}(s)h_{L}^{-1}(s)$
with $h_{L}(s)$ and $\tilde{h}_{2}(s)$ given by

$\displaystyle h_{L}(s)$ $\displaystyle=$
$\displaystyle\frac{2+s}{4}+\frac{1}{2(1+s)},$ $\displaystyle\tilde{h}_{2}(s)$
$\displaystyle=$ $\displaystyle 1-\frac{2}{1+s}-2H_{s}.$ (17)

The inverse Laplace transform of the coefficients $h_{i}$ in Eq. (16), defined
by the kernels
$\widehat{J}_{i}(\upsilon){\equiv}{\mathcal{L}^{-1}}[h_{i}(s);\upsilon]$, is
straightforward. We find that

$\displaystyle\widehat{J}_{1}(\upsilon)$ $\displaystyle=$ $\displaystyle
4e^{-\frac{3}{2}\upsilon}\bigg{[}\cos{\bigg{(}}\frac{\sqrt{7}}{2}\upsilon{\bigg{)}}-\frac{\sqrt{7}}{7}\sin{\bigg{(}}\frac{\sqrt{7}}{2}\upsilon{\bigg{)}}\bigg{]},$
(18)

and

$\displaystyle\widehat{J}_{2}(\upsilon)$ $\displaystyle=$
$\displaystyle-4C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}e^{-\frac{3}{2}\upsilon}\bigg{[}(1+2\psi(1))\cos{\bigg{(}}\frac{\sqrt{7}}{2}\upsilon{\bigg{)}}-(5+2\psi(1))\frac{\sqrt{7}}{7}\sin{\bigg{(}}\frac{\sqrt{7}}{2}\upsilon{\bigg{)}}\bigg{]}$
(19)
$\displaystyle+8C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\widehat{f}(\upsilon),$

where

$\displaystyle\widehat{f}(\upsilon)$ $\displaystyle=$
$\displaystyle{\mathcal{L}^{-1}}\bigg{[}\frac{(s+1)\psi(s+1)}{4+3s+s^{2}},s;\upsilon\bigg{]}$
(20)
$\displaystyle=-2e^{-\frac{3}{2}\upsilon}\bigg{[}0.1704\cos\bigg{(}\frac{\sqrt{7}}{2}\upsilon\bigg{)}+1.211\sin\bigg{(}\frac{\sqrt{7}}{2}\upsilon\bigg{)}\bigg{]}+\sum_{k=1}^{\infty}\frac{k}{(k+1)^{2}-3(k+1)+4}e^{-(k+1)\upsilon},$

as shown in the Appendix. Transforming back in to $x$ space, the longitudinal
structure function $F_{L}(x,Q^{2})$ is given by

$\displaystyle F_{L}(x,Q^{2})$ $\displaystyle=$ $\displaystyle
4\int_{x}^{1}\frac{d{F}_{2}(z,Q^{2})}{d{\ln}Q^{2}}(\frac{x}{z})^{3/2}\bigg{[}\cos{\bigg{(}}\frac{\sqrt{7}}{2}{\ln}\frac{z}{x}{\bigg{)}}-\frac{\sqrt{7}}{7}\sin{\bigg{(}}\frac{\sqrt{7}}{2}{\ln}\frac{z}{x}{\bigg{)}}\bigg{]}\frac{dz}{z}$
(21)
$\displaystyle-4C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\int_{x}^{1}F_{2}(z,Q^{2})(\frac{x}{z})^{3/2}\bigg{[}(1.6817+2\psi(1))\cos{\bigg{(}}\frac{\sqrt{7}}{2}{\ln}\frac{z}{x}{\bigg{)}}+(2.9542-2\frac{\sqrt{7}}{7}\psi(1))\sin{\bigg{(}}\frac{\sqrt{7}}{2}{\ln}\frac{z}{x}{\bigg{)}}\bigg{]}\frac{dz}{z}$
$\displaystyle+8C_{F}\frac{\alpha_{s}(Q^{2})}{2\pi}\sum_{k=1}^{\infty}\frac{k}{(k+1)^{2}-3(k+1)+4}\int_{x}^{1}F_{2}(z,Q^{2})(\frac{x}{z})^{k+1}\frac{dz}{z}.$

##  IV Results and Discussions

With the explicit form of the longitudinal structure function (i.e., Eq.
(21)), we begin to extract the numerical results at small $x$ in a wide range
of $Q^{2}$ values, using the parametrization of $F_{2}(x,Q^{2})$ given by Eq.
(3). The QCD parameter $\Lambda$ for four numbers of active flavor has been
extracted [10] due to $\alpha_{s}(M_{z}^{2})=0.1166$ with respect to the LO
form of $\alpha_{s}(Q^{2})$ with $\Lambda=136.8~{}\mathrm{MeV}$. In order to
make the effect of production threshold for charm quark with
$m_{c}=1.29^{+0.077}_{-0.053}~{}\mathrm{GeV}$ [42, 43], the rescaling variable
$\chi$ is defined by the form $\chi=x(1+4\frac{m_{c}^{2}}{Q^{2}})$ where
reduced to the Bjorken variable $x$ at high $Q^{2}$ [43].

Figure 1: The function $\widehat{f}(\upsilon)$ [i.e., Eq.(20)] is plotted in a
wide range of $\upsilon$. The solid black, dashed red, and dot-dashed green
curves are shown according to the expansion of the second term in Eq.(20) in
the ranges $k=1-2$, $k=1-10$, and $k=1-50$, respectively. The function
$\widehat{f}(\upsilon)$ diverges as ln $\upsilon$ for $\upsilon\to 0$.

In Fig. 1, we analyze the function $\widehat{f}(\upsilon)$ in a wide range of
$\upsilon$ according to the expansion of the second term in Eq. (20) in the
ranges $k=1-2$, $k=1-10$, and $k=1-50$, respectively. We observe that the
function $\widehat{f}(\upsilon)$ is very small as would be expected from the
decreasing exponential factor in $\upsilon$ in Eq. (20). It is nearly
independent of the cutoff in the expansion for $\upsilon\geq 1$, but the
expansion must be carried to large $k$ for $\upsilon$ small. In the following
we choose the value of $k=50$ in the numerical results.

In Fig. 2, using the parametrization for $F_{2}(x,Q^{2})$ in Eq. (3), we have
plotted the longitudinal structure function $F_{L}(x,Q^{2})$ from Eq. (21) as
a function of $x$ ($5{\times}10^{-5}<x<10^{-2}$) for the values of
$Q^{2}=5,15,25$, and $45~{}\mathrm{GeV}^{2}$, respectively. In the figure, the
H1 Collaboration data (H1 2011 [36] and H1 2014 [42]) accompanied by total
errors for $Q^{2}=5,15,25$, and $45~{}\mathrm{GeV}^{2}$ are also shown.

In Fig. 3, we have separated our analysis of the longitudinal structure
function at any fixed $Q^{2}$ and compare it with the results in [8] and [10]
at the LO approximation as a function of $x$. The longitudinal structure
function extracted at $Q^{2}$ values are in good agreement with experimental
data in comparison with those in [8] at the LO approximation, as the
mathematical structure of Eq. (8) in momentum space differs from the DGLAP
equations for structure functions.

Figure 2: Plots of the longitudinal structure function $F_{L}(x,Q^{2})$ versus
$x$ for $Q^{2}=5$ (solid brown), $15$ (dash green), $25$ (dash-dot turquoise),
and $45~{}\mathrm{GeV}^{2}$ (dash-dot-dot purple), respectively. H1
Collaboration data are selected from [36](H1 2011) and [42](H1 2014).  Figure
3: The longitudinal structure function $F_{L}(x,Q^{2})$ (solid curves) plotted
at fixed $Q^{2}$ as a function of $x$ variable, compared with the Mellin
transforms method [8] (dot curves) at the LO approximation. Experimental data
(up-triangle H1 2014, down-triangle H1 2011) are from the H1-Collaboration
[42, 36] as accompanied with total errors.  Figure 4: The central values of
the longitudinal structure function are plotted at low $x$ as accompanied by
the LHeC total errors [13]. H1 Collaboration data [36, 42] are collected at
$Q^{2}=8.5$ and $35~{}\mathrm{GeV}^{2}$ with total errors.

In Fig. 4, the longitudinal structure functions in momentum space at selected
$x$ and $Q^{2}$ are associated with the LHeC simulated uncertainties [13]. We
observe that the longitudinal structure functions (central values) are
determined owing to Eq. (21) for the $Q^{2}$ values (4.5, 8.5, 18, and
35$~{}\mathrm{GeV}^{2}$) and accompanied with the simulated uncertainties
reported by the LHeC study group [13]. The H1 collaboration data with total
errors for $Q^{2}=8.5$ and $35~{}\mathrm{GeV}^{2}$ are shown in Fig.4.

In Fig. 5, we show a comparison between the longitudinal structure functions
in momentum space with the H1 Collaboration data at a fixed value of the
invariant mass $W$ (i.e. $W=230~{}\mathrm{GeV}$) at low values of $x$. Figure
5 clearly demonstrates that the Laplace transform method in momentum space
provides correct behaviors of the extracted longitudinal structure function in
comparison with the LO and NLO analysis reported in [10]. As can be seen in
this figure, the results are comparable with the H1 data and the NLO
corrections to the Mellin transform method at all $Q^{2}$ values.

Figure 5: The extracted longitudinal structure function (solid brown curve) in
momentum space at fixed value of the invariant mass $W$
($W=230~{}\mathrm{GeV}$) compared with the H1 Collaboration data [42] as
accompanied with total errors and the results in Refs. [8] and [10] at the LO
(dotted black curve) and NLO (dot-dashed green curve) approximation.

##  V Conclusions

We have presented a method based on the Laplace transform method to determine
the longitudinal structure function at the LO approximation in momentum space.
This method relies on the parametrization of the function $F_{2}(x,Q^{2})$ and
its derivative $dF_{2}(x,Q^{2})/d{\ln}Q^{2}$ within a kinematical region
characterized by low values of the Bjorken variable $x$. The $x$ dependence of
$F_{2}(x,Q^{2})$ and its evolution with $Q^{2}$ are determined much better by
the data than $F_{L}(x,Q^{2})$, so this method provides both a direct check on
$F_{L}(x,Q^{2})$ where measured, and a way of extending $F_{L}(x,Q^{2})$ into
regions of $x$ and $Q^{2}$ where there are currently no data. We find that the
Laplace transform method in momentum space provides correct behaviors of the
extracted longitudinal structure function $F_{L}(x,Q^{2})$ and that our
results for $F_{L}(x,Q^{2})$ demonstrate comparability with data from the H1
Collaboration and other results obtained using the Mellin transform method.

##  VI ACKNOWLEDGMENTS

G.R.Boroun thanks M.Klein and N.Armesto for allowing access to data related to
simulated errors of the longitudinal structure function at the Large Hadron
Electron Collider (LHeC). Phuoc Ha would like to thank Professor Loyal Durand
for useful comments and invaluable support.

##  VII Appendix

The inverse Laplace transform

$\widehat{f}(\upsilon)={\mathcal{L}^{-1}}\bigg{[}f(s),s;\upsilon\bigg{]}=\frac{1}{2\pi
i}\int_{-i\infty}^{i\infty}dsf(s)e^{s\upsilon}$ (22)

of the function

$f(s)=\frac{(s+1)\psi(s+1)}{4+3s+s^{2}}$ (23)

can be evaluated analytically in terms of an infinite series, rapidly
convergent except for $\upsilon$ near zero where it grows logarithmically. The
denominator in $f(s)$ has zeros at $s_{\pm}=-\frac{3}{2}\pm
i\frac{\sqrt{7}}{2}$. which lead to simple poles in $f(s)$ at those points.
The function $\psi(s+1)$ has simple poles with residue $-1$ at
$s+1=0,-1,-2,...$ [44] . There are no other sigularities in the integrand
which decreases exponentially rapidly for $s\rightarrow-\infty$ and
Re$(\upsilon)>0$. We can therefore close the integration contour in the left-
half $s$ plane and evaluate the integral as the sum of the residues at the
poles multiplied by $2\pi i$ by Cauchy’s residue theorem.

This gives

$\widehat{f}(\upsilon)=\frac{s_{+}+1}{s_{+}-s_{-}}\psi(s_{+}+1)e^{s_{+}\upsilon}-\frac{s_{-}+1}{s_{+}-s_{-}}\psi(s_{-}+1)e^{s_{-}\upsilon}+\sum_{k=1}^{\infty}\frac{k}{(k+1)^{2}-3(k+1)+4}e^{-(k+1)\upsilon}.$
(24)

Since

$\frac{s_{+}+1}{s_{+}-s_{-}}\psi(s_{+}+1)=-0.1704+1.211i\,,\,\,\,\,\,\frac{s_{-}+1}{s_{+}-s_{-}}\psi(s_{-}+1)=0.1704+1.211i.$
(25)

we find

$\widehat{f}(\upsilon)=-2e^{-\frac{3}{2}\upsilon}\bigg{[}0.1704\cos\big{(}\frac{\sqrt{7}}{2}\upsilon\big{)}+1.211\sin\big{(}\frac{\sqrt{7}}{2}\upsilon\big{)}\bigg{]}+\sum_{k=1}^{\infty}\frac{k}{(k+1)^{2}-3(k+1)+4}e^{-(k+1)\upsilon}.$
(26)

For $\upsilon\ll 1$, $k$ large, the series is approximately

$\sum_{k}\frac{1}{k}e^{-(k+1)\upsilon}\approx e^{-\upsilon}\,{\rm
ln}\,(1-e^{-\upsilon}),$ (27)

so diverges as ln $\upsilon$ for $\upsilon\to 0$.

Table 1:  The effective parameters [7] at low $x$ for $0.15~{}\mathrm{GeV}^{2}<Q^{2}<3000~{}\mathrm{GeV}^{2}$. parameters value |  |
---|---|---
$a_{0}$ |  $8.205\times 10^{-4}~{}~{}\pm 4.62\times 10^{-4}$ |  |
$a_{1}$ |  $-5.148\times 10^{-2}\pm 8.19\times 10^{-3}$ |  |
$a_{2}$ |  $-4.725\times 10^{-3}\pm 1.01\times 10^{-3}$ |  |
$b_{0}$ |  $2.217\times 10^{-3}\pm 1.42\times 10^{-4}$ |  |
$b_{1}$ |  $1.244\times 10^{-2}\pm 8.56\times 10^{-4}$ |  |
$b_{2}$ |  $5.958\times 10^{-4}\pm 2.32\times 10^{-4}$ |  |
$c_{1}$ |  $1.475\times 10^{-1}~{}\pm 3.025\times 10^{-2}$ |  |
$n$ |  $11.49\pm 0.99$ |  |
$\lambda$ |  $2.430~{}\pm 0.153$ |  |
$M^{2}$ |  $0.753\pm 0.068~{}\mathrm{GeV}^{2}$ |  |
$\mu^{2}$ |  $2.82\pm 0.290~{}\mathrm{GeV}^{2}$ |  |
$c_{0}$ |  $0.255\pm 0.016$ |  |
$\chi^{2}(\mathrm{goodness~{}of~{}fit})$ |  $0.95$ |  |

## References

  * [1] T. Lappi, H. Mantysaari, H. Paukkunen and M.Tevio, Eur.Phys.J.C 84, 84 (2024).
  * [2] V. N. Gribov and L. N. Lipatov, Sov.J.Nucl.Phys. 15, 438 (1972).
  * [3] G. Altarelli and G. Parisi, Ncul.Phys.B 126, 298 (1977).
  * [4] Y. L. Dokshitzer, Sov.Phys.JETP 46, 641 (1977).
  * [5] B. Badelek, et al., J.Phys.G 22, 815 (1996).
  * [6] S. Catani, F. Hautmann, Nucl.Phys.B 427, 475 (1994).
  * [7] M. M. Block, L. Durand and P. Ha, Phys.Rev.D 89, 094027 (2014).
  * [8] L.P. Kaptari, et al., JETP Lett. 109, 281 (2019).
  * [9] G.R. Boroun, Phys.Rev.C 97, 015206 (2018).
  * [10] L.P. Kaptari, et al., Phys.Rev.D 99, 096019 (2019).
  * [11] G.R.Boroun and B.Rezaei, Phys.Lett.B 816, 136274 (2021).
  * [12] M. Froissart, Phys. Rev. 123, 1053 (1961).
  * [13] P.Agostini et al. [LHeC Collaboration and FCC-he Study Group], J. Phys. G: Nucl. Part. Phys. 48, 110501 (2021).
  * [14] R. Abdul Khalek et al., Snowmass 2021 White Paper, arXiv [hep-ph]:2203.13199.
  * [15] R.Abir et al., The case for an EIC Theory Alliance: Theoretical Challenges of the EIC, arXiv [hep-ph]:2305.14572.
  * [16] M. Gluck, E. Reya, and A. Vogt, Z. Phys. C 48, 471 (1990).
  * [17] W. Furmanski and R. Petronzio, Nucl. Phys. B 195, 237 (1982).
  * [18] R. Toldra, Comput. Phys. Commun. 143, 287 (2002).
  * [19] N. Cabibbo and R. Petronzio, Nucl. Phys. B 137, 395 (1978).
  * [20] J.Rausch, V.Guzey and M.Klasen, Phys.Rev.D 107, 054003 (2023).
  * [21] G.R.Boroun and B.Rezaei, Phys.Rev.D 105, 034002 (2022).
  * [22] H.Khanpour, A.Mirjalili and S.Atashbar Tehrani, Phys.Rev.C 95, 035201 (2017).
  * [23] G.R.Boroun, Eur.Phys.J.Plus 137, 32 (2022).
  * [24] S.Atashbar Tehrani, F.Taghavi-Shahri, A.Mirjalili, M.M.Yazdanpanah, Phys.Rev.D 87, 114012 (2013).
  * [25] S.Dadfar and S.Zarrin, Eur.Phys.J.C 80, 319 (2020).
  * [26] M.Mottaghizadeh, P.Eslami and F.Taghavi-Shahri, Int.J.Mod.Phys.A 32, 1750065 (2017).
  * [27] N.Olanj, M.Lotfi Parsa and L.Asgari, Phys.Lett.B 834, 137472 (2022).
  * [28] R. D. Ball and S. Forte, Phys. Lett. B 336, 77 (1994).
  * [29] A. V. Kotikov and G. Parente, Nucl. Phys. B 549, 242 (1999).
  * [30] L. Mankiewicz, A. Saalfeld, and T. Weigl, Phys. Lett. B 393, 175 (1997).
  * [31] M. Markovych and A. Tandogan, arXiv [hep-ph]:2304.10458.
  * [32] A. Simonelli, arXiv [hep-ph]: 2401.13663.
  * [33] M. M. Block and F. Halzen, Phys. Rev. D 70, 091901 (2004).
  * [34] H. Abramowicz et al. (H1 and ZEUS Collaborations), Eur. Phys. J. C 78, 473 (2018).
  * [35] V. Andreev et al. (H1 Collaboration), Eur. Phys. J. C 74, 2814 (2014).
  * [36] F. D. Aaron et al. (H1 Collaboration), Eur. Phys. J. C 71, 1579 (2011).
  * [37] Martin M. Block, Loyal Durand and Douglas W. McKay, Phys.Rev.D 79, 014031 (2009).
  * [38] Martin M. Block, Loyal Durand, Phuoc Ha and Douglas W. McKay, Phys.Rev.D 83, 054009 (2011).
  * [39] Martin M. Block, Loyal Durand, Phuoc Ha and Douglas W. McKay, Phys.Rev.D 84, 094010 (2011).
  * [40] Martin M. Block, Loyal Durand, Phuoc Ha and Douglas W. McKay, Phys.Rev.D 88, 014006 (2013).
  * [41] H1 and ZEUS Collaborations (H. Abramowicz et al.), Eur. Phys. J. C 78, 473 (2018).
  * [42] H1 Collab. (V.Andreev et al.), Eur.Phys.J.C 74, 2814(2014).
  * [43] M.A.G.Aivazis et al., Phys.Rev.D 50, 3102 (1994).
  * [44] ”NIST Digital Library of Mathematical Functions” edited by F. W. J. Olver, A. B. Olde Daalhuis, D. W. Lozier, B. I. Schneider, R. F. Boisvert, C. W. Clark, B. R. Miller, B. V. Saunders, H. S. Cohl, and M. A. McClain, Sec. 5.2(i), https://dlmf.nist.gov/, Release 1.1.12 of 2023-12-15.


# People Counting using Radio Irregularity in Wireless Sensor Networks – An
Experimental Study

Wei-Chuan Lin, Winston K.G. Seah and Wei Li  Wei-Chuan Lin, and Winston K.G.
Seah are with the School of Engineering and Computer Science, Victoria
University of Wellington, P.O. Box 600, Wellington 6140, New Zealand. E-mail:
{wei-chuan.lin<EMAIL_ADDRESS>Wei Li is with the Department of
Electrical and Computer Engineering, University of Victoria, P.O. Box 3055,
Victoria, B.C. V8W 3P6, Canada. Email<EMAIL_ADDRESS>

###### Abstract

The Internet has grown into a large cyber-physical system centered that
connects not just computer systems but a plethora of systems, devices, and
objects, collectively referred to as “Things”, giving rise to the term
“Internet of Things” (IoT). It encompasses technologies for identification and
tracking, sensing and actuation, both wired and wireless communications, and
also, intelligence and cognition. Wireless communications, which is an
integral part of IoT, suffers from radio irregularity – a phenomenon referring
to radio waves being selectively absorbed, reflected or scattered by objects
in their paths, e.g., human bodies that comprises liquid, bone and flesh.
Radio irregularity is often regarded as a problem in wireless communications
but, with the envisioned pervasiveness of IoT, we aim to exploit radio
irregularity as a means to detect and estimate the number of people. We
demonstrate how radio signal fluctuations arising from radio irregularity,
combined with discriminant analysis, can be used to provide a simple low-cost
alternative to dedicated sensing systems for indoor people counting.

###### Index Terms:

People Counting; Radio Irregularity; Internet of Things; Wireless Sensor
Networks; Discriminant Analysis;

##  I Introduction

The Internet has grown beyond connecting computer systems and platforms that
run applications to meet endusers’ computing and communication needs to
connecting a plethora of systems, appliances, devices, objects, etc.,
collectively referred to as “Things”, turning it into a large cyber-physical
system and giving rise to a new paradigm known as the “Internet of Things”
(IoT) [1]. Likewise, the technologies that the IoT encompasses extend beyond
computation and communication, to identifcation and tracking, sensing and
actuation, and even intelligence and cognition. Ensuring connectivity in the
IoT is increasingly reliant on wireless communications as connected devices
become more ubiquitous and embedded into our daily living space.

The key technology in wireless communications is radio frequency (RF)
transmission. When an RF signal propagates within a medium, it may be
reflected, diffracted, and scattered. Each effect occurs to a different extent
in various media, depending on factors such as wavelength and intensity of the
wave, thickness and physical composition (permittivity and permeability) of
the medium. The human body comprises liquid, bone and flesh, which selectively
absorb, reflect or scatter RF signals, leading to the phenomenon known as
radio irregularity. Consequently, in the presence of human activity within a
network, the radio irregularity phenomenon is seen as signal strength
fluctuations at the receiver, and the degree of signal fluctuation exhibits a
significant level of correlation to the level of human activity in the
network. Researchers have exploited this phenomenon in intrusion detection [2]
and device-free localization applications [3, 4].

Applications like automated people counting cannot tolerate false positives
that result in overcounting, giving inaccurate data that are used for
forecasting and resource allocation. People counting is extensively used in
different industries, including retail (stores, malls and shopping centres),
colleges and universities, government facilities, government non-profits
organizations, visitor centres, libraries, museums and art galleries. In the
retail industry, it is a form of intelligence-gathering that helps a retailer
determine the percentage of visitors who actually make purchases. This is a
key performance indicator of a store’s performance as compared to just looking
at the sales data. It also helps the management to optimize the usage of staff
resources, e.g. deploy more staff during peak periods and cutting down during
lull periods in order to save wages. For building automation and management
purposes, people counting is used to optimize the use of resources as well as
ensure that the safe level of occupancy is maintained.

With the emergence of IoT leading to pervasive wireless communication devices,
radio irregularity which has often been viewed as a problem can instead be
exploited for automated people counting with minimal additional hardware and
installation costs. The goal of our approach is to provide a simple low-cost
system for people counting that can be used independently or as a supplement
to other more sophisticated methods. In the next section, we examine the
related research on automated people counting with a focus on indoor use
cases. We then present our approach to indoor automated people counting based
on the signal fluctuations arising from radio irregularity. This is followed
by the discussion of the experimental study and the application of
discriminant analysis to process the results obtained from tests carried out
indoors within a building before concluding the paper.

##  II Related Work

###  II-A People Counting Methods

Widespread usage of people counting has fed the growth of commercially
available automated people counting technology, among which infrared (IR) beam
counters, thermal counters and video/CCTV cameras are most often used. The
simplest and possibly cheapest approach is a single-beam IR counter placed
across an entrance. However, such a counter suffers from numerous drawbacks
and is only suitable detecting someone passing, e.g. entering/leaving a shop.
These commonly used counters have very high percentage of errors when multiple
persons cross their monitoring area at a time. When multiple (IR) beams or
other forms of boundary sensors are deployed with careful placements
strategies and coupled with artificial intelligence and/or analytical
techniques for processing, a more accurate and versatile people counting
system can be realized [5, 6].

People counters that use thermal imaging are typically mounted overhead and
have the ability to simultaneously maintain separate counts for multiple
people moving in two directions (in and/or out). The IR images captured by the
heat detectors are then processed to determine the number of people [7].
Video-based people counters work on video streams obtained through video/CCTV
camera which are then run through intelligent video-processing techniques to
identify and count the people in the video [8]. The accuracy of such
approaches can vary according to the level of ambient lighting and background
colour contrasts [9]. Hybrid approaches combining IR and video cameras,
together with neural networks, have been proposed to improve the accuracy of
visual-based automated people counting [10].

###  II-B Radio-based Detection and Counting Methods

It was first reported in [11] that the shadowing effect caused by an object
moving between two communicating wireless devices can be used for detection
purposes. In particular, a human body comprises liquid, bone and flesh, that
selectively absorb, reflect or scatter RF signals, leading to the phenomenon
known as radio irregularity. This phenomenon has been extensively used for
device-free localization in wireless networks [3, 4]. Radio Tomographic
Imaging [12] measures the attenuation of signals across wireless links between
many pairs of nodes in a wireless network to create “images” of objects moving
within the network area. The variance of measured received signal strength
(RSS) on the links in a network has also been used to infer the locations of
people or objects moving in the network deployment area [13]. This approach
utilizes a statistical model for the RSS variance as a function of a person’s
position with respect to the transmitter and receiver locations. The approach
adopted by [11] has also been extended in [14] for outdoor people counting by
measuring the RSS level measured at the receiver. The reliance on (absolute)
RSS values, however, has a drawback during deployment, which is the need to
take into consideration the channel model and other related factors like the
impact of path loss and fading. These approaches also require complex signal
processing techniques and a calibration phase for each deployment environment,
and this can significantly affect the ease of deployment. Our targetted
application differs from the proposed device-free localization schemes [3, 4]
as we are not aiming to localize or track individual objects.

It has been observed in [2] that human movement through the path of the radio
signal causes the histogram of the absolute RSS values to become more spread;
this is manifested quantitatively as higher standard deviation. However, the
standard deviation varies significantly across environments, making it
difficult to define a universal threshold to detect movement in terms of these
first order statistics. While also exploiting the RSS spread caused by human
movement, the approach adopted focused on the fluctuation in signal strength
instead, in order to reduce the impact of channel models and other
environmental factors. However, there are false positives reported in their
results which are deemed to be acceptable for intrusion detection
applications. This has since been improved by eliminating the false positives
and extended to detect two persons crossing a narrow corridor [15].

##  III Counting People using RSSI Fluctuations

###  III-A RSSI fluctuations caused by human activity

Most, if not all, the approaches that use the variance of RSS caused by human
motion across the signal transmission paths require complex signal processing
techniques and a calibration phase for each deployment environment, and this
can significantly affect the ease of deployment. In our scheme, we adopt a
network-oriented approach that relies on RSS information of received packets
which can be easily obtained from device drivers of wireless network
interfaces when the packets are received. A key goal of our approach is to be
able to easily utilize the existing wireless transmitters and receivers
already deployed in the environment without the need for accurate channel
models nor complex signal processing techniques. We extend the method of using
Received Signal Strength Indicator (RSSI) fluctuations proposed in [2] where
the RSSI fluctuation for a packet $p_{i}$, denoted by $F(p_{i})$, is given by
$F(p_{i})=RSSI(p_{i})-RSSI(p_{i-1})$. It has been shown that two consistent
patterns of RSSI fluctuations can be observed for two key scenarios of
interest to us, namely, without human movement and with human movement across
the signal transmission path. Given the absolute RSSI readings for packets
recorded at the receiver over time (Fig. 1), the histogram of RSSI
fluctuations shows narrower distribution when there is no human movement
across the signal path, i.e., less fluctuation across RSSI readings (Fig. 2a)
and, conversely, signals fluctuate more in the presence of human movement
resulting in the spread out distribution shown in Fig. 2b.

Figure 1: Absolute RSSI reading

(a) RSSI Fluctuation without movement

(b) RSSI Fluctuation with movement

Figure 2: RSSI Fluctuation Patterns [2]

###  III-B Human Detection

Our proposed algorithms computes the fluctuation between the RSSI of packets
received at a receiver since the absolute RSSI largely depends on external
environment. From the absolute RSSI readings for packets recorded at the
receiver over a period of time as shown in Fig. 1, the fluctuation of RSSI
readings is calculated, as shown in Fig. 3.

Figure 3: RSSI Fluctuation

We define a sliding window of $n$ samples, where $n$ is a parameter that can
be tuned to achieve the desired accuracy for the target environment. The
desired value for $n$ would be one that achieved the lowest error detection
rate, which is calculated by dividing the number of incorrect detections by
the total number of samples. We studied the accuracy of detection using
different sliding window sizes, from $n=5$ to $n=100$, across different
environments, and also comparing with other schemes, namely, the probability
density function approach [15] and the overcomplete dictionary approach [16],
to obtain the results shown in Fig. 4.

Figure 4: Optimal Window Size

Based on the results shown in Fig. 4, we use a sliding window of size $n=10$
to observe the behaviour of RSSI fluctuation. Therefore, a window of RSSI
fluctuations at sample 200 is shown in Fig. 5.

Figure 5: RSSI fluctuations over a window size of 10

At sample 200, using the window of 10 previous readings, the mean and standard
deviation are computed as 0.2727 and 4.6280 respectively. We then map the RSSI
fluctuations into the normal distribution with the mean and standard deviation
for that window, i.e. $\mu=0.2727$ and $\sigma=4.6280$, as shown in Fig. 6a
representing the case where the signal has been subjected to interference by
human movement across its path. Similarly, the normal distribution of RSSI
fluctuation at sample 600, where there is no movement, is shown for comparison
in Fig. 6b.

(a) Sample 200 (movement)

(b) Sample 600 (no movement)

Figure 6: Normal distribution showing probability in fluctuation range [-1,1]

From the graphs, we compute the probability of the RSSI fluctuation falling
within the range [-1,1] (i.e. area under the curve from -1 to 1) to be 0.17078
for the case where there is movement across the signal path (i.e. sample 200)
and 0.84303 for the case where there is no movement (sample 600). For the
dataset shown in Fig. 1, we compute the probability of falling with the
fluctuate range [-1,1] and plot the results in Fig. 7. As shown, the
probability of fluctuations falling in the range of [-1, 1] is below 0.3 in
the presence of human movement. Hence, a probability value that is higher than
0.3 implies no human movement. Based on this threshold, we then infer from the
results whether or not there has been human movement across the signal path,
and the results are shown in Fig. 8.

Figure 7: Probability of fluctuation within [-1,1] for RSSI readings in Fig.1
Figure 8: Inferred presence of human movement using RSSI Fluctuations

The approach used in [2] has resulted in false positives as shown in Fig. 9a.
We applied our approach to the dataset used by the detection algorithm [2]
that produced the results shown in Fig. 9a, and confirmed that our algorithm
is able to achieve better accuracy in eliminating false positives, as shown in
Fig. 9b.

(a) False positives by detection method in [2]

(b) Movement detected by our method

Figure 9: Detection results using dataset of [2]

##  IV Pedestrian Traffic Monitoring

Next, we extend the detection algorithm to monitor pedestrian traffic with the
aim of counting the number of people passing through a particular location of
interest.

###  IV-A Single transmitter-single receiver configuration

First, a series of experiments were conducted to observe the precision of the
detection algorithm in a realistic indoor environment, namely, a corridor in a
university building, as shown in Fig. 10, where the two red dots indicated by
the arrows refer to the transmitter/receiver pair using IEEE802.15.4
technology. The devices are spaced 1.5m apart (width of corridor) and placed
at a height of 1.1m, on a ledge. Each data collection duration was 300 seconds
with inter-packet interval time of 0.15 seconds, during which the number of
people who have walked past the devices were recorded and tagged with the
time. Fig. 11 shows the results for one data collection period, during which
nine persons walked through individually and two pairs of people past while
walking close to each other, at the sample index of 484 and 925. In the
detection results, shown in Fig. 11, a total of 11 movements were detected. It
is clear that detecting two people walking side by side is a major challenge
as the RSSI fluctuations arising from one and two persons passing are quite
indistinguishable.

Figure 10: Deployment along corridor of building in university Figure 11:
Detection of pedestrian traffic along corridor

###  IV-B Single transmitter-multiple receiver configuration

A cyber-physical system like IoT is made up of numerous small wireless
devices, from which a conceptual deployment scenario like that shown in Fig.
12a can be assumed; we first look at a subset configuration of one-transmitter
and two-receivers as shown in Fig. 12b [15]. As people walk along the path
between T and the two receivers in the direction shown in Fig. 12b, they first
cross the T-R2 signal transmission path, followed by the T-R1 signal path. A
key point to note is the different signal interference zones that result from
the movement of the people.

(a) Conceptual Configuration

(b) One-Transmitter Two-Receiver

Figure 12: Single Transmitter with Multiple Receivers

Using the one-transmitter two-receiver configuration, the transmitter
broadcasts packets at a rate of one packet every 0.15 seconds. Receiver R1 is
1.5m from transmitter T (d1 = 1.5m) and R2 is 1.5m from R1 (d2 = 1.5m). As two
persons walk along the path between T and the two receivers in the direction
shown in Fig. 12b, they first cross the T-R2 signal transmission path,
followed by the T-R1 signal path.

First, we collected data for one person walking across the signal transmission
path, passing first R2 then R1 to be used as the reference case. The detection
results correctly show that one person passed at around the time of sample 100
and another at around sample 200, as shown in Fig. 13. Intuitively, the
detection result at sample 100 is more logical since the person passed R2
first, then R1. However, as the two receivers at very close to each other,
having the two receivers showing signal fluctuations at almost the same
instant is also likely especially when the person is walking fast.

Next, we collected data for the case of two persons walking side-by-side in
the direction of R2 to R1 as shown in Fig. 12b. We expect that the detection
duration of T-R2 should be longer than T-R1. This is because the T-R2 signal
experienced a longer duration of interference than the T-R1 signal. The
detection result of two people walking from R2 to R1 shown in Fig. 14 confirms
our hypothesis. However, we also observed a false positive detection at sample
64. As the two receivers are placed closed to each other, 1.5m apart, we can
assume that it is unlikely for a moving object to be detected by one receiver
but not the other. Therefore, by comparing and matching the data from both
receivers, we can perform a simple optimization process to remove such false
positive detections, to achieve the desired results as shown in Fig. 15.

Figure 13: One person walking in the direction of R2 to R1 Figure 14: Two
people walking in the direction of R2 to R1 Figure 15: Optimised Result of
multiple receivers

##  V Enhanced Accuracy for People Counting

###  V-A Standard Deviation of RSSI fluctuation Approach

Using the probability density function of RSSI fluctuations falling within the
range [-1,+1] as the threshold to signify no human movement has eliminated the
occurrence of false positives [15], but it disregards information from the
distribution of RSSI fluctuations that lie outside the region between -1 and
+1. This is undesirable as the distribution of RSSI fluctuations has been
shown to be a good indication of the size or crowd density of moving objects
[17].

Figure 16: Standard Deviation of fluctuation for RSSI readings in Fig.1 Figure
17: Inferred presence of human movement using standard deviation

According to Fig. 4, we use a sliding window of size $n=10$ to observe the
behaviour of RSSI fluctuation. Signal interference due to human motion causes
rapid RSSI fluctuations which results in an increased standard deviation. For
the dataset shown in Fig. 1, we compute the standard deviation of RSSI
fluctuation and the computed results are shown in Fig. 16; at sample 200
(human presence), the standard deviation (of the most recent 10 RSSI
fluctuation readings) is 4.6280 and at sample 600 (no human presence), the
standard deviation is 0.6325. The standard deviation of RSSI fluctuations
being higher than 2 implies the presence of human movement, which provides a
clearer threshold than the approach presented in [15]. We infer the existence
of human movement based on the standard deviation threshold and the results
are shown in Fig. 17. More importantly, additional information can be derived
from the data, such as, the peak of standard deviation.

###  V-B Multiple Object Detection Setup

In order to exploit the effect of RSSI fluctuation caused by signal
interference by a few people, we conducted a series of experiments with a
variable group size (ranging from one to five people) under a controlled
indoor environment. Sensor motes using IEEE802.15.4 wireless technology were
deployed in a $6m{\times}8m$ room in a one-transmitter two-receiver
configuration (Fig. 12b.) The motes were placed at a height of 0.9m. Receiver
R1 was placed 3.5m from transmitter T (d1 = 3.5m) and R2 at 2m from R1(d2 =
2m). The transmitter, T, broadcasts packets continuously in time intervals of
0.15s. The absolute RSSI values were recorded upon packet reception at the
receivers. Then, groups comprising one to five people made five consecutive
round trips between the transmitter and two receivers. Five experiments were
conducted for each group.

###  V-C Data Analysis

We use discriminant analysis, a method to find the linear combination of
measurements which characterize two or more groups [18], to analyze the
collected data. The key concept of the discriminant analysis approach is to
classify the number of people based on the difference in the influence on RSSI
readings due to the size of interference zone. We let a finite number $g$
denote the distinct number of groups which in our case is five (i.e. $g$=5,
denoting groups of 1–5 persons). We refer to the $G_{i}$ as groups where
$i=1,\ldots,5$. There are two phases of discriminant analysis which are
training and classification. The system identifies differences in RSSI
fluctuations caused by different number of people as signatures in the
training phase. It is important to note that this training is only required
once to identify the fluctuations caused by the different number of people,
and not for each deployment scenario. During the training phase, a total of
$g-1$ orthogonal discriminant functions are constructed such that the groups
differ as much as possible on discriminant score $D$. The form of the linear
discriminant function is:

$D=v_{1}X_{1}+v_{2}X_{2}+...+v_{i}X_{i}+c$ (1)

where

$\displaystyle v$ $\displaystyle=\mathrm{the\ discriminant\ coefficient}$
$\displaystyle X$ $\displaystyle=\mathrm{the\ value\ of\ independent\
variable}$ $\displaystyle c$ $\displaystyle=\mathrm{a\ constant}$
$\displaystyle i$ $\displaystyle=\mathrm{the\ number\ of\ independent\
variables}$

Once the discriminant functions are constructed, the discriminant analysis
enters the second phase which is classification. Fisher’s linear discriminant
analysis [19] is used for data classification; the purpose of Fisher’s
technique is to find the line of projection that separates different groups
[20]. We use the standard deviation of RSSI fluctuations of detected movement
as the primary dataset. For instance, the standard deviation of RSSI
fluctuations of an experiment involving one person is shown in Fig. 18. We
utilize the information from each positive detection, particularly mean,
standard deviation (std), coefficient of variation (CV), detection duration
and area under the curve, to be the independent variables in order to achieve
high discrimination between groups. The more statistical information we can
extract from these positive detections, the greater the ability to
discriminate between the different size groups.

For a total of 50 samples of each group, the mean of each independent variable
is plotted in Fig. 19. Duration and area under the curve are the two most
significant independent variables. The detection duration of R1 is stable
throughout all groups while the detection duration of R2 steadily increases as
the number of people increase. This trend is as expected for more time is

Figure 18: Experimental Data for one person in a single experiment Figure 19:
Mean of independent variables of each group

needed to pass the T-R2 transmission path than T-R1. The area under the curve
follows the same trends as longer detection duration makes the area under the
curve larger.

Next, these independent variables are taken to construct the discriminant
functions that maximize the separation between each group. In Table I, it
provides very high F values and lower Wilks’ Lambda as evidence of significant
difference in R2 detection duration than any other independent variables. In
addition, in statistical significance testing, the null hypothesis is rejected
when the $p$-value is smaller than the significance level $\alpha$ which is
0.05. The results are considered to be statistically significant when null
hypothesis is rejected. The information of each discriminant function is shown
in Table II. There are five groups, namely ‘one person’ to ‘five people’ and
as a result four discriminant functions are produced. With high Eigenvalue and
percentage of variance, ‘Function 1’ covers most total statistical population.

Independent Variable | Wilks’ Lambda | F |  $p$-value
---|---|---|---
R1 mean | 0.773 | 17.986 | 0.001
R2 mean | 0.992 | 5.207 | 0.000
R1 std | 0.640 | 34.437 | 0.000
R2 std | 0.692 | 27.311 | 0.000
R1 CV | 0.485 | 65.034 | 0.000
R2 CV | 0.436 | 79.137 | 0.000
R1 duration | 0.633 | 35.572 | 0.000
R2 duration | 0.256 | 178.008 | 0.000
R1 area | 0.713 | 24.605 | 0.000
R2 area | 0.603 | 40.256 | 0.000
TABLE I: Test of Equaltiy of Group Means Functions | Eigenvalue | % of variance | Canonical Correlation
---|---|---|---
1 | 5.717 | 82.9 | 0.923
2 | 0.952 | 13.8 | 0.698
3 | 0.146 | 2.1 | 0.357
4 | 0.084 | 1.2 | 0.279
TABLE II: Table of Eigenvalues

The relative importance of each independent variable in each discriminant
function can be found by analyzing the structure matrix table (Table IV).
These ‘Pearson’ coefficients are discriminant loadings which act like factor
loadings in factor analysis [21]. Generally, an absolute value of 0.3 is used
as the threshold to separate a significant from an insignificant variable. For
example, we have five independent variables in discriminant function $D_{1}$,
namely, ‘R1 CV’, ‘R2 CV’, ‘R1 duration’, ‘R2 duration’, and ‘R2 area’ that
discriminates between groups. The discriminant function coefficient shows the
contribution of each independent variable to the discriminant function (Table
IV). It operates like a regression equation.

| Function
---|---
| 1 | 2 | 3 | 4
R1 mean | 0.153 | 0.384 | -0.047 | 0.477
R2 mean | 0.094 | -0.088 | 0.243 | 0.464
R1 std | 0.275 | 0.361 | -0.051 | 0.235
R2 std | 0.269 | -0.109 | 0.049 | 0.486
R1 CV | 0.417 | 0.219 | 0.387 | -0.162
R2 CV | 0.468 | -0.093 | 0.070 | 0.603
R1 duration | 0.306 | 0.111 | 0.417 | -0.296
R2 duration | 0.692 | -0.411 | -0.206 | -0.098
R1 area | 0.191 | 0.438 | 0.020 | 0.372
R2 area | 0.326 | -0.212 | 0.094 | 0.288
TABLE III: Table of Eigenvalues

| Function
---|---
| 1 | 2 | 3 | 4
R1 mean | 0.260 | 0.863 | 0.611 | 0.084
R2 mean | -0.837 | -1.044 | 0.843 | 0.460
R1 std | -1.033 | -1.023 | -4.612 | -0.471
R2 std | 0.999 | 1.046 | -1.165 | -0.934
R1 CV | 17.142 | 17.026 | 32.236 | -12.682
R2 CV | 3.964 | -6.788 | 2.522 | 12.729
R1 duration | 0.298 | 0.225 | 0.514 | -0.358
R2 duration | 0.306 | -0.203 | -0.140 | -0.146
R1 area | 0.011 | 0.017 | -0.001 | 0.018
R2 area | 0.002 | -0.003 | 0.005 | 0.012
(Constant) | -10.805 | -2.414 | -10.974 | 2.054
TABLE IV: Canonical Discriminant Function Coefficient

For example, the discriminant function $D_{1}$ and $D_{2}$ are shown as
follows:

$\begin{split}D_{1}=&\left(0.26\times R_{1}\,mean)+(-0.837\times
R_{2}\,mean)+(-1.033\times R_{1}\,std)+(0.999\times R_{2}\,std)\right.\\\
&\left.+\,(17.142\times R_{1}\,CV)+(3.964\times R_{2}\,CV)+(0.298\times
R_{1}\,duration)\right.\\\ &\left.+\,(0.306\times
R_{2}\,duration)+(0.011\times R_{1}\,area)+(0.002\times
R_{2}\,area)-10.805\right.\end{split}$ (2)
$\begin{split}D_{2}=&\left(0.863\times R_{1}\,mean)+(-1.044\times
R_{2}\,mean)+(-1.023\times R_{1}\,std)+(1.046\times
R_{2}\,std)\\!\\!\right.\\\ &\left.+\,(17.026\times R_{1}\,CV)+(-6.788\times
R_{2}\,CV)+(0.225\times R_{1}\,duration)\right.\\\ &\left.+\,(-0.203\times
R_{2}\,duration)+(0.017\times R_{1}\,area)+(-0.003\times
R_{2}\,area)-2.414\right.\end{split}$ (3)

Fig. 20 plots all samples of discriminant ‘Function 1’ and ‘Function 2’ which
cover 96.7% of variance. The group centroids are plotted from left to right as
the number of people increase. This is a good indication that the groups are
well discriminated by functions $D_{1}$ and $D_{2}$.

Figure 20: Combined Group Plots of Discriminant Functions 1 and 2

Finally, the classification is performed based on the discriminant score. The
classification results are shown in Table V. All correct prediction cases lie
on the diagonal. The classification results show 81.6% overall accuracy in
detecting the number of people comprising a given group. Further, an overall
accuracy of 97.9% was achieved in predicting individual head counts. For
example, the actual head count of 250 samples is 750 and the predicted head
count was 734.

| NPeople | Predicted Group Membership |
---|---|---|---
|  | 1 | 2 | 3 | 4 | 5 | Total
Count | 1 | 47 | 3 | 0 | 0 | 0 | 50
| 2 | 3 | 46 | 1 | 0 | 0 | 50
| 3 | 0 | 3 | 44 | 2 | 1 | 50
| 4 | 0 | 3 | 7 | 31 | 9 | 50
| 5 | 0 | 0 | 0 | 14 | 36 | 50
% | 1 | 94 | 6 | 0 | 0 | 0 | 100
| 2 | 6 | 92 | 2 | 0 | 0 | 100
| 3 | 0 | 6 | 88 | 4 | 2 | 100
| 4 | 0 | 12 | 14 | 62 | 18 | 100
| 5 | 0 | 0 | 0 | 28 | 72 | 100
TABLE V: Cassification Results

##  VI Multiple Object Detection under Uncontrolled Environment

In order to evaluate the accuracy in counting the number of people with no
specific formation using the discriminant function that is calculated in
previous section, we deployed the wireless sensor nodes in an uncontrolled
environment. We deployed a configuration of one transmitter and three
receivers, as shown in Fig. 21a where the darker areas represent corridors and
walkways. Three receivers were placed on one side of a wide corridor with one
transmitter on the opposite side. The distance of T-R2 is 3m which allows up
to five people walking past the corridor side-by-side. The distance of R1-R2
is 2.5m is the same as the distance of R2-R3. This sensor placement
configuration can be considered as two sets of one-transmitter two-receiver
configuration in symmetry. Each experiment was conducted for 15 minutes. The
total number of people walking past was recorded on video to compare for
verification of the detection result. A video frame of an experiment in
progress is shown in Fig.21b.

(a) Layout of Sensing Zone

(b) Experiment in progress with passerbys crossing

Figure 21: Test Environment Figure 22: Inferring Results under Uncontrolled
Environment

For convenience in data analysis and comparison, the data collection was
separated into 3-minute periods. Fig. 22 shows the results for one data
collection period, during which our proposed detection system has been able to
correctly infer every single occurrence of passerbys crossing the sensing
zone. The recorded video shows the sequence of people in groups walking past
the corridor as {2, 1, 1, 3, 1, 1, 1, 1, 3, 1}. The detection duration
difference between receivers increases for 1st, 4th and 9th positive
detections. As discussed previously, this is caused by the larger interference
zone when more than one person passes, as shown in Fig. 12b.

Using the one-transmitter three-receiver configuration, two sets of
statistical information of positive detections can be extracted since it
contains two sets of one-transmitter two-receiver configuration. When applying
new cases, the discriminant score will be calculated using the discriminant
functions that are constructed during training phase. All new samples’
discriminant scores of discriminant ‘Function 1’ and ‘Function 2’ for R1-R2
and R3-R2 are plotted in Fig. 23 and Fig. 24 respectively. A sample can be
classified when its discriminant score lies on the diagonal to the group
centroid which is shown as a filled marker.

We select a sample of the experiments for discussion, where a total of 39
people walked past the sensing zone in groups ranging from one to three
persons; as this was an uncontrolled environment, we had not been able to
detect groups of larger sizes although some of the instances could have been
classified as a larger group.

Figure 23: Combined Group Plots of Discriminant Functions 1 & 2 for R1-R2
Figure 24: Combined Group Plots of Discriminant Functions 1 & 2 for R3-R2

First, we present the classification results for R1-R2 which are summarized in
Table VI . From the R1-R2 dataset, the accuracy of predicting the number of
people comprising a given group is 78.6%, which is slightly lower than that
for the controlled environment. This is expected as a few people walking in
close proximity can be grouped in different ways from different angles and
their relative positions with the group can change dynamically. More
importantly, we aim to estimate the number of people in total, and that we
have been able to achieve with an accuracy of around 97.4% for this
experiment, and with similar accuracies for other experiments. This is very
close to the accuracy achieved in the controlled environment.

| NPeople | Predicted Group Membership |
---|---|---|---
|  | 1 | 2 | 3 | 4 | 5 | Total
Count | 1 | 17 | 1 | 2 | 0 | 0 | 20
| 2 | 2 | 3 | 0 | 0 | 0 | 5
| 3 | 1 | 0 | 2 | 0 | 0 | 3
| 4 | 0 | 0 | 0 | 0 | 0 | 0
| 5 | 0 | 0 | 0 | 0 | 0 | 0
% | 1 | 85 | 5 | 10 | 0 | 0 | 100
| 2 | 40 | 60 | 0 | 0 | 0 | 100
| 3 | 33.3 | 0 | 66.7 | 0 | 0 | 100
| 4 | 0 | 0 | 0 | 0 | 0 | 0
| 5 | 0 | 0 | 0 | 0 | 0 | 0
TABLE VI: R1-R2 Prediction Results under Uncontrolled Environment

Next, we present the data from R3-R2, tabulated in Table VII. The accuracy of
group size estimation dropped to 60.7% for the same reasons as the R1-R2 case.
The accuracy for number of people detected is 71.8% with a noticeable degree
of over-counting. Upon careful analysis of the video, we identified a possible
cause of the over-counting that led to the higher estimation error. From Fig.
21a, we note that R3 is located next to a stairwell on the left. That is a
highly utilized staircase and, as a result, many people pass close to R3 but
not necessarily continue right along the corridor across the R3-R2 sensing
zone. This has resulted in two instances of over-estimating the group size by
more than one, where a group of two people has been predicted as a group of
four people. E.g., if a group four people come down the stairs and split into
two groups of two, with one group continuing right along the corridor across
the R3-R2 sensing zone while the other, moving in the opposite direction or
straight across to the walkway on the lower left of Fig. 21a, then there is a
higher probability of the group size prediction to be affected.

| NPeople | Predicted Group Membership |
---|---|---|---
|  | 1 | 2 | 3 | 4 | 5 | Total
Count | 1 | 14 | 4 | 2 | 0 | 0 | 20
| 2 | 1 | 2 | 0 | 2 | 0 | 5
| 3 | 0 | 1 | 1 | 1 | 0 | 3
| 4 | 0 | 0 | 0 | 0 | 0 | 0
| 5 | 0 | 0 | 0 | 0 | 0 | 0
% | 1 | 70 | 20 | 10 | 0 | 0 | 100
| 2 | 20 | 40 | 0 | 40 | 0 | 100
| 3 | 0 | 33.3 | 33.3 | 33.3 | 0 | 100
| 4 | 0 | 0 | 0 | 0 | 0 | 0
| 5 | 0 | 0 | 0 | 0 | 0 | 0
TABLE VII: R3-R2 Prediction Results under Uncontrolled Environment

Based on the assumption that we intend to utilize pre-deployed wireless
communication devices, the R3-R2 case highlights the need for careful
selection of transmitter-receiver combinations in order to minimize prediction
errors. For the R1-R2 case, although R1 is also near the building
entrance/exit on the right, its position is sufficiently far away such that
people moving in and out of the entrance/exit will not have a significant
effect, if at all, on its prediction.

##  VII Conclusion

The use of radio irregularity resulting from the movement of human objects
crossing the path of a radio signal to detect human presence has been
demonstrated previously and applied to intrusion detection [2, 3, 15].
However, the ability to detect more than one person remains a challenge if we
rely on the characteristics of one signal’s fluctuations. With pervasive
networking brought about by large cyber-physical systems like the Internet of
Things, the presence of numerous wireless communication devices allow us to
study the fluctuations of multiple signals in close proximity of one another
as a result of human interference and deduce the number of human objects that
have crossed the paths of these signals.

In this paper, we have proposed a network-oriented approach that utilizes
received signal strength information (RSSI) of received packets to detect and
count people when they cross the signal transmission paths. This information
can be easily obtained from device drivers of wireless network interfaces when
the packets are received and the goal of our approach is to be able to easily
utilize the existing wireless transmitters and receivers already deployed in
the environment. Our approach which is based on the RSSI fluctuations between
consecutive packets does not require accurate channel models nor complex
signal processing techniques. It only needs to be trained to detect the RSSI
fluctuation patterns associated with the objects of interest, e.g. groups of
different number of people; no additional tuning is required during
deployment.

Using a simple configuration of two receivers deployed in close proximity to
each other, we have first demonstrated the ability to detect two persons
walking side-by-side along a typical 1.5m wide corridor using a
straightforward approach based on the difference in the periods of
fluctuations experienced by the two signals paths as the two human subjects
pass. We then extended our scheme to detect more human subjects using the same
two-receiver configuration together with discriminant analysis to process the
signal fluctuation data, we validated our scheme in a controlled environment
and showed that it is able to accurately detect and count up to five persons
with an accuracy of almost 98%. Next, we deployed our scheme in a public area
without the ability to control the mobility patterns and group structure of
passerbys, and achieved comparable accuracy in counting people. However, we
also note that a poorly located receiver can induce high estimation errors and
significantly reduce the accuracy of the system. Ideally, a robust placement
strategy should be chosen [22]. However, in our target scenario, achieving the
best device placements from among the already deployed wireless devices may
not always be possible. Additional wireless devices need to be deployed to
complement the existing topology, if the critical voids in detection need to
be covered.

From this study, we aim to show that a large cyber-physical system, like the
Internet of Things, can be exploited for applications like people counting
without the need for specialized hardware. However, our method is not aimed to
completely replace the specialized hardware for automated people counting but
more as a complementary technology. While the scheme in its current form
requires further work to enhance its capabilities further, it presents an
exciting opportunity to turn an existing wireless communications network into
a sensing system for automated people counting.

As our ongoing and future work, we are extending the scheme for automated
people counting in outdoor environments, e.g. to count visitors in public
parks [14], crowd size estimation [23], etc. Such technologies are
increasingly being deployed for crowd size estimation to assist in crowd
control and prevent any potential problems arising from loss of control over
crowd size. Often, agencies involved in crowd safety and management require
quick estimates to assist personnel on the ground, but most of the available
technology rely on image and video processing which are complex and expensive.

## References

  * [1] O. Vermesan _et al._ , _Internet of Things Strategic Research Roadmap_. European Commission - Information Society and Media DG, 15 Sep 2009.
  * [2] P. Lee, W. K. G. Seah, H. P. Tan, and Z. Yao, “Wireless Sensing Without Sensors – An Experimental Study of Motion/Intrusion Detection using RF Irregularity,” _Measurement Science and Technology, Special Issue on Wireless Sensor Networks: designing for real-world deployment and deployment experiences_ , vol. 21, no. 12, Nov 2010.
  * [3] M. Youssef, M. Mah, and A. Agrawala, “Challenges: device-free passive localization for wireless environments,” in _Proceedings of the 13th Annual ACM International Conference on Mobile Computing and Networking_ , Montréal, QC, Canada, 9-14 Sep 2007, pp. 222–229. [Online]. Available: http://doi.acm.org/10.1145/1287853.1287880
  * [4] N. Patwari and J. Wilson, “RF Sensor Networks for Device-Free Localization: Measurements, Models, and Algorithms,” _Proceedings of the IEEE_ , vol. 98, no. 11, pp. 1961–1973, Nov 2010.
  * [5] E. Mathews and A. Poigne, “An echo state network based pedestrian counting system using wireless sensor networks,” in _Proceedings of the International Workshop on Intelligent Solutions in Embedded Systems_ , Regensburg, Germany, 10-11 Jul 2008, pp. 1–14.
  * [6] Y. Zheng, D. J. Brady, and P. K. Agarwal, “Localization using boundary sensors: An analysis based on graph theory,” _ACM Transactions on Sensor Networks_ , vol. 3, no. 4, Oct 2007. [Online]. Available: http://doi.acm.org/10.1145/1281492.1281496
  * [7] J. Nakamura, Y. Tomita, and S. Honda, “Algorithm for counting the number of persons using IR image processing,” in _Proceedings of the IEEE Instrumentation and Measurement Technology Conference (IMTC)_ , Hamamatsu, Japan, 10-12 May 1994, pp. 220–223.
  * [8] T. Teixeira and A. Savvides, “Lightweight People Counting and Localizing for Easily Deployable Indoors WSNs,” _IEEE Journal of Selected Topics in Signal Processing_ , vol. 2, no. 4, pp. 493–502, Aug 2008.
  * [9] H. Celik, A. Hanjalic, and E. Hendriks, “Towards a robust solution to people counting,” in _Proceedings of the IEEE International Conference on Image Processing (ICIP)_ , Atlanta, GA, USA, 8-11 Oct 2006, pp. 2401–2404.
  * [10] I. Amin, A. Taylor, F. Junejo, A. Al-Habaibeh, and R. Parkin, “Automated people-counting by using low-resolution infrared and visual cameras,” _Measurement_ , vol. 41, no. 6, pp. 589 – 599, 2008.
  * [11] K. Woyach, D. Puccinelli, and M. Haenggi, “Sensorless sensing in wireless networks: Implementation and measurements,” in _Proceedings of the 2nd International Workshop on Wireless Network Measurement (WiNMee)_ , Boston, MA, USA, Apr 2006.
  * [12] J. Wilson and N. Patwari, “Radio tomographic imaging with wireless networks,” _IEEE Transactions on Mobile Computing_ , vol. 9, pp. 621–632, 2010.
  * [13] N. Patwari and J. Wilson, “Spatial models for human motion-induced signal strength variance on static links,” _IEEE Transactions on Information Forensics and Security_ , vol. 6, no. 3, pp. 791–802, Sep 2011.
  * [14] D. Puccinelli, A. Foerster, A. Puiatti, and S. Giordano, “Radio-based trail usage monitoring with low-end motes,” in _Proceedings of the 7th IEEE International Workshop on Sensor Networks and Systems for Pervasive Computing (PerSeNS)_ , Seattle, WA, USA, 21 Mar 2011.
  * [15] W.-C. Lin, W. K. G. Seah, and W. Li, “Exploiting Radio Irregularity in the Internet of Things for Automated People Counting,” in _Proceedings of the IEEE 22nd International Symposium on Personal Indoor and Mobile Radio Communications (PIMRC)_ , Toronto, ON, Canada, 11-14 Sep 2011, pp. 1015–1019.
  * [16] W. C. Lin, W. Li, W. K. G. Seah, and T. A. Gulliver, “Inferring human activity using overcomplete dictionary based pattern recognition,” in _Proceedings of the IEEE Pacific Rim Conference on Communications, Computers and Signal Processing (PacRim)_ , Victoria, BC, Canada, 23-26 Aug 2011, pp. 155–159.
  * [17] M. Nakatsuka, H. Iwatani, and J. Katto, “A study on passive crowd density estimation using wireless sensors,” in _Proceedings of 4th International Conference on Mobile Computing and Ubiquitous Networking (ICMU)_ , Tokyo, Japan, 11-13 Jun 2008.
  * [18] G. J. McLachlan, _Discriminant Analysis and Statistical Pattern Recognition_. John Wiley & Sons, Inc., 2004.
  * [19] R. A. Fisher, “The Use of Multiple Measurements in Taxonomic Problems,” _Annals of Human Genetics_ , vol. 7, no. 2, pp. 179–188, 1936. [Online]. Available: http://dx.doi.org/10.1111/j.1469-1809.1936.tb02137.x
  * [20] M. Welling, “Fisher Linear Discriminant Analysis,” 2005, Class Notes, University of California, Irvine, CA, USA. [Online]. Available: www.ics.uci.edu/~welling/classnotes/papers_class/Fisher-LDA.pdf
  * [21] R. B. Burns and R. A. Burns, _Business Research Methods and Statistics Using SPSS_. SAGE publication Ltd., 2008\.
  * [22] A. Krause, C. Guestrin, A. Gupta, and J. Kleinberg, “Robust sensor placements at informative and communication-efficient locations,” _ACM Transactions on Sensor Networks_ , vol. 7, no. 4, pp. 31:1–31:33, Feb 2011. [Online]. Available: http://doi.acm.org/10.1145/1921621.1921625
  * [23] Y.-L. Hou and G. Pang, “Automated people counting at a mass site,” in _Proceedings of the IEEE International Conference on Automation and Logistics (ICAL)_ , Qingdao, China, 1-3 September 2008, pp. 464–469.