text:
# 1 Introduction

## 1 Introduction

The Sudbury Neutrino Observatory (SNO) collaboration has presented results from 241 days running of their solar neutrino experiment , a reminder that it is still early days yet in this fruitful area of astro-particle physics. From the data on charged current (CC) deuteron disintegration ($`\nu _e+dp+p+e^{}`$) they observe a $`\nu _e`$ flux which is $`3.3\sigma `$ less than the neutrino flux measured by the Super-Kamiokande collaboration (SK) via $`\nu e^{}`$ scattering (ES) 111SNO also has results from ES events  which agree with the SK results  but is of poorer statistics currently.. This further strengthens the oscillation interpretation of the solar neutrino problem and, at least in a two neutrino scenario, disfavours the mixing of the $`\nu _e`$ with a sterile species.

Violation of the Equivalence Principle (VEP) has been examined in the literature as a possible solution to the solar neutrino problem , one of the alternatives to the mass-mixing vacuum oscillation and the matter induced Mikheyev-Smirnov-Wolfenstein (MSW) approaches . In this note we revisit the VEP solution in a two-flavor picture with the $`\nu _e`$ oscillating to another active species both of which are assumed to be massless.

Neutrino oscillation in VEP can occur if (a) the weak equivalence principle is not satisfied, i.e., the coupling of neutrinos to the gravitational field is nonuniversal and (b) the flavor eigenstates are not identical to the states that couple to gravity . It does not require neutrinos to carry a non-zero mass. The weak equivalence principle demands the coupling of particles to the ambient gravitational potential $`\varphi `$ to be of the form $`f\varphi E`$, where $`E`$ is the particle energy, and $`f`$ a universal coupling constant. If the latter varies from one neutrino species to another then that would constitute a violation of the equivalence principle. If $`f_1f_2`$ in a two-neutrino framework, then these states define a basis in the two-dimensional space which, in general, could be different from the flavor basis.

$$\nu _e=\nu _1\mathrm{cos}\theta +\nu _2\mathrm{sin}\theta ;\nu _x=\nu _1\mathrm{sin}\theta +\nu _2\mathrm{cos}\theta ,$$
(1)
where $`x=\mu `$, $`\tau `$. The effect due to a small splitting $`\mathrm{\Delta }f`$ will manifest itself in the form of flavor oscillations, the wavelength going to infinity as $`\mathrm{\Delta }f`$ tends to zero. The essential difference between this approach and the mass-mixing solution appears in the energy dependence of the survival probability. For a two-neutrino picture, the general expression for the survival probability for an initial $`\nu _e`$ after propagation through a distance $`L`$ is given by:

$$P_{ee}(E_\nu ,L)=1\mathrm{sin}^22\theta \mathrm{sin}^2\left(\frac{\pi L}{\lambda }\right),$$
(2)
where, for the VEP induced case, the oscillation length $`\lambda `$ is:

$$\lambda =\frac{2\pi }{E_\nu \varphi \mathrm{\Delta }f}.$$
(3)
$`E_\nu `$ is the neutrino energy. Thus, unlike the mass-mixing case where $`\lambda E_\nu `$, for VEP $`\lambda 1/E_\nu `$. Due to the different energy dependences of the survival probability in the mass-mixing and the VEP alternatives, their predictions can be quite different for solar neutrinos of different energies. We follow the prevalent practice of choosing the gravitational potential, $`\varphi `$, to be a constant over the neutrino path. This is the case if the potential due to the Great Attractor  dominates over that due to the sun and other heavenly bodies in our neighborhood. In such an event, writing $`\mathrm{\Delta }F=\varphi \mathrm{\Delta }f/2`$, the expression for the oscillation wavelength, Eq. (3), becomes

$$\lambda =\frac{6.20\times 10^{13}\mathrm{m}}{\mathrm{\Delta }F}\left(\frac{1\mathrm{M}\mathrm{e}\mathrm{V}}{E_\nu }\right).$$
(4)
An alternate physics scenario which also leads to neutrino oscillations with $`\lambda `$ $`1/E_\nu `$ is a picture of violation of special relativity (VSR) . If special relativity is violated, the maximum attainable speed of a particle in vacuo need not universally be the speed of light $`c`$. In particular, if the maximum possible velocities of the two types of neutrinos be $`v_1`$ and $`v_2`$ and these velocity eigenstate neutrinos be related to the $`\nu _e`$ and $`\nu _x`$ through a mixing angle $`\theta `$ (see Eq.(1)) then in this case the expression for $`\lambda `$ is:

$$\lambda =\frac{2\pi }{E_\nu \mathrm{\Delta }v},$$
(5)
where $`\mathrm{\Delta }v`$ is the velocity difference for the neutrinos $`\nu _1`$ and $`\nu _2`$. Comparing Eqs. (5) and (3) one finds that the energy dependence of the oscillation length is identical in the two cases222It has been shown that inclusion of CPT-violating interactions in addition to Violation of Special Relativity can lead to more general energy dependences involving $`1/E_\nu `$, $`E_\nu `$, and constant terms  and the role of $`\mathrm{\Delta }v`$ in the VSR case is the same as that of $`\varphi \mathrm{\Delta }f`$ in the VEP formalism. Here, we use the terminology of the VEP mechanism but the results can be taken over mutatis mutandis to the VSR situation.

In this work we use the BBP2000 calculation  of the solar neutrino flux as the SSM reference. In addition to this, we explore the possibility of the absolute normalizations of the solar 8B- and $`hep`$-neutrino spectra333The SNO and SK experiments are sensitive only to the 8B and $`hep`$ neutrinos. being different from their BBP2000 SSM predictions. If $`X_B`$ and $`X_{hep}`$ denote the factors by which the absolute normalizations are multiplied, we use the data to find the best-fit values for these. We find that all the fits are improved in a noteworthy manner when $`X_B`$ is permitted to be less than unity.

## 2 Update of VEP for the SK ES results

In contrast to other recent work on VEP  based on the 825-day SK spectral data, in a previous publication , the SK ES results from 1117 days of running  were considered. It was found that the predictions were by-and-large robust. One trend, which was remarked upon, was that the softening of the highest energy SK ES data had resulted in a reduction of $`X_{hep}`$. In fact, it was noted that the best-fit to the day-night spectrum of the 1117-day SK data preferred a negligibly small value of $`X_{hep}`$.

For a further update, we have performed a similar analysis with the SK 1258-day data. The SK results are presented in the form of number of events (normalized to the SSM expectation) in 19 electron recoil energy bins of width 0.5 MeV in the range 5 MeV to 14 MeV and a 19th bin which covers the events in the range 14 to 20 MeV .

The definition of $`\chi ^2`$, the error correlations, etc. are chosen in the same manner as in . Suffice it to say that we have included the statistical error, the uncorrelated systematic errors, and the energy-bin-correlated experimental errors  as well as those from the calculation of the shape of the expected spectrum . When we allow the normalizations of the 8B and $`hep`$ fluxes to vary, we do not include their astrophysical uncertainties separately. The results are presented in Table 1.

| Fitted | Case | $`\mathrm{sin}^22\theta `$ | $`\mathrm{\Delta }F`$ | $`X_B`$ | $`X_{hep}`$ | $`\chi _{min}^2`$/d.o.f. | g.o.f |
| --- | --- | --- | --- | --- | --- | --- | --- |
| Data |  |  | $`(10^{24})`$ |  |  |  | (%) |
| SK1258 ES | 1a | 0.72 | 0.22 | 1.0 | 1.0 | 29.10/17 | 3.37 |
| Day-Night averaged |  |  |  | (fixed) | (fixed) |  |  |
| Spectrum | 1b | 1.0 | 0.99 | 0.77 | $`2.4\times 10^5`$ | 11.38/15 | 72.52 |
| SK1258 ES | 1c | 0.72 | 0.21 | 1.0 | 1.0 | 43.71/36 | 17.67 |
| separate Day & Night |  |  |  | (fixed) | (fixed) |  |  |
| Spectra | 1d | 1.0 | 0.98 | 0.77 | $`9.3\times 10^8`$ | 24.68/34 | 87.92 |

Table 1: The best-fit values of the parameters, $`\mathrm{sin}^22\theta `$, $`\mathrm{\Delta }F`$, $`X_B`$, $`X_{hep}`$, the $`\chi _{min}^2`$, and the goodness of fit (g.o.f.) for fits to the SK1258 ES spectra.

We find that the fits are by-and-large robust; using the separate day and night data or the day-night averaged data affects them only marginally. Further, the best-fit parameters have not changed much and they are well within the 90% C.L. allowed regions of the previous analysis using the 1117-day data. For example, for the fit to the the day-night averaged data sample, in the SSM case ($`X_B=X_{hep}=1`$), the best-fit values of $`\mathrm{sin}^22\theta `$ and $`\mathrm{\Delta }F`$ have changed $`<`$ (5 - 10)% but the quality of fit has significantly dropped (g.o.f. = 76.7% for the fit to the 1117-day data). When $`X_B`$ and $`X_{hep}`$ are allowed to vary, the best-fit value of (a) $`X_{hep}`$ is tiny, and (b) $`\mathrm{\Delta }F`$ is larger than earlier and closer to the value for the SSM fit, a welcome feature. The goodness of fit is also a bit lower. For the fit to the separate day and night data, the differences between the fits to the 1258-day and 1117-day data broadly follow the same pattern as above.

## 3 VEP Fits to the SNO CC and SK ES spectra

The SNO collaboration has published the first measurement of the spectrum of electron kinetic energy produced from CC deuteron disintegration by solar neutrinos. The results are presented in eleven electron kinetic energy bins ranging from 6.5 MeV to 13.0 MeV. This provides a new check on solutions to the solar neutrino puzzle. The best-fit values of the VEP parameters obtained from fitting this data sample are shown in Table 2.

When the SNO CC and SK 1258-day ES data (separate day and night or day-night averaged) are simultaneously addressed within the VEP picture, the best-fit values are remarkably close to each other in the SSM fits (2a), (2c), and (2e) though the goodness of fit is low. When $`X_B`$ and $`X_{hep}`$ are allowed to vary, then the fits to the combined SK ES and SNO CC spectra improve. These results are also presented in Table 2.

| Fitted | Case | $`\mathrm{sin}^22\theta `$ | $`\mathrm{\Delta }F`$ | $`X_B`$ | $`X_{hep}`$ | $`\chi _{min}^2`$/d.o.f. | g.o.f |
| --- | --- | --- | --- | --- | --- | --- | --- |
| Data |  |  | $`(10^{24})`$ |  |  |  | (%) |
|  | 2a | 0.70 | 0.20 | 1.0 | 1.0 | 9.74/9 | 37.19 |
| SNO CC |  |  |  | (fixed) | (fixed) |  |  |
| Spectrum | 2b | 0.72 | 0.20 | 1.09 | $`6.6\times 10^7`$ | 9.70/7 | 20.62 |
| alone |  |  |  |  |  |  |  |
|  | 2c | 0.72 | 0.21 | 1.0 | 1.0 | 42.04/28 | 4.30 |
| SNO CC & SK1258 |  |  |  | (fixed) | (fixed) |  |  |
| (Day-Night averaged) | 2d | 0.66 | 0.99 | 0.62 | $`1.9\times 10^{10}`$ | 24.60/26 | 54.17 |
| Spectra |  |  |  |  |  |  |  |
|  | 2e | 0.72 | 0.21 | 1.0 | 1.0 | 56.60/47 | 15.93 |
| SNO CC & SK1258 |  |  |  | (fixed) | (fixed) |  |  |
| (Separate Day & Night) | 2f | 0.66 | 0.99 | 0.61 | $`9.3\times 10^{10}`$ | 37.79/45 | 76.84 |
| Spectra |  |  |  |  |  |  |  |

Table 2: The best-fit values of the parameters, $`\mathrm{sin}^22\theta `$, $`\mathrm{\Delta }F`$, $`X_B`$, $`X_{hep}`$, the $`\chi _{min}^2`$, and the g.o.f. for fits to SNO CC and SK1258 ES spectra.

As seen in the fit to the SK ES spectra alone (Table 1), here again the fits do not change whether the separate day and night data or the day-night averaged data are used. In view of this, in the next section, we use only the day-night averaged data.

We show in Fig. 1, the allowed regions for the parameters $`\mathrm{sin}^22\theta `$ and $`\mathrm{\Delta }F`$ at 90% C.L. obtained by fitting the SK 1258-day ES spectrum (dot-dashed line), the SNO CC spectrum (small-dashed line), and for a combined fit to the SNO CC and the SK ES spectra (solid – for the SSM case – and dotted – when $`X_B`$ and $`X_{hep}`$ are allowed to vary). It is noteworthy that for the SSM, the allowed regions from the different fits overlap very nicely.

## 4 VEP Fits to total rates and the SNO CC and SK ES spectra

The SK ES and SNO CC spectrum measurements are excellent filters for choosing between alternative solutions of the solar neutrino problem. But it has to be borne in mind that they are both sensitive only to the highest energy 8B neutrinos. Further constraints on any putative solution come from the total solar neutrino rates measured by the Chlorine and Gallium experiments. In order to be acceptable, the VEP solution must also be confronted with these results, to which we now turn. The data used in the $`\chi ^2`$ analysis of total rates are given in Table 3.

| Experiment | Chlorine | Gallium | SK ES | SNO CC |
| --- | --- | --- | --- | --- |
| $`\frac{\mathrm{Observed}\mathrm{Rate}}{\mathrm{SSM}\mathrm{Prediction}}`$ | $`0.337\pm 0.030`$ | $`0.584\pm 0.039`$ | $`0.4594\pm 0.016`$ | $`0.347\pm 0.029`$ |

Table 3: The ratio of the observed solar neutrino rates to the corresponding BBP2000 SSM predictions used in this analysis. The results are from Refs. .The Gallium rate is the weighted average of the Gallex, SAGE, and GNO results.

The best-fit parameters obtained by fitting the total rates are presented in Table 4 (sets 4a and 4b). Here the goodness of fits are not high. In order to ascertain how the total rate measurements mesh with the SK ES and SNO CC measurements, we have performed combined fits to the rates and spectral data. These results are also shown in Table 4. Note that the SSM is disfavoured by this analysis.

| Fitted | Case | $`\mathrm{sin}^22\theta `$ | $`\mathrm{\Delta }F`$ | $`X_B`$ | $`X_{hep}`$ | $`\chi _{min}^2`$/d.o.f. | g.o.f |
| --- | --- | --- | --- | --- | --- | --- | --- |
| Data |  |  | $`(10^{24})`$ |  |  |  | (%) |
|  | 4a | 1.0 | 1.66 | 1.0 | 1.0 | 5.50/2 | 6.39 |
| Total Rates |  |  |  | (fixed) | (fixed) |  |  |
| alone | 4b | 1.0 | 1.66 | 0.74 | 1.0 | 2.40/1 | 12.13 |
|  |  |  |  |  | (fixed) |  |  |
| SNO CC, SK1258 | 4c | 0.72 | 0.21 | 1.0 | 1.0 | 117.62/32 | $`1.0\times 10^9`$ |
| ES spectra & Total |  |  |  | (fixed) | (fixed) |  |  |
| Rates (incl. | 4d | 1.0 | 1.60 | 0.77 | $`6.8\times 10^8`$ | 33.54/30 | 29.96 |
| SNO & SK) |  |  |  |  |  |  |  |
| SNO CC, SK1258 | 4e | 0.72 | 0.21 | 1.0 | 1.0 | 92.43/30 | $`2.8\times 10^6`$ |
| ES spectra & Total |  |  |  | (fixed) | (fixed) |  |  |
| Rates (excl. | 4f | 1.0 | 1.62 | 0.77 | $`1.8\times 10^5`$ | 31.04/28 | 31.53 |
| SNO & SK) |  |  |  |  |  |  |  |

Table 4: The best-fit values of the parameters, $`\mathrm{sin}^22\theta `$, $`\mathrm{\Delta }F`$, $`X_B`$, $`X_{hep}`$, the $`\chi _{min}^2`$, and the g.o.f. for fits to the total rates along with the SNO CC and SK ES spectra.

When we fit the total rates and the SNO CC and SK ES spectra together, the SNO and SK total rates and the corresponding spectra cannot be considered to be statistically independent pieces of data. However, their exclusion from the fits would also not be entirely correct as the fitted spectral ratios (observed/SSM) do not completely determine the total rates. We have therefore chosen to examine two extreme alternatives; i.e., in the fits to the spectra and the total rates together, (a) in one case we include the SK and SNO total rates in the fit, and (b) in the other case we have excluded them. It is gratifying that the best-fit values of the parameters of the two cases (4c, 4e) or (4d, 4f) are not very different, a reflection of the comparative energy independence of the (observed/SSM) spectra.

The best-fit values of the combined total rates + spectra fits for the SSM case (4c, 4e) are within the 90% C.L. allowed regions obtained from fits to the spectra alone (see Fig. 1).

## 5 SNO NC expectation, Comparison of different fits

SNO will soon start taking data for a calorimetric measurement of neutral current (NC) deuteron disintegration ($`\nu +d\nu +p+n`$) to which all active neutrinos contribute equally. A measure of neutrino oscillations is provided by the ratio of the NC and CC rates, $`R_{NC}`$ and $`R_{CC}`$, which is somewhat less sensitive to theoretical uncertainties than the rates themselves. We define

$$R_{SNO}=\frac{R_{NC}}{R_{CC}}.$$
(6)
For the best fit values of parameters obtained in the earlier sections, we present the predicted values of $`R_{SNO}`$ in Table 5. For an assessment of the various fits, for every set of best-fit parameters we also show in Table 5 the following quantities: the calculated rates for the Chlorine, Gallium, SK ES and SNO CC observations, and also the $`\chi ^2`$ obtained for the calculated SK ES and SNO CC electron spectra.

| Case | Cl | Ga | SK | | SNO | | |
| --- | --- | --- | --- | --- | --- | --- | --- |
|  |  |  | Rate | $`\chi ^2`$ | Rate | $`\chi ^2`$ | $`R_{SNO}`$ |
|  | ($`0.337\pm 0.030`$) | ($`0.584\pm 0.039`$) | ($`0.4594\pm 0.016`$) | (spect) | ($`0.347\pm 0.029`$) | (spect) |  |
| 1a | 0.516 | 0.940 | 0.479 | 29.10 | 0.353 | 11.65 | 1.00 |
| 1b | 0.437 | 0.796 | 0.462 | 11.38 | 0.380 | 20.96 | 0.71 |
| 1c | 0.508 | 0.940 | 0.476 | 29.40 | 0.341 | 10.73 | 1.04 |
| 1d | 0.436 | 0.798 | 0.458 | 11.50 | 0.380 | 21.62 | 0.71 |
| 2a | 0.522 | 0.941 | 0.502 | 45.50 | 0.356 | 9.74 | 0.99 |
| 2b | 0.523 | 0.942 | 0.522 | 78.60 | 0.356 | 9.70 | 1.07 |
| 2c | 0.510 | 0.940 | 0.481 | 29.96 | 0.343 | 10.35 | 1.03 |
| 2d | 0.496 | 0.858 | 0.461 | 13.14 | 0.418 | 12.53 | 0.51 |
| 2e | 0.506 | 0.940 | 0.479 | 30.25 | 0.336 | 10.15 | 1.05 |
| 2f | 0.495 | 0.859 | 0.459 | 13.20 | 0.416 | 12.56 | 0.51 |
| 4a | 0.443 | 0.609 | 0.605 | 355.6 | 0.520 | 23.18 | 0.67 |
| 4b | 0.339 | 0.596 | 0.448 | 23.60 | 0.385 | 12.37 | 0.67 |
| 4c | 0.511 | 0.940 | 0.482 | 30.26 | 0.344 | 10.29 | 1.03 |
| 4d | 0.352 | 0.613 | 0.456 | 18.40 | 0.400 | 12.94 | 0.67 |
| 4e | 0.511 | 0.940 | 0.480 | 29.75 | 0.344 | 10.44 | 1.03 |
| 4f | 0.350 | 0.608 | 0.456 | 18.40 | 0.400 | 12.97 | 0.67 |

Table 5: The total rates for the different experiments obtained by using the best-fit values of the VEP parameters presented in Tables 1, 2, and 4. Also presented are the $`\chi ^2`$ for fits to the observed SK ES and SNO CC spectra using these values of the parameters. The predictions for $`R_{SNO}`$ are also shown.

From Table 5 it can be seen that fits to the spectra by themselves give best-fit parameters grossly disfavoured by the total rates measured in the Chlorine and Gallium experiments. The most acceptable results are obtained in cases (4d, 4f) – combined fits to the SK ES spectra, SNO CC spectra, and the Chlorine and Gallium total rates measurement with $`X_B`$ = 0.77 (see Table 4). These values of parameters will be further strengthened if $`R_{SNO}`$ is measured around 0.67.

## 6 Conclusions and Discussions

In this work, we have considered the VEP oscillation explanation of the solar neutrino problem within a two-flavor (active) scenario in the light of the SNO CC measurement, the SK 1258-day ES data, and the total rates from the Chlorine and Gallium experiments.

We find that the VEP fits to the SK ES data are robust. The best-fit parameters have changed only marginally from the earlier fit to the 1117-day data.

In the literature there are many variants of the data-fitting procedure for solar neutrinos. The solution preferred by Nature should fit all the solar neutrino results collectively and also individually. A good fit to one (or some) of the results may do rather poorly when confronted with another piece of datum. A good fit to all the data taken together may mask bad fits to some of the results taken in isolation. This caveat is overlooked in much of the current literature. We have attempted to fit the SK ES and SNO CC spectra individually and jointly, and also the total rates by themselves and jointly with the spectral data in an attempt to find a region of parameter space that is acceptable in every respect. These results are given in Table 5 and in Tables 1, 2, and 4. We find that the best-fit parameters of the cases (4d) and (4f) (see Tables 4 & 5) best meet the requirements alluded to above. Note that these fits correspond to $`X_B0.77`$ and predict $`R_{SNO}`$ 0.67.

During the passage of the neutrinos from their point of production to the solar surface, interactions with the ambient matter, responsible for the MSW effect, become important. In the presence of VEP and this MSW contribution, the effective neutrino mass matrix in flavor space takes the form

$$M=\frac{1}{2}\left|\begin{array}{cc}E_\nu \mathrm{\Delta }F\mathrm{cos}2\theta 2\sqrt{2}G_Fn_e(r)& E_\nu \mathrm{\Delta }F\mathrm{sin}2\theta \\ E_\nu \mathrm{\Delta }F\mathrm{sin}2\theta & E_\nu \mathrm{\Delta }F\mathrm{cos}2\theta \end{array}\right|,$$
(7)
where we have dropped an irrelevant part proportional to the identity matrix. The MSW contribution in (7) inside the sun turns out to be several orders of magnitude larger than the terms due to VEP that we have discussed in this work. Recall that we have assumed the neutrinos to be massless. For maximal mixing (e.g., the preferred 4d and 4f cases), there is no resonance effect and, in fact, till such time that the neutrino emerges from the sun, the MSW contribution controls the masses in (7). Inside the sun, the $`\nu _e`$ is, therefore, a mass eigenstate to a very good approximation. The effect of VEP oscillations begins to manifest itself only from then onwards.

The VEP oscillation wavelengths favored by the data are comparable to the distance of the sun to the earth. Therefore, seasonal effects are expected in this picture. When more seasonal data from SK accumulates, it will provide a good check on this scenario.

Since the publication of the first results by the SNO collaboration, a number of analyses of the data in the mass-mixing vacuum oscillation and MSW approaches have appeared . The quality of these fits are comparable to those found in the VEP picture. Here too, the parameters which best fit the rates alone are somewhat different from those that best fit the rates together with the spectra . Further sharpening of the data from new and ongoing experiments, it is hoped, will help to distinguish between these alternatives.

Acknowledgements

The authors are grateful to the Calcutta University Computer Centre for the use of their Origin 2000 computer. They thank Debasish Majumdar for his help. The research of A.R. is supported in part by a grant from CSIR, India while A.S. enjoys a fellowship from UGC, India.


format: text

text:
# References

On the Erraticity in Random-Cascading $`\alpha `$ Model111Supported in part by the NSFC under Grant No 19975021

ZHOU Yi-Fei LIU Qin TANG Ying CHENG Chun222The authors are students of the State-Level Personnel Training Base for Research and Teaching in Fundamental Sciences (Physics).

(Institute of Particle Physics, Huazhong Normal University Wuhan 430079)

Abstract

The erraticity in the random-cascading $`\alpha `$ model is revisited. It is found that in contrary to the previous expectation, even in the pure single-$`\alpha `$ random-cascading model without putting in any particle there exists erraticity behavior and the corresponding entropy indices do not vanish. This means that the dynamical fluctuations in a pure single-$`\alpha `$ model already fluctuate event-by-event. Models with multiple-$`\alpha `$ strengthen this fluctuation. Taking double-$`\alpha `$ model as example, the variation of the event-space fluctuation strength with the mixing ratio of the two $`\alpha `$’s is studied in some detail. The influence of particle number on the results when particles are putted into the system is also investigated.

PACS: 13.85Hd

Keywords: Dynamical fluctuation random cascading $`\alpha `$ model erraticity

entropy index

Since the observation of large local fluctuations in a high multiplicity nucleus-nucleus collision event induced by cosmic ray in 1983, similar phenomenon has also been observed in the accelerator experiments. Since then, extensive attentions are paid to the non-linear phenomena in high energy collisions.

Due to the finiteness of particle number, the particle distribution fluctuates with respect to the probability distribution obtained from theoretical calculation. This is the so called statistical fluctuations. For a long time the popular idea was that the probability distribution obtained from theoretical calculation is smooth and the actual particle distribution fluctuates statistically on top of the smooth probability distribution. However, is it really the case? Is it possible that there are non-linear dynamics which causes the dynamical probability distribution itself to be non-smooth? This question is still to be answered. Although a systematic quantum non-linear theory has not been constructed yet and a definite reply to this quesition is still impossible by now, non-linear phenomena have already been observed in high energy experiments\[4–6\] as stated above.

In studying the non-linear dynamical fluctuations, besides the factorial moments averaged over event sample, the importance of the fluctuations of single-event moments in the event samples has also been noticed. It is shown that these fluctuations are associated with the system’s chaotic behavior. A parameter $`\mu `$, called entropy index, is introduced to describe the chaotic behavior. An entropy index greater than zero is the signal of chaos. This method is called “erraticity analysis”.

The random-cascading $`\alpha `$ model is a simple model which can simulate the non-linear dynamical fluctuations. It has been used in Ref.\[8–11\] to study the dynamical fluctuations in event space, i.e. erraticity. Since the value of parameter $`\alpha `$ determines the strength of dynamical fluctuations, it was argued in Ref. that there is no erraticity in the model with a single parameter $`\alpha `$, and therefore a model with Gaussian-distributed $`\alpha `$ was proposed instead.

In this letter the erraticity of $`\alpha `$ model is revisited. It is found that there is erraticity already in a pure single-$`\alpha `$ model. The corresponding entropy index is small but does not vanish. It increases with the increasing of $`\alpha `$. The mixing of different $`\alpha `$’s makes the entropy index increases further. The dependence of entropy index on the mixing ratio is studied. It is shown through putting in particles that in all cases the entropy index calculated by particle number is bigger than that in pure $`\alpha `$ model and tends to the latter only when the particle number tends to infinity. This reflects the effect of statistical fluctuations.

Let us first briefly recall the definition of erraticity and entropy index.

The single-event factorial moments $`F_q^{(e)}`$ and probabilitic moments $`C_q^{(e)}`$ are defined as

$$F_q^{(e)}=\frac{\frac{1}{M}\underset{m=1}{\overset{M}{}}n_m(n_m1)\mathrm{}(n_mq+1)}{\left(\frac{1}{M}\underset{m=1}{\overset{M}{}}n_m\right)^q},$$
(1)
$$C_q^{(e)}=M^{q1}\underset{m=1}{\overset{M}{}}\left(p_m^{(e)}\right)^q,$$
(2)
respectively, where a phase space region is divided into $`M`$ equivalent sub-regions, $`p_m`$ and $`n_m`$ are correspondingly the probability and particle numbers in the $`m`$th sub-region.

The fluctuation of moments $`F_q^{(e)}`$ and $`C_q^{(e)}`$ in the event space can be characterized by their event-space moments

$$C_{p,q}^{(F)}=F_q^{(e)p}/F_q^{(e)}^p,C_{p,q}^{(C)}=C_q^{(e)p}/C_q^{(e)}^p.$$
(3)
If $`C_{p,q}^{(F)}`$ or $`C_{p,q}^{(C)}`$ has power law behavior

$$C_{p,q}^{(F)}M^{\psi _q^{(F)}(p)},C_{p,q}^{(C)}M^{\psi _q^{(C)}(p)},$$
(4)
when $`M`$ is big then we say that there is erraticity in the system. It can be characterized quantitatively by the entropy index

$$\mu _q=\frac{\mathrm{d}}{\mathrm{dp}}\psi _q(p)|_{p=1}.$$
(5)
Another method for calculating entropy index is to express $`C_{p,q}`$ as

$$C_{p,q}=\mathrm{\Phi }_q^{(e)p},$$
(6)
where

$$\mathrm{\Phi }_q^{(e)}=F_q^{(e)}/F_q^{(e)},\mathrm{\Phi }_q^{(e)}=C_q^{(e)}/C_q^{(e)}.$$
(7)
Define

$$\mathrm{\Sigma }_q=\mathrm{\Phi }_q^{(e)}\mathrm{ln}\mathrm{\Phi }_q^{(e)},$$
(8)
then in the scaling region, i.e. in the region where $`\mathrm{\Sigma }_q`$ depends on $`\mathrm{ln}M`$ linearly, we can write

$$\mu _q=\frac{\mathrm{\Sigma }_q}{\mathrm{ln}M}.$$
(9)
In the following we will study the erraticity in the random-cascading $`\alpha `$ model. In this model the $`M`$ partition of the phase space region $`\mathrm{\Delta }`$ is realized in $`\nu `$ steps. At the first step, it is divided into two equal parts; at the second step, each part is further divided into two parts, …, and so on. The steps are repeated until $`M=\mathrm{\Delta }Y/\delta y=2^\nu `$. In this step-by-step partition how particles are distributed between the two parts of a given phase space cell is determined by the value of independent random variable $`\omega _{\nu j_\nu }`$, where $`j_\nu `$is the position of the sub-cell ($`(1j_\nu 2^\nu )`$), $`\nu `$ is the number of steps. The value of the random variable $`\omega _{\nu j_\nu }`$ is given by:

$$\omega _{\nu ,2j1}=\frac{1+\alpha r}{2},\omega _{\nu ,2j}=\frac{1\alpha r}{2},j=1,\mathrm{},2^{\nu 1}$$
(10)
where $`r`$ is a random number distributed uniformly in the interval \[-1,1\], $`\alpha `$ is a positive number less than or equal to unity, which determines the value-region of the random variable $`\omega `$, and thus describes the strength of dynamical fluctuations in the model. After $`\nu `$ steps, the probability in the $`m`$th sub-cell ($`m=1,\mathrm{},M`$) is

$$p_m=\omega _{1j_1}\omega _{2j_2}\mathrm{}\omega _{\nu j_\nu }.$$
(11)
Using these probabilities the probability moment $`C_q^{(e)}`$ in each event is calculated according to (2), and the entropy index $`\mu _q`$ is obtained through Eq’s.(2)-(9).

Besides the single-$`\alpha `$ model described above, the double-$`\alpha `$-mixing model will also be used later. It is constructed by letting $`r_1𝒩`$ events among the total $`𝒩`$ events have $`\alpha =\alpha _1`$, and the remaining $`(1r_1)𝒩`$ events have $`\alpha =\alpha _2`$.

In order to study the influence of limited particle number, we put a certain number $`N`$ of particles into the $`M`$ subcells according to the Bernauli distribution

$$B(n_1,\mathrm{},n_M|p_1,\mathrm{},p_M)=\frac{N!}{n_1!\mathrm{}n_M!}p_1^{n_1}\mathrm{}p_M^{n_M}.$$
(12)
The particle numbers $`n_m(m=1,\mathrm{},M)`$ in every window are thus obtained and the event factorial moments $`F_q^{(\mathrm{e})}`$ can be calculated. The entropy index $`\mu _q`$ is then obtained through Eq’s. (2)–(9).

Eq.’s (5) and (9) are two equivalent methods for getting entropy index. We now discuss which one is more suitable for our purpose.

What we are interested in is the case when the partition number $`M`$ of phase space is very big. When the number of particles is small and that of sub-cells is big, it is probable that the particle number in every window is 0 or 1. It can be seen from Eq.(1) that this will result in $`F_q^{(e)}=0`$ for all $`q>1`$, and then lead to $`\mathrm{\Phi }_q^{(e)}`$ and $`_q^{(e)}`$ both vanish. Therefore, the second method fails in this case. So we choose the first method to get $`\mu _q`$.

When using the first method, it can be seen from Eq.(5) that $`\mu _q`$ is obtained through partially differentiating $`\psi _q^{(\mathrm{e})}`$ with respect to $`p`$ near $`p=1`$. Then how to choose the range

of $`p`$ value near $`p=1`$ to fit the partial derivative is also a question which we must answer. A too small range will lead to a wrong result due to the limitation of computer precision. To make the choise we show in Fig.1 the $`\psi `$-$`p`$ plot. It can be seen from the figure that within the range $`0.5p2`$ $`\psi `$ versus $`p`$ is apprximately linear. Therefore, we take it as the right range for fitting the partial derivative $`\psi /p`$.

Fig.1 $`\psi p`$ relation

In real calculation when $`\alpha `$ is very small (e.g. $`\alpha `$=0.1) the value of single event probability moment $`C_q^{(e)}`$ is close to zero. This results in big errors in numerical calculation. In order to solve this problem, we multiply the value of single event probability moments by a large number, e.g. 1000. This extra factor will appear in both the denominator and the numerator in the formula for calculating $`C_{p,q}`$, Eq. (3), and will not effect the final results.

The calculating results of single $`\alpha `$ model show that there exists erraticity. As example, in Fig.2 is shown the variation of $`C_{p,q}`$ with $`M`$ for $`\alpha `$=0.3. It shows a typical characteristic of erraticity. The entropy indices $`\mu _2`$ in this model for various values of $`\alpha `$ are plotted at the right side of Fig.3. It can be seen that $`\mu _2>0`$, indicating the existence of erraticity. The value of entropy index $`\mu _2`$ increases with the increasing of $`\alpha `$.

Fig.2 The variation of $`C_{p,q}M`$ Fig.3 Entropy index for single

for $`\alpha =0.3`$ $`\alpha `$ model

On the left of Fig. 3 is shown the variation of entropy index with the number of particles when particles are putted into the model. It can be seen that when the number of particles is small the corresponding entropy index is much bigger than that of the pure $`\alpha `$ model without particle. The entropy index decreases with the increasing of the number of particles, tends to the value of the pure $`\alpha `$ model when the number of particles tends to infinity.

The reason why erraticity exists also in the single-$`\alpha `$ model is that, in the expression Eq. (10) for the elementary probability in the model, there is a random number $`r`$, causing the window probability to fluctuate inspite of the fact that the values of $`\alpha `$ are the same in different events.

It is clear that if the $`\alpha `$’s are different in different events, the fluctuation will be strengthened. We take the double-$`\alpha `$ model as example to study this situation.

The two $`\alpha `$’s are taken to be $`\alpha _1=0.2,\alpha _2=0.8`$. There are in total $`𝒩`$ events amoung which $`r_1𝒩`$ events have $`\alpha _1=0.2`$ and the remaining $`(1r_1)𝒩`$ events have $`\alpha _1=0.8`$.

In Fig.4 are shown the log$`{}_{10}{}^{}\mu `$ versus the number of particles $`n`$ for $`r_1=0,0.1,0.3,0.5,0.7,`$ $`0.9,0.95,0.98,1`$, repectively. It can be seen that the entropy index $`\mu `$ firstly increases with the increasing of $`r_1`$ starting from $`r_1=0`$, but when $`r_1`$ is big (near to unity) it turns to

decreasing. The reason is that $`r_1=0`$ and $`r_1=1`$ correspond to single-$`\alpha `$ model with $`\alpha =0.2`$ and 0.8 respectively. The entropy indices for these two cases are in general smaller than that of the cases with the two $`\alpha `$ values mixed. The only exception is when only a very little fraction ( $`<0.1`$) of big $`\alpha `$ (0.8) is mixed in. In this case, the double-$`\alpha `$-mixing model will have a smaller entropy index than the single-$`\alpha `$ model with $`\alpha =0.8`$.

Fig.4 Entropy index versus

multiplicity for different

double-$`\alpha `$ model

The variation of entropy index in double-$`\alpha `$-mixing model with the mixing ratio $`r_1`$ is plotted in Fig.5 for $`\alpha _1,\alpha _2=0.1,0.9;0.2,0.8;0.3,0.7;0.4,0.6`$, respectively. The increasing of $`\mu `$ with $`r_1`$ when the latter is not very big and the decreasing of $`\mu `$ when $`r_1`$ near to unity can be seen in all the cases.

In order to show the fine structure near $`r_1=1`$ the curves for $`r_1,r_2=`$ 0.1,0.9 (case-A) and 0.4,0.6 (case-B) are magnified and shown in the small figure inside Fig.5. The crossing of these two curves is understandable. This is because $`\alpha _1^{(\mathrm{A})}<\alpha _1^{(\mathrm{B})}`$ while $`\alpha _2^{(\mathrm{A})}>\alpha _2^{(\mathrm{B})}`$.

Fig.5 Entropy index in double-$`\alpha `$

model as function of mixing ratio

Therefore, when $`r_1`$ is small (1-$`r_1`$ is big) $`\mu ^{(\mathrm{A})}>\mu ^{(\mathrm{B})}`$, but when $`r_1`$ is big (1-$`r_1`$ is small) $`\mu ^{(\mathrm{B})}>\mu ^{(\mathrm{A})}`$.

In summary, it is found in this letter that even in the pure single-$`\alpha `$ random-cascading model without putting in any particle there exists erraticity behavior and the corresponding entropy indices do not vanish. This observation is in contrary to the previous expectation. It means that the dynamical fluctuations in a pure single-$`\alpha `$ model already fluctuate event-by-event.

This fluctuation is strengthened when different values of $`\alpha `$ are mixed in the model. Taking double-$`\alpha `$-mixing model as example, the variation of the entropy index with the mixing ratio of the two $`\alpha `$’s is studied. It is found that the entropy index $`\mu `$ firstly increases with the increasing of $`r_1`$ and then decreases. The maximum is located near $`r_1=1`$ when $`\alpha _1<\alpha _2`$.

The influence of particle number on the results when particles are putted into the system is also investigated. In all cases the entropy index calculated by particle number is bigger than that in pure $`\alpha `$ model and tends to the latter only when the particle number tends to infinity.

Acknowledgement

The authors show sincere thanks to the patient direction given by Liu Lianshou and the warm help given by Liu Fuming.


text:
# Optimal detection of sources on a homogeneous & isotropic background

## 1 Introduction

One of the challenges in the analysis of 1D spectra, 2D images or 3D volumes in Astrophysics and Cosmology is to overcome the problem of separating a localized signal from a background. In particular, we are interested in localized sources with central symmetry and a background that we shall assume with the properties of homogeneity and isotropy and it will be characterized by a power spectrum. Typical cases in the 1D case include: a) the spectra of QSOs, how to detect/extract absorption lines from a noisy spectrum once a profile is assumed, b) time series analysis where a localized emission is superposed on a background, c) point and extended sources to be detected in time-ordered data where the scanning strategy (for satellites like Planck) is affected by the well-known $`1/f`$ noise. In the 2D case, we mention as typical cases: a) cleaning of images to detect typical astrophysical sources on a white noise background, b) the detection/extraction of point sources and extended sources (clusters of galaxies) on microwave/IR maps where the background is dominated by white noise or intrinsic CMB signal or galactic emission. In the 3D case, we remark as an example: a) the detection of clusters and large-scale structure in 3D catalogs.

The classical treatment to remove the background has been filtering. Low and high-pass filters reduced the contribution from the high and low frequencies present in the spectrum or image. In general, this process is extremely inefficient in dealing with very localized sources. The reason for that is that a very localized source can be decomposed in Fourier space but many waves are needed (infinite for a delta distribution!), so if low/high-pass filters are applied then at the end many artefacts (rings in 2D images) usually appear surrounding the sources.

A very important application of these principles is the detection of sources in two-dimensional astronomical images. Several packages are commonly used for this task, such as DAOfind (Stetson 1992), used to find stellar objects in astronomical images, and SExtractor (Bertin & Arnouts 1996). When trying to detect sources, two different problems have to be solved: first, it is necessary to take account for the background variations across the image. In addition, if instrumental noise (i. e. white noise) appears, it should be removed as far as possible in order to increase the SNR ratio. SExtractor estimates the local background on a grid of size set by the user and then approximates it by a low-order polynomial. Once such a background is defined the detection of sources at a certain level can be established through sets of connected pixels above certain threshold. DAOfind implicitly assumes that the background is smooth and its characteristic scale is very much larger than the scale of the stars. But in the case where the characteristic scale of variation of the background is approximately the scale of the structures the previous schemes obviously fail. An incorrect estimation of the background leads to a biased estimation of the amplitude of the sources. An example is a typical image of the cosmic microwave background radiation (CMB) at a resolution of several arcmin. If the intrinsic signal is due to a cold dark matter model (i. e. the characteristic scale of the background is of the order of $`10`$ arcmin), then the separation of point sources with the same characteristic scale becomes a very difficult task.

To deal with the instrumental noise the traditional procedure is filtering (e. g. Gaussian window). DAOfind filters with an elliptical Gaussian that mimics the stellar psf and then it performs the detection looking for peaks above certain threshold. SExtractor includes the possibility of filtering with several kind of filters (Top Hat, Gaussian, Mexican Hat and even self-made filters for every particular situation). An obvious advantage of this procedures (background estimation plus filtering) is that no a priori information on the structures is needed. A serious drawback is that the choice of the filter will have a great influence on the final result. The choice of filter depends on many factors, including in most cases personal preferences.

In this context, it is necessary to find a systematic way to determine the optimal filter for every case.

Other methods have been used to separate different components given several microwave maps: Wiener filtering (WF) and maximum entropy methods (MEM). Regarding point sources, WF assumes a mean spectral flux dependence together with other assumptions for the other components (Tegmark & Efstathiou, 1996; Bouchet et al. 1999) whereas MEM assumes that the point sources are distributed like noise (Hobson et al. 1999). These methods are powerful regarding extended sources like clusters of galaxies because they use the concrete spectral dependence for the Sunyaev-Zeldovich effect. It is clear that the unknown spectral dependence for point sources remark the inefficiency of the previous methods.

A possible solution to overcome this problem came with the usage of wavelets. These are localized bases that allow a representation of a local object due, in general, to their basic properties: spatial and frequency localization (as opposed to the trigonometric functions appearing in Fourier decomposition). We remark at this point the success of such a technique dealing with simulated microwave maps: the ”Mexican Hat” wavelet can be used in a nice way (no extra assumptions about the background) to detect/extract point sources and create a simulated catalog (Cayón et al. 2000, Vielva et al. 2000). Two advantages emerge: on the one hand, one localizes the structures associated to the maxima of the wavelet coefficients and, what is more remarkable, we gain in the detection (as compared to real space) due to the amplification effect because at the scale of the source the background is not contributing to the dispersion. One relevant question concerns the possibility to find optimal filters. Tegmark & Oliveira-Costa (1998) introduced a filter that minimizes the variance in the map. With this method one can identify a big number of point sources in CMB maps. However, they failed to introduce the appropriate constraints in the minimization problem, i. e. the fact that we have a maximum at the source position at the scale defined by the source in order not to have spurious identifications.

Thus, this type of analysis lead us to the following questions: is there an optimal filter (or better pseudo-filter) given the source profile and the power spectrum of the background?, Is the ”Mexican Hat” wavelet the optimal pseudo-filter dealing with point sources? In order to answer these questions, we will assume that the sources can be approximated by structures with central symmetry given by a profile $`\tau (x)`$, $`x|\stackrel{}{x}|`$ with a characteristic scale (e.g. a single maximum at its center and rapid decay at large distances). If the nD image contains different types of sources, a number of pseudo-filters adapted to each profile must be used to detect them. A possible generalization to include more general profiles is under study. The background will be modelled by a homogeneous and isotropic random field given in terms of the power spectrum $`P(q)`$, $`q|\stackrel{}{q}|`$. In particular, we shall explore a scale-free spectrum $`P(q)q^\gamma `$ that includes the cases of white noise ($`\gamma =0`$), $`1/f`$ noise ($`\gamma =1`$), etc. Moreover, any spectrum of physical interest often can be locally approximated by a power-law. If the characteristics of the noise are not known a priori it can be always estimated directly from the nD image. We consider the n-dimensional case and make special emphasis on the analysis of spectra ($`n=1`$), 2D images ($`n=2`$) and 3D volumes ($`n=3`$). In all the calculations we assume that the overlapping by nearby sources is negligible and also that their contribution to the total power spectrum is also negligible. All of this is a very good approximation at least above a certain flux level.

An analytical approach to get the optimal pseudo-filter is presented in section 2. Section 3 deals with properties of the optimal pseudo-filters on real space. Section 4 introduces the concepts of detection level and gain. Sections 5 and 6 are dedicated to the important cases of sources with profiles described by a Gaussian and an exponential, respectively. An example of the performance of optimal pseudo-filters applied to simulated one-dimensional data is presented in section 7. Section 8 deals with the extraction of the sources. Conclusions are summarized in section 9.

## 2 Optimal pseudo-filter

Let us consider an n-dimensional ($`n=1,2,3`$) image with data values defined by

$$y(\stackrel{}{x})=s(x)+n(\stackrel{}{x}),$$
(1)
where $`x`$ is the spatial coordinate (in the 1D case $`x`$ can be also time, when we are dealing with time-ordered data sets) and $`s(x),x|\stackrel{}{x}|`$, represents a source with central symmetry placed at the origin with a characteristic scale (e. g. a single maximum at its center and rapid decay at large distances) and $`n(\stackrel{}{x})`$ a homogeneous & isotropic background (random field) with mean value $`<n(\stackrel{}{x})>=0`$ and characterized by the power spectrum $`P(q),q|\stackrel{}{q}|`$ (this can represent instrumental noise and/or a real background),

$$<n(\stackrel{}{q})n^{}(\stackrel{}{q^{}})>=P(q)\delta ^n(\stackrel{}{q}\stackrel{}{q^{}}),$$
(2)
where $`n(\stackrel{}{q})`$ is the nD Fourier transform ($`n(\stackrel{}{q})=(2\pi )^{n/2}𝑑\stackrel{}{x}e^{i\stackrel{}{q}\stackrel{}{x}}n(\stackrel{}{x})`$), symbol $`n^{}`$ represents the complex conjugate of $`n`$, $`\stackrel{}{q}`$ is the wave vector and $`\delta ^n`$ is the nD Dirac distribution.

Let us introduce a spherical (centrally symmetric) pseudo-filter, $`\mathrm{\Psi }(\stackrel{}{x};\stackrel{}{b},R)`$, dependent on $`n+1`$ parameters .

$$\mathrm{\Psi }(\stackrel{}{x};R,\stackrel{}{b})=\frac{1}{R^n}\psi (\frac{|\stackrel{}{x}\stackrel{}{b}|}{R}),$$
(3)
where $`\stackrel{}{b}`$ defines a translation whereas $`R`$ defines a scaling. Then, we define the pseudo-filtered field $`w(R,\stackrel{}{b})`$

$$w(R,\stackrel{}{b})=𝑑\stackrel{}{x}y(\stackrel{}{x})\mathrm{\Psi }(\stackrel{}{x};\stackrel{}{b},R).$$
(4)
We do not assume a priori the positiveness of $`\mathrm{\Psi }`$. The previous convolution can be written as a product in Fourier space, in the form

$$w(R,\stackrel{}{b})=𝑑\stackrel{}{q}e^{i\stackrel{}{q}\stackrel{}{b}}y(\stackrel{}{q})\psi (Rq).$$
(5)
where $`y(\stackrel{}{q})`$ and $`\psi (q)`$ are the Fourier transforms of $`y(\stackrel{}{x})`$ and $`\psi (\stackrel{}{x})`$, respectively. Because of the central symmetry assumed for the pseudo-filter, $`\psi (q)`$ depends only on the modulus of $`\stackrel{}{q}`$. A simple calculation -taking into account eqs. (1) and (2)- gives the average at the origin $`\stackrel{}{b}=0`$, $`<w(R,\stackrel{}{0})>`$, and the variance, $`\sigma _w^2(R)=<w^2(R,\stackrel{}{b})><w(R,\stackrel{}{b})>^2`$, of the pseudo-filtered field

$$<w(R,\stackrel{}{0})>=\alpha 𝑑qq^{n1}s(q)\psi (Rq),\sigma _w^2(R)=\alpha 𝑑qq^{n1}P(q)\psi ^2(Rq).$$
(6)
where $`q=|\stackrel{}{q}|`$, $`\alpha =2,2\pi ,4\pi `$ for $`n=1,2,3`$, respectively, (for n-dimensions $`\alpha =2\pi ^{n/2}\mathrm{\Gamma }^1(n/2)`$) and the limits in the integrals go from $`0`$ to $`\mathrm{}`$.

Now, we are going to express the conditions in order to obtain an optimal pseudo-filter for the detection of the source $`s(x)`$ at the origin. One basic idea is to find a pseudo-filter such that when the original image is filtered with a scale $`R_oR_s`$ -being $`R_s`$ the characteristic scale of the source- one obtains the maximum detection level $`𝒟_w`$

$$𝒟_w\frac{<w(R,\stackrel{}{0})>}{\sigma _w(R)}.$$
(7)
Taking into account the fact that the source is characterized by a single scale $`R_s`$, other basic idea is to generate a filter giving the maximum contribution at the center of the source at a filtering scale $`R_oR_s`$. Finally, we would like to estimate directly the amplitude of the source by the previous number. Therefore, taking into account these basic ideas we will introduce from the mathematical point of view the optimal pseudo-filters.

By definition a pseudo-filter will be called optimal if the following conditions are satisfied:

i) there exists a scale $`R_o`$ such that $`<w(R,\stackrel{}{0})>`$ has a maximum at that scale, ii) $`<w(R_o,\stackrel{}{0})>=s(0)A`$, i. e. $`w(R,\stackrel{}{0})`$ is an unbiased estimator of the amplitude of the source and iii) the variance of $`w(R,\stackrel{}{b})`$ has a minimum at the scale $`R_o`$, i. e. we have an efficient estimator. As a by-product, the ratio given by eq. (7) will be maximum. We remark that no other information about the source profile is assumed, so “optimal” must be understood in the previous sense.

By introducing the profile $`\tau (x)`$ of the source, $`s(x)=A\tau (x)`$, the condition ii) and the equation (6) give the constraint

$$𝑑qq^{n1}\tau (q)\psi (R_oq)=\frac{1}{\alpha },$$
(8)
whereas the condition i) gives the constraint

$$𝑑qq^{n1}\tau (q)\psi (R_oq)[n+\frac{dln\tau }{dlnq}]=0.$$
(9)
So, the problem is reduced to the functional minimization (with respect to $`\psi `$) of $`\sigma _w^2(R)`$ given by equation (6) with the constrains given by equations (8) and (9). This minimization incorporate these constraints through a couple of Lagrangian multipliers. The solution (optimal pseudo-filter) is found to be

$$\stackrel{~}{\psi }(q)\psi (R_oq)=\frac{1}{\alpha }\frac{\tau (q)}{P(q)}\frac{1}{\mathrm{\Delta }}[nb+c(na+b)\frac{dln\tau }{dlnq}],\mathrm{\Delta }=acb^2,$$
(10)
$$a𝑑qq^{n1}\frac{\tau ^2}{P},b𝑑qq^{n1}\frac{\tau }{P}\frac{d\tau }{dlnq},c𝑑qq^{n1}\frac{1}{P}[\frac{d\tau }{dlnq}]^2.$$
(11)
Therefore, we have obtained analytically the functional form of the pseudo-filter (its shape and characteristic scale are associated to the source profile and power spectrum). It is clear that assuming an adimensional dependence $`\tau (x/R_s)`$, where $`R_s`$ is the characteristic scale of the source, then such scale will appear explicitly in the form $`\stackrel{~}{\psi }(qR_s)`$. Obviously, we assume all the differentiable and regularity conditions at $`q0,\mathrm{}`$ for $`\tau `$ and $`P(q)`$ in order to have finite expressions for $`a,b,c`$. Generically, $`\psi `$ is a pseudo-filter, i. e. it is not positive (filter). Let us remark that if we assume the behavior $`\frac{\tau (q)}{P(q)}0`$, $`\frac{1}{P}\frac{d\tau }{dlnq}0`$ for $`q0`$ then $`\psi (q)0`$ and $`\mathrm{\Psi }`$ is a ”compensated” filter, i. e. $`𝑑\stackrel{}{x}\mathrm{\Psi }=0`$. Strictly speaking there is another condition to be satisfied to get the reconstruction of the image and thus to have a wavelet: $`𝑑qq^1\psi ^2(q)<\mathrm{}`$ (the admissibility condition).

Taking into account eq. (5) the amplitude will be estimated as

$$A=w(R_o,\stackrel{}{0})=𝑑\stackrel{}{q}y(\stackrel{}{q})\stackrel{~}{\psi }(q),$$
(12)
where $`\stackrel{~}{\psi }`$ is given by eq.(10).

## 3 Optimal pseudo-filter on real space

The equation (10) can be written on real space as follows

$$\mathrm{\Psi }(\stackrel{}{x};R_o,\stackrel{}{0})=\frac{1}{R_o^n}\psi (\frac{|\stackrel{}{x}|}{R_o})=\frac{1}{\alpha \mathrm{\Delta }}[(nb+c)F(x)(na+b)G(x)],$$
(13)
where $`F`$ and $`G`$ are the inverse Fourier transform of $`\frac{\tau }{P}`$ and $`\frac{1}{P}\frac{d\tau }{dlnq}`$, respectively.

For a flat background, i. e. $`P=constant`$, and assuming the behaviour $`q^{\frac{1}{2}(n+1)}\tau 0`$ when $`q\mathrm{}`$, one obtains

$$b=\frac{n}{2}a,F=\frac{\tau (x)}{P},G(x)=\frac{1}{P}[n\tau +\frac{d\tau }{dlnx}].$$
(14)
If we also assume a Gaussian profile, i. e. $`\tau =e^{\frac{x^2}{2\theta ^2}}`$, one finds the pseudo-filter

$$\mathrm{\Psi }(\stackrel{}{x};\theta ,\stackrel{}{0})=\frac{1}{\pi ^{n/2}\theta ^n}e^{\frac{x^2}{2\theta ^2}}[1+\frac{n}{2}(\frac{x}{\theta })^2],$$
(15)
that is a useful formula to be used for the detection of nD-Gaussian structures on nD-images (e.g. spectra, 2D images or 3D volumes) on a flat background.

On the other hand, if one assumes a Gaussian profile but a non-flat spectrum one can easily find

$$\mathrm{\Psi }(\stackrel{}{x};\theta ,\stackrel{}{0})=\frac{1}{\alpha \mathrm{\Delta }}[(nb+c)F(x)(na+b)\theta ^2^2F(x)],$$
(16)
being $`F(x)`$ the Fourier transform of $`\frac{\tau }{P}`$.

## 4 Detection level, gain and reliability

Taking into account the previous expression (10), one can calculate the detection level (see equation(7))

$$𝒟_w=A[\frac{\alpha \mathrm{\Delta }}{n^2a+2nb+c}]^{1/2}.$$
(17)
On the other hand, we can calculate the dispersion of the field $`n(\stackrel{}{x})`$

$$\sigma _b^2=\alpha 𝑑qq^{n1}P(q),$$
(18)
that allows to define a detection level on real space as

$$𝒟=\frac{A}{\sigma _b}.$$
(19)
The gain going from real space to pseudo-filter space is defined by

$$g\frac{𝒟_w}{𝒟}=\frac{\sigma _b}{\sigma _w(R_o)}.$$
(20)
If the background has a characteristic scale different from the scale of the structures (sources) to be detected, it is obvious that $`\sigma _w<\sigma _b`$, so that we have a real gain going from real space to pseudo-filter space.

The identification of sources as peaks above a high threshold (e. g. $`3\sigma _w`$) in pseudo-filter space gives a low probability of false detections (reliability) because if the background has a characteristic scale different from the sources then everything detected with our method is real, but if both scales are comparable one can give an estimate based on the fluctuations of the background. For instance, in the case of a Gaussian background, false detections above $`3\sigma _w`$, due to the Gaussian background, appear with a probabilty $`1.5\times 10^3`$ based on the formula

$$Pr(w>A)=\frac{1}{2}erfc(\frac{A}{2^{1/2}\sigma _w}).$$
(21)
To select the false detections from the real ones one can study the pseudo-filter profile nearby any real source (see the last paragraph of section 8).

Regarding the completeness (i. e. how many real sources we miss with our method), this is a complicated topic because the background can slightly modify the location of the peaks and their amplitude. We will address this problem via numerical simulations (see section 7).

## 5 Gaussian source on a background

In many physical applications the standard response of the instruments can be approximated by a Gaussian function. In particular, the Point Spread Function (PSF) for many telescopes is of Gaussian type. Other more specific applications are related to the Cosmic Microwave Background (CMB), where the antennas used are well approximated by a Gaussian. Dealing with absorption systems associated to QSOs, if the absorption line is not saturated and is dominated by thermal motions, then the line is usually approximated by an inverted Gaussian inserted in a continuum plus noise.

Let us assume that the source and the background can be represented by

$$s(x)=Ae^{\frac{x^2}{2\theta ^2}},P=Dq^\gamma .$$
(22)
The structure to be detected could have an intrinsic Gaussian profile or it could be a point source in n-dimensions observed with an instrument that can be modelled through a Gaussian pattern with a beam size $`\theta `$. The background can be described by a scale-free spectrum. In this case: $`\tau (q)=\theta ^ne^{\frac{1}{2}(q\theta )^2}`$, and equations (11) give

$$a=\frac{\theta ^{n\gamma }}{2D}\mathrm{\Gamma }(m),b=ma,c=m(1+m)a,m\frac{n+\gamma }{2}$$
(23)
and the pseudo-filter is

$$\stackrel{~}{\psi }(q)=\frac{1}{\alpha \mathrm{\Gamma }(m)}(q\theta )^\gamma e^{\frac{1}{2}(q\theta )^2}[2+\gamma n+\frac{n\gamma }{m}(q\theta )^2].$$
(24)
Taking into account the $`q`$ behaviour in this formula, one obtains a compensated filter (i. e. $`\psi (q=0)=0`$) if $`\gamma >0`$ or $`\gamma =n2`$. In figure Optimal detection of sources  on a homogeneous & isotropic background appear the optimal pseudo-filters for the 1D, 2D and 3D cases, respectively, and scale-free power spectrum with indexes $`\gamma =0,1,2,3`$. There is a degeneration in the case $`n=3`$, where the $`\gamma =1`$ line overlaps with the $`\gamma =3`$ line, and in the case $`n=2`$, where the $`\gamma =0`$ line overlaps with the $`\gamma =2`$ one. This degeneration can be deduced directly from equation (23).

On the other hand, the detection level in pseudo-filter space is given by equation (18)

$$𝒟_w=A[\frac{2\alpha }{D}\frac{\mathrm{\Gamma }(1+m)}{4m+(n\gamma )^2}]^{1/2}\theta ^{\frac{n\gamma }{2}}.$$
(25)
Finally, it is interesting to remark that the cases $`\gamma =n`$ and $`\gamma =n2`$ give the same pseudo-filter

$$\stackrel{~}{\psi }(q)=\frac{2}{\alpha \mathrm{\Gamma }(n)}(q\theta )^ne^{\frac{x^2}{2\theta ^2}}.$$
(26)
Therefore, in the cases $`n=2,\gamma =0,2`$ the Mexican Hat is found to be the optimal pseudo-filter. This justify the use of this wavelet to detect point sources in CMB maps (Cayón 2000, Vielva 2000).

a) Gaussian source and white noise:

This subcase corresponds to $`P=constant`$ or $`\gamma =0`$, and the pseudo-filter is

$$\stackrel{~}{\psi }(q)=\frac{2}{\alpha \mathrm{\Gamma }(\frac{n}{2})}e^{\frac{1}{2}(q\theta )^2}[1\frac{n}{2}+(q\theta )^2],$$
(27)
that gives a pseudo-filter for the analysis in the different dimensions except for $`n=2`$ that gives the Mexican Hat wavelet ($`\psi =\frac{1}{\pi }e^{\frac{1}{2}(q\theta )^2}(q\theta )^2`$).

The detection level is given by equation (24):

$$𝒟_w=A[\frac{\alpha \mathrm{\Gamma }(\frac{n}{2})}{2D(1+\frac{n}{2})}]^{1/2}\theta ^{n/2}.$$
(28)
We have calculated the contribution of sources to the power spectrum in order to estimate their influence in calculating the pseudo-filter. We arrive to the conclusion that if the signal/noise ratio (i. e. dispersion associated to the sources over dispersion associated to the background) is $`\frac{\sigma _s}{\sigma _b}<0.6(\frac{l_p}{\theta })`$, being $`l_p`$ the pixel scale and $`\theta `$ the width of the source, then the extra contribution to the pseudo-filter coefficients $`a,b,c`$ is less than a $`10\%`$.

b) Gaussian source and $`1/f`$ noise:

Let us assume a source with a Gaussian profile and a background represented by $`1/f`$ noise, i.e. $`P=Dq^1`$ or $`\gamma =1`$. In this case: $`\tau (q)=\theta ^ne^{\frac{1}{2}(q\theta )^2}`$, and equation (22) gives the pseudo-filter

$$\stackrel{~}{\psi }(q)=\frac{1}{\alpha \mathrm{\Gamma }(\frac{n+1}{2})}e^{\frac{1}{2}(q\theta )^2}(q\theta )[3n+2\frac{n1}{n+1}(q\theta )^2],$$
(29)
For instance, in the case $`n=1`$ one has the wavelet $`\psi =(q\theta )e^{\frac{1}{2}(q\theta )^2}`$, that is the optimal pseudo-filter to be used to detect a Gaussian signal on 1D spectra. In this case the detection level is (see equation(24)) $`𝒟_w=AD^{1/2}`$.

## 6 Exponential source on a background

Typical example in astrophysics is the exponential disk associated to spiral galaxies. Another interesting application could be in some areas of physics where the profile expected for the signal associated to the detection of some particles could be of exponential type.

Let us assume that the source and background can be represented by

$$s(x)=e^{\frac{x}{\lambda }},P=Dq^\gamma .$$
(30)
In this case:

$$\tau (q)=\beta \lambda ^n[1+(q\lambda )^2]^{\frac{n+1}{2}},$$
(31)
where $`\beta =(\frac{2}{\pi })^{1/2},1,2(\frac{2}{\pi })^{1/2}`$ for $`n=1,2,3`$, respectively, and equations (11) give

$$a=\frac{\beta ^2\lambda ^{n\gamma }}{2D}\frac{\mathrm{\Gamma }(m)\mathrm{\Gamma }(1+\frac{n\gamma }{2})}{\mathrm{\Gamma }(n+1)},$$
(32)
$$b=ma,c=\frac{n+1}{n+2}m(1+m)a,m\frac{n+\gamma }{2}$$
(33)
and the pseudo-filter is

$$\stackrel{~}{\psi }(q)=\frac{2L}{\alpha \beta }(q\lambda )^\gamma [1+(q\lambda )^2]^{\frac{n+1}{2}}[1+\frac{\gamma n}{2}(n+1)+M\frac{(q\lambda )^2}{1+(q\lambda )^2}],$$
(34)
$$L\frac{\mathrm{\Gamma }(n+1)}{\mathrm{\Gamma }(\frac{n+\gamma }{2})\mathrm{\Gamma }(2+\frac{n\gamma }{2})},M\frac{n\gamma }{n+\gamma }(n+1)(n+2).$$
(35)
In figure Optimal detection of sources  on a homogeneous & isotropic background appear the optimal pseudo-filters for the 1D, 2D and 3D cases, respectively, and power spectrum with indexes $`\gamma =0,1,2,3`$. The filter profiles are more extended than in the Gaussian source, as one can expect from the more gentle fall of the exponential source.

An interesting case is $`\gamma =n`$, then the pseudo-filter is

$$\stackrel{~}{\psi }(q)=\frac{2n}{\alpha \beta }(q\lambda )^n[1+(q\lambda )^2]^{\frac{n+1}{2}}.$$
(36)
a) Exponential source and white noise:

For this subcase $`\gamma =0`$, then equation (36) leads to the pseudo-filter

$$\stackrel{~}{\psi }(q)=\frac{2}{\alpha \beta }\frac{\mathrm{\Gamma }(n+1)}{\mathrm{\Gamma }(\frac{n}{2})\mathrm{\Gamma }(2+\frac{n}{2})}[1+(q\lambda )^2]^{\frac{n+1}{2}}[1+\frac{n}{2}(n+1)+\frac{(n+1)(n+2)(q\lambda )^2}{1+(q\lambda )^2}].$$
(37)
For instance, for an exponential structure to be optimally detected in a 1D spectrum, we must use

$$\stackrel{~}{\psi }(q)=\frac{8}{(2\pi )^{1/2}}\frac{(q\lambda )^2}{[1+(q\lambda )^2]^2}.$$
(38)
b) Exponential source and $`1/f`$ noise:

An interesting case is $`n=\gamma =1`$, equation (36) gives

$$\stackrel{~}{\psi }(q)=(\frac{\pi }{2})^{1/2}\frac{(q\lambda )}{1+(q\lambda )^2}.$$
(39)
## 7 Simulations of one-dimensional Gaussian sources and $`1/f`$ noise

In order to test some of the ideas proposed in previous sections, we simulated the case of one-dimensional Gaussian sources on a background. The kind of background simulated is the well-known $`1/f`$ noise. This kind of noise appears very often in many devices in experimental physics. Further simulations with 2-dimensional data and realistic realizations of noise will be carried on in future work.

For the sake of simplicity, all the simulated sources have the same amplitude and that is set to be 1 (in arbitrary units). 100 of these sources were deployed over a 32768 pixel ’field’. The number of sources and the size of the field were selected in order to have enough sources for statistical studies, to avoid (as far as possible) the overlapping of the sources and to minimize the contribution of the sources to the total dispersion of the simulations. The width of the Gaussian profiles was chosen to be $`\theta 3\theta _p`$: this is the case for a pixel of $`1.5^{}`$ and a Gaussian with a FWHM of $`5^{}`$. Noise was added so that the signal-to-noise ratio of the sources, defined as the ratio $`A/\sigma _b`$ (where $`A`$ is the amplitude of the source and $`\sigma _b`$ is the standard deviation of the noise), assumes values of $`2,3,4`$ and $`5`$. Finally, the optimal filter, given by eq. (29) with $`n=1`$, was applied to the image.

To compare with a more traditional filtering scheme, we filtered the images also with a Gaussian of width equal to $`\theta `$ and a Mexican Hat wavelet of width equal to $`\theta `$. This is a rather naive usage of the Mexican Hat wavelet and the Gaussian source because the optimal width for these filters in the general case is not the source scale (Cayón 2000, Vielva 2000), but it serves us well because what we intend is to compare how do filters work when we have no further information about the data (i.e., the optimal scale, which is different for each background). The result of these simulations is shown in tables 1 and 2. Table 1 refers to the original simulations. It shows the original signal-to-noise ratio of each simulation as well as statistical quantities of interest such as the dispersion of the map and the mean amplitude of the sources in it. Finally, it shows the number of sources directly detected from the simulations above $`3\sigma `$ and $`5\sigma `$ thresholds and the number of spurious detections above these tresholds. As expected, only a few sources are detected, except for the most favorable cases (high original signal-to-noise rato and low detection threshold). The small bias in the mean measured amplitude is due to pixelization effects. Table 2 refers to the simulations in table 1 after filtering with a Gaussian of width $`\theta `$, a Mexican Hat wavelet of width $`\theta `$ and the optimal pseudo-filter. Each row in table 2 corresponds to the same row in table 1.

A Gaussian filter smoothes the image, removing small-scale noise. It also smoothes the source peaks, thus lowering the amplitude of detected sources. For the case of $`1/f`$ noise the dominant fluctuations appear at large scales and are not affected by the Gaussian filter. The large-scale features may contribute to contamination in two different ways: they can conceal sources in large ’valleys’ and can produce spurious peaks. None of these effects can be avoided with a Gaussian filter. On the other hand, the smoothing effect of the Gaussian filter takes place normally and lowers the amplitude of the sources. Therefore, the number of true detections is smaller than in the non-filtered image, and the spurious detections are not removed even in the highest $`s/n`$ case. The gains, indicated in column 6, clearly reflects this situation ($`g<1`$).

The Mexican Hat wavelet has a better performance under $`1/f`$ conditions. The Mexican Hat removes large-scale fluctuations, allowing the ’hidden’ sources to arise above the detection threshold. For example, in the case of original $`s/n=2.95`$ there were 47 sources above the $`3\sigma `$ level and only 1 above the $`5\sigma `$ level. After filtering with the Mexican Hat, there are 93 detections above $`3\sigma `$ level and 64 above the $`5\sigma `$ level, a significant improvement. The number of spurious sources remains almost untouched.

The optimal pseudo-filter also deals with the large-scale structure. It is constructed to enhance all fluctuations in the source scale, while removing fluctuations that arise in other scales. In addition, it is required to be unbiased with respect to the amplitude. In practice, the amplitude is slightly underestimated due to the propagation of pixelization effects. This small bias is lower than a $`10\%`$ and can be calibrated in any case. In the $`1/f`$ case the number of true detections is higher than in the Mexican Hat case and the number of spurious sources is comparable or slightly reduced. Only in the case of low initial signal-to-noise ratio the number of spurious detections is greater. This is due to the fact that this pseudo-filter enhances all fluctuations in the source scale. Future work will take care of this weakness of the method including more information about the shape of the sources. For the case of initial signal-to-noise ratio of 2.95 we find 94 sources (of 100) and 8 spurious detections (a reliability close to $`10\%`$) above the $`3\sigma `$ level, a result very similar to the obtained with the Mexican Hat. Above the $`5\sigma `$ detection level the optimal pseudo-filter finds 79 sources where the Mexican Hat found only 64. The number of spurious sources have not increased significantly (from 4 to 5). For higher initial signal-to-noise ratios the completeness and reliability quickly improve.

The gain obtained with the optimal pseudo-filter is greater than the one obtained with the Mexican Hat. It can be analytically calculated, using eqs. (6,17) for the pseudo-filter and its equivalents for the Mexican Hat, leading to:

$$g_{op}/g_{mh}=\left[\frac{4}{\pi }\frac{\mathrm{\Gamma }(\frac{1+\gamma }{2})\mathrm{\Gamma }(\frac{5\gamma }{2})}{1+\frac{(1\gamma )^2}{2(1+\gamma )}}\right]^{1/2}$$
(40)
for the one-dimensional case. This formula holds while $`\gamma 1`$. According to equation (40), the ratio $`g_{op}/g_{mh}`$ is 1.41 in the case $`\gamma =0`$ and 1.13 in the case $`\gamma =1`$. The mean observed ratio in the simulations is 1.31 and 1.08 respectively, and fits well with our expectations. As a conclusion we have that optimal filter gives higher gains than the classical Mexican Hat filter.

In figure Optimal detection of sources  on a homogeneous & isotropic background an example of the simulations is shown. In the top panel there is a 500 pixel wide subsection of the $`s/n=3`$, $`\gamma =1`$ simulation. This subsection corresponds to a region in which the large-scale noise has a positive value. Four sources are present in this area, all of them arising above the $`3\sigma `$ level (indicated with a dotted line). The position of the sources are marked with an asterisk in the lower panel. Additionally, there are three peaks, corresponding to background fluctuations, that arise above the $`3\sigma `$ level. The second panel from the top shows the image after filtering with the optimal filter. The large-scale features have been removed and also the small-scale noise is reduced. The sources have been amplified with respect to the original map and now all of them reach the $`3\sigma `$ level but the spurious peaks have been removed. The amplitudes of the sources remain unbiased and close to the true value of 1. In the third panel from the top there is the image after filtering with a Gaussian. The whole image has been smoothed and now one of the sources barely reaches the $`3\sigma `$ level. The large-scale fluctuations remain untouched and all the spurious peaks remain in the filtered image. In the bottom panel we see the image after filtering with the Mexican Hat. The large-scale fluctuations are also removed as well as the small-scale noise, as in the case of the optimal pseudo-filter. Nonetheless, the gain is lower and only three of the sources reach the $`3\sigma `$ level. Additionally, it is found that the small-scale noise removal is less efficient in the case of the Mexican Hat.

## 8 Extraction of sources

The optimal pseudo-filter gives the position and an unbiased estimator of the amplitude of the source. We propose to make the extraction of the source on real space, i. e. one subtracts the function $`A\tau (|\stackrel{}{x}\stackrel{}{x_o}|)`$, being $`\tau `$ the given profile and $`A`$ the estimated amplitude, at the position of the source $`\stackrel{}{x_o}`$.

From the practical point of view, in order to select the appropriate sources (with a given scale and avoiding to select spurious detections if the background and/or noise are manifest at scales comparable to the sources) we can operate with the optimal pseudo-filter at other different scales $`R`$ as given by equation (10) but with $`\stackrel{~}{\psi }(qx),xR/R_o`$. If the scale that gives the maximum do not correspond to the scale we are looking for then this is a spurious source (or another type of source with a different scale). As an additional check, we can calculate the source profile in the pseudo-filter space nearby any real source, e. g. for a Gaussian profile the behaviour around the maximum $`R_oR_s`$ must be

$$<w(R)>=Ax^\gamma (\frac{2}{1+x^2})^m[1+\frac{n\gamma }{2}\frac{1x^2}{1+x^2}],x\frac{R}{R_o},m\frac{n+\gamma }{2},$$
(41)
and an analogous behaviour can be found for the exponential profile. If a detected source do not follows such a behaviour then it would be consider as a false detection and must be deleted from the initial catalog.

## 9 Conclusions

We have introduced for the first time (as far as the authors know) the concept of optimal pseudo-filter to detect/extract spherical sources on a background modelled by a (homogeneous & isotropic) random field characterized by its power spectrum. We have obtained a generic analytical formula that allows to calculate such a pseudo-filter either in Fourier or real space as a function of the source profile and power spectrum of the background. The psesudofilter is an unbiased an efficient estimator of the amplitude of the sources.

We have applied the previous formula to the cases of a Gaussian and an exponential profile and studied scale-free spectra. In particular, we have remarked the interesting cases of white noise and $`1/f`$ noise. We have calculated the detection level for the physically interesting cases of spectra, images and volumes. For some particular cases, the optimal pseudo-filters are wavelets (e. g. a Gaussian source embedded in white noise in the 2D case).

We have simulated Gaussian sources embedded in a $`1/f`$ noise in order to see the performance of the optimal filter against the Mexican Hat wavelet. In the last case the gain is lower, the noise removal is less efficient and the number of real detections is smaller. We also remark that filtering with a Gaussian window is not the optimal procedure.

The extraction of the sources identified at a certain scale is proposed to be done directly on real space. At the location of the source $`\stackrel{}{x_o}`$ one subtracts the function $`A\tau (|\stackrel{}{x}\stackrel{}{x_o}|)`$, being $`\tau `$ the given profile.

All the calculations assume that the overlapping of nearby sources is negligible and the contribution of the sources to the background is also negligible. This is a very good approximation in many cases of interest at least above a certain flux level.

We remark the advantages of our method over traditional ones (DAOfind, SExtractor): we do not need to assume a smooth background and/or some filters (e.g. Gaussian) in order to detect the sources. For some astrophysical cases (CMB) the background can be complex so a smooth surface could be not a reasonable assumption. However, we need to assume the profile of the source and statistical properties of the background in order to get the optimal filter. The main advantage of our method is the amplification effect (gain) going to pseudo-filter space.

Generalization of these studies, considering different kind of sources (including non-centrally symmetric ones), are now being undertaken. The applications of this type of methodology is without any doubt of interest not only for Astrophysics/Cosmology but for other sciences.

This work has been supported by the Comision Conjunta Hispano-Norteamericana de Cooperación Científica y Tecnológica ref. 98138, Spanish DGESIC Project no. PB98-0531-c02-01, FEDER project no. 1FD97-1769-c04-01 and the EEC project INTAS-OPEN-97-1992 for partial financial support. J. L. S. acknowledges partial financial support from Spanish MEC and thanks CfPA and Berkeley Astronomy Dept. hospitality during year 1999. D. H. acknowledges a Spanish M.E.C. PhD. scholarship.