text:
# Novel infrared generation in low-dimensional semiconductor heterostructures via quantum coherence

## Abstract

A new scheme for infrared generation without population inversion between subbands in quantum-well and quantum-dot lasers is presented and documented by detailed calculations. The scheme is based on the simultaneous generation at three frequencies: optical lasing at the two interband transitions which take place simultaneously, in the same active region, and serve as the coherent drive for the IR field. This mechanism for frequency down-conversion does not rely upon any ad hoc assumptions of long-lived coherences in the semiconductor active medium. And it should work efficiently at room temperature with injection current pumping. For optimized waveguide and cavity parameters, the intrinsic efficiency of the down-conversion process can reach the limiting quantum value corresponding to one infrared photon per one optical photon. Due to the parametric nature of IR generation, the proposed inversionless scheme is especially promising for long-wavelength (far- infrared) operation.

-

Low-dimensional semiconductor heterostructures are almost ideally suited for generation in the mid- to far-infrared range (denoted below as IR for brevity), because the spacing between levels of dimensional quantization can be conveniently manipulated over the region from several to hundreds of microns, and injection pumping is possible. There is however a major problem: strong non-resonant losses of the IR field due to free-carrier absorption and diffraction, which become increasingly important at longer far-IR wavelengths. Due to very short lifetime of excited states, it is difficult to maintain a large enough population inversion and high gain at the intersubband transitions, necessary to overcome losses. There were many suggestions to solve this problem by rapid depletion of the lower lasing state using, e.g., the resonant tunneling to adjacent semiconductor layers or transition to yet lower subbands due to phonon emission , or even stimulated interband recombination ; see  for recent reviews. The successful culmination of these studies is the realization of quantum cascade lasers , in which the lower lasing state is depopulated either by tunneling in the superlattice or due to transition to lower-lying levels separated from the lasing state by nearly the energy of a LO-phonon; see e.g.  and references therein.

We here put forward another possibility, allowing us to achieve IR generation without population inversion at the intraband transition. This becomes possible with the aid of laser fields simultaneously generated at the interband transitions (called optical fields for brevity), which serve as the coherent drive for the frequency down- conversion to the IR. Employing self-generated optical lasing fields provides the possibility of injection current pumping and also removes the problems associated with external drive (beam overlap, drive absorption, spatial inhomogeneity), which were inherent in previous works on parametric down-conversion in semiconductors, see  for a review. The second important feature is the great enhancement of nonlinear wave mixing near resonance with intersubband transitions. The third feature is the possibility of canceling the resonance one-photon absorption for the generated IR field due to coherence effects provided by self-generated driving optical fields . The processes incorporating these three features constitute an important field of research with a variety of physical effects and promising applications. We consider just one example of such processes of IR generation and nonlinear mixing with self-generated optical fields in semiconductor heterostructures.

To avoid misunderstanding, we note that the mechanism of IR generation discussed in this paper is different from the approach in which the resonant tunneling and Fano-type interference are used to establish a large coherence at intersubband transitions . These quantum interference ideas imply usually the presence of a long- lived coherence at the intraband transitions. Our approach does not require long dephasing times. We here focus on nonlinear wave mixing phenomenon, which is greatly enhanced near the resonance with intersubband transition. Note that this process has some common features with recently observed generation of coherent IR emission in rubidium vapor in four-wave mixing experiments . Also, coherent microwave generation at the difference frequency under the action of two resonant external optical fields was observed in cesium vapor .

Generic three-level scheme. As the simplest case, consider the situation when only three levels of dimensional quantization are involved in generation: one (lowest-lying) heavy-hole level, and two electron levels; see Fig. 1. Of course, this scheme also describes the situation when there are two hole levels and one electron level involved. We need all three transitions to be allowed by selection rules. In a quantum well (QW) this will generally require using asymmetric structures, e.g. rectangular well with different barrier heights. For example, in a Al0.3Ga0.7As/GaAs/Al0.2Ga0.8As QW parametric IR generation is possible either between two electron subbands $`1e`$ and $`2e`$ separated by 98 meV ($`\lambda 13`$ $`\mu `$m) or between two lowest heavy-hole subbands ($`\lambda 60`$ $`\mu `$m).

Symmetric QWs can also be employed. e.g. in the case of a strong coupling between different subbands of heavy and light holes. A typical example is Al0.3Ga0.7As/GaAs/Al0.3Ga0.7As QW in which the second subband of heavy holes $`2hh`$ and the first subband of light holes $`1lh`$ happen to be very close to each other (within homogeneous linewidth) in the $`\mathrm{\Gamma }`$-point for a wide range of thicknesses $`58`$ nm, and therefore are strongly mixed. In this case two optical fields correspond to $`1e1hh`$ and $`1e1lh`$ transitions. And IR field is generated via the $`1hh2hh`$ transition.

A third configuration of interest is the quantum dot (QD). For example in a self-assembled InAs/GaAs quantum dot (QD) the three-level scheme can be easily realized with all three transitions allowed .

When the injection current density reaches the threshold value $`j_{\mathrm{th}}`$, optical generation starts due to recombination transitions between ground electron and hole states. Upon increasing the pumping current, optical generation can start also from excited states and the laser can be completely switched to lasing from the excited-state which has higher maximum gain due to a larger density of states. The effect of excited-state lasing was studied both in QW and QD lasers . It was found that with optimized laser parameters the region of simultaneous ground-state and excited- state lasing can be around $`j2j_{\mathrm{th}}`$ . In order to have the region of two-wavelength lasing sufficiently broad ($`\mathrm{\Delta }j(0.10.2)j_{\mathrm{th}}`$), gains for the two wavelengths should be close to each other.

The presence of one or two strong optical driving fields in the cavity gives rise to a rich variety of resonant coupling mechanisms by which the IR field can be produced. Here we will concentrate on one such scheme in which the two coherent optical fields having frequencies $`\omega _1`$ and $`\omega _2`$ excite the induce electronic oscillation at the difference frequency $`\omega _2\omega _1`$. It is important to note that the coherent IR polarization is parametrically excited independent of the sign of population difference at the IR transition.

The resulting output intensity of IR radiation depends on the coupling coefficient between the IR polarization and the cavity modes. It is clear that polarization wave has longitudinal wavenumber $`k_x`$ equal to the difference $`k_{2x}k_{1x}`$ of longitudinal wavenumbers of the two optical fields. Therefore, only the mode having the above wavenumber is efficiently excited. The field intensity is maximized when the frequency of this mode is equal to the difference frequency of optical fields. This requires special waveguide design since refractive indices of bulk semiconductor materials for optical and IR frequencies are different. For far-IR generation, there is more flexibility due to efficient manipulation of the refractive index by a slight doping.

Basic model. To quantify the above ideas, we have calculated the excited IR polarization and field by solving the coupled electronic density-matrix equations and electromagnetic Maxwell equations for the three fields, assuming steady-state. It is convenient to expand all fields in an orthonormal set of cavity modes $`𝐅_\lambda `$ and to introduce slowly varying complex amplitudes of fields and polarizations. For example, for the IR field at the $`32`$ transition we can write

$$𝐄(𝐫,t)=\underset{\lambda }{}\frac{1}{2}(t)𝐅_\lambda (𝐫)\mathrm{exp}(i\omega t)+\mathrm{c}.\mathrm{c}.$$
(1)
We will assume that the mode has a simple $`\mathrm{exp}(\pm ik_xx)`$ dependence in the propagation direction $`x`$ with the refractive index $`\mu =k_xc/\omega `$ and the transverse structure defined by a specific waveguide.

After introducing complex Rabi frequency $`e(t)=d(t)/2\mathrm{}`$, the wave equation can be written as

$$\frac{de}{dt}+(\kappa +i(\omega _c\omega ))e=\frac{2\pi i\omega d^2N}{\mathrm{}\mu ^2}_{V_c}\underset{j}{}\sigma _{32}^jF_\lambda (𝐫)d𝐫.$$
(2)
Here $`d`$ is the dipole moment of the IR transition, $`N`$ is the total volume density of electron states in the active region, $`\omega _{32}`$ the central frequency of $`32`$ transition, $`\omega _c`$ the frequency of the IR cavity mode with given $`k_x`$, $`\kappa `$ the cavity losses, $`V_c`$ the cavity volume. The variable $`\sigma _{32}^j`$ is the slowly varying amplitude of the element $`\rho _{32}^j`$ of density matrix, index $`j`$ labels different electron states contributing to the inhomogeneously broadened line. For a system of QDs, $`j`$ is simply the dot label. In QWs, index $`j`$ labels different $`𝐤_{}`$ states with respect to longitudinal quasimomentum. Only $`𝐤_{}`$-conserving transitions are considered, where $`𝐤_{}`$ lies in the plane of the wells ($`xy`$-plane).

The same representation is assumed for the two optical fields $`𝐄_{1,2}`$, with off-diagonal density matrix elements $`\sigma _{21}^j`$, $`\sigma _{31}^j`$ on the right-hand side, and corresponding parameters $`\mu `$, $`d`$, $`\omega `$, $`\omega _c`$, $`\kappa `$, $`k_x`$, $`g`$ having index 1 or 2. Note that the optical fields can have arbitrary polarization with respect to the parametrically excited IR field. In particular, for a QW laser the IR field should be $`z`$-polarized in the general case, while the optical fields are preferentially $`y`$-polarized.

Expressions for $`\sigma _{ik}`$ and population differences $`n_{ik}=\rho _{ii}\rho _{kk}`$, $`i,k=1,2,3`$, are found from the density matrix equations with phenomenological rates of relaxation and pumping. For simplicity, we will assume bipolar injection with equal injection rates of electrons to level 3 and holes to level 1, so the total particle density is conserved. In the limit $`|e|<\gamma _{ik}`$ the resulting amplitude of the IR field is found to be

$$e\frac{ig^2e_1^{}e_2}{\kappa }\underset{j}{}\left(\frac{n_{12}(\nu _j)}{\mathrm{\Gamma }_{21}^{}\stackrel{~}{\mathrm{\Gamma }}_{32}}+\frac{n_{13}(\nu _j)}{\mathrm{\Gamma }_{31}\stackrel{~}{\mathrm{\Gamma }}_{32}}\right),$$
(3)
where $`g^2=2\pi \omega d^2NG/(\mathrm{}\mu ^2)`$,

$`\mathrm{\Gamma }_{21}`$ $`=`$ $`\gamma _{21}+i(\omega _{21}+\nu _j\omega _1),`$ (4)
$`\mathrm{\Gamma }_{31}`$ $`=`$ $`\gamma _{31}+i(\omega _{31}+\nu _j\omega _2),`$ (5)
$`\mathrm{\Gamma }_{32}`$ $`=`$ $`\gamma _{32}+i(\omega _{32}+\nu _j\omega _2+\omega _1),`$ (6)
$`\stackrel{~}{\mathrm{\Gamma }}_{32}`$ $`=`$ $`\mathrm{\Gamma }_{32}+|e_1|^2/\mathrm{\Gamma }_{31}+|e_2|^2/\mathrm{\Gamma }_{21}^{}.`$ (7)
Here $`G`$ is the optical confinement factor for the IR field, $`\nu _j`$ is the difference between the transition frequency for a given $`j`$th state and the central frequency $`\omega _{21}`$, $`\omega _{31}`$, or $`\omega _{32}`$.

The resonance $`\omega _c=\omega _2\omega _1`$ with a given cavity mode having the wavenumber $`k_xk_{2x}k_{1x}`$ is assumed.

The value of the dipole moment at the IR transition is typically $`d(13)`$ nm, while it is $`(0.31)`$ nm for the optical transitions . The relaxation rates $`\gamma _{ik}`$ of both optical and IR polarizations are of order 5-10 meV in QW lasers at room temperature (i.e., the relaxation time $`0.1`$ ps)  and can be several times lower in self- assembled QDs (1 ps). In the IR range the cavity losses are mainly due to free carrier absorption. At the high carrier densities $`N10^{18}`$ cm-3 necessary for the excited-state lasing in QW’s the material losses in a bulk active medium are of order 100 cm-1 at $`\lambda 6`$ $`\mu `$m, and grow as $`\lambda ^2`$ or $`\lambda ^3`$ depending on the dominant scattering mechanism . Thick cladding layers have smaller doping density of order $`4\times 10^{16}10^{17}`$ cm-3, but can contribute significantly to the losses due to the large overlap factor $`G1`$. Large losses are a major problem in all proposed far-IR lasing schemes involving free carriers in semiconductors. Note, however, that increasing losses is not a principal limitation in the present parametric scheme since they do not prohibit generation itself, but rather decrease the IR field intensity. This possibility of far-infrared operation is an important advantage of the proposed mechanism. Note also that in QD lasers the excited state lasing can be achieved at lower carrier densities due to the state filling effect, and the intrinsic losses in active medium are lower.

Homogeneous broadening. This case is relevant for QWs at low temperatures. For room temperature, the quasi-Fermi energy for electrons is expected to be larger than the homogeneous bandwidth of 7-10 meV, and inhomogeneous broadening cannot be neglected. As for the self-assembled QDs, present-day structures have a large inhomogeneous broadening $`20`$ meV associated with the spread of dot sizes, which is definitely much larger than the homogeneous linewidth.

The resulting expression for the IR field is

$$|e|\frac{|e_1||e_2|}{\gamma _{32}}\left(\frac{\omega }{\omega _1}\frac{d^2}{d_1^2}\frac{\kappa _1}{G_1}\frac{G}{\kappa }+\frac{\omega }{\omega _2}\frac{d^2}{d_2^2}\frac{\kappa _2}{G_2}\frac{G}{\kappa }\right),$$
(8)
where $`G`$, $`G_1`$ and $`G_2`$ are the IR and optical confinement factors. At a given wavelength, the crucial parameter in Eq. (8), which governs the efficiency of down-conversion, is $`\eta =(\kappa _{1,2}/G_{1,2})(G/\kappa )`$. The main source of IR losses is free-carrier absorption. For the optical fields the ratio $`\kappa _{1,2}/G_{1,2}`$ in which the material gain at the optical transition, which is of order $`10^3`$ cm-1 in QWs and $`10^4`$ cm-1 in QDs. Therefore, even for very high IR material losses of order $`10^310^4`$ cm-1, the $`\eta `$ parameter can still be close to unity. As we already mentioned, IR losses can be dominated by absorption in doped cladding layers. Let us take as an example $`2\kappa 150`$ cm-1 as measured in quantum cascade lasers at 17 $`\mu `$m wavelength. If $`2\kappa _1/G_11500`$ cm-1, then we have $`\eta 1`$ for $`G0.1`$. Evidently, the value of $`\eta `$ decreases rapidly with increasing wavelength due to growing $`\kappa `$ and decreasing $`G`$.

Inhomogeneous broadening. We will assume here that the inhomogeneous widths $`u_{ik}`$ of all transitions are much larger than all homogeneous bandwidths $`\gamma _{ik}`$ and Rabi frequencies $`|e_{1,2}|`$ of optical laser fields that are at exact resonance with the centers of inhomogeneously broadened lines: $`\omega _1=\omega _{21}`$ and $`\omega _2=\omega _{31}`$. We consider two different situations in which explicit formulas can be obtained: when the optical field intensities are (i) much smaller and (ii) much greater than the saturation values.

(i). In this case the populations have the spectral distributions as supported by pumping in the absence of generated fields (no spectral hole burning). When $`u_{ik}\gamma _{ik}`$, the precise shape of the inhomogeneous line is not important and simple expression for the IR field can be obtained:

$$|e|\frac{2|e_1||e_2|}{(\gamma _{32}+\gamma _{21})}\frac{u_{21}}{u_{32}}\frac{\omega }{\omega _1}\frac{d^2}{d_1^2}\frac{\kappa _1}{G_1}\frac{G}{\kappa },$$
(9)
or, in terms of intensities,

$$||^2|_1|^2\frac{|_2|^2}{|_2|_s^2}\left(\frac{2\gamma _{32}}{(\gamma _{32}+\gamma _{21})}\frac{d}{d_1}\frac{\omega }{\omega _1}\frac{u_{21}}{u_{32}}\frac{\kappa _1}{G_1}\frac{G}{\kappa }\right)^2.$$
(10)
Here $`|_2|_s^2=\mathrm{}^2\gamma _{32}^2/d_2^2`$. Equation (9) is similar to (8) except for the ratios of inhomogeneous linewidths that are expected to be of the order of unity. As we see, the ratio of the IR to optical intensity does not change much as compared with homogeneously broadened transitions, but, of course, the threshold conditions do change, so the value of injection current is now much larger.

(ii). In this case the population distributions have narrow spectral holes burned at the spectral positions of optical field modes. This can be relevant for QD lasers or powerful pulsed multiple QW lasers. To simplify calculations, let us assume that the absolute values of the Rabi frequencies of two optical fields are equal, $`|e_1|=|e_2|`$, all $`\gamma _{ik}`$ and $`r_i`$ are equal to the same value $`\gamma `$, and also $`u_{21}=u_{31}=u`$. After somewhat lengthy calculations, we arrive at a surprisingly simple result:

$$|e|\frac{|e_1||e_2|u}{\gamma u_{32}}\left|0.9\frac{\omega }{\omega _1}\frac{d^2}{d_1^2}\frac{\kappa _1}{G_1}\frac{G}{\kappa }0.1\frac{\omega }{\omega _2}\frac{d^2}{d_2^2}\frac{\kappa _2}{G_2}\frac{G}{\kappa }\right|.$$
(11)
The requirement $`|e|<\gamma `$ employed in the derivation of Eq. (11) means that $`|e||e_1|`$ in the asymptotics (10), since the latter was obtained in the limit $`|e_1|\gamma `$.

Expressions (8)-(11) predict the rapid growth of the IR intensity, proportional to the product of optical field intensities. Of course, this tendency holds only until the IR intensity reaches the saturation value. Above this value, the IR field begins to deplete the electron populations, and its growth becomes nonlinearly saturated. In the optimal case, the maximum efficiency of down-conversion approaches the limiting quantum value corresponding to one IR photon per one optical photon. For the mid-IR range $`510`$ $`\mu `$m the maximum IR power is $`1050`$ mW if we take the value of $`100`$ cm-1 for IR losses. Here we assumed the laser with parameter $`\eta 0.1`$; see the discussion after Eq. (8). Beyond the restsrahlen region of strong phonon dispersion ($`\lambda 50`$ $`\mu `$m for AlGaAs/GaAs structure) the expected IR power is $`1`$ mW due to rapidly growing losses. We note, however, that progress in fabricating low-loss IR waveguides for semiconductor lasers is impressive. For example, in  waveguide losses as low as 14 cm-1 at 9.5 $`\mu `$m wavelength were reported.

Our calculations demonstrate the generation of coherent IR emission at intersubband transitions due to nonlinear wave mixing in standard multiple QW or QD laser diodes. The prerequisite for this is simultaneous lasing at two optical wavelengths which provide the necessary drive fields. This mechanism does not require population inversion at the IR transition, and its threshold current is determined by the minimum injection current necessary for the excited-state lasing.

It is important to note that the proposed parametric scheme seems to be viable in the far-IR region ($`\lambda 10100`$ $`\mu `$m) using transitions between hole subbands in the valence band. Indeed, in conventional lasers the increase of IR losses with wavelength is a significant problem for generation since it makes it difficult to reach the lasing threshold. In our case, due to the parametric nature of IR generation, the increase in losses will only lead to a decrease in output intensity; not a complete loss of IR generation as is the case in the usual below threshold laser behavior. For a practical device, the goal is to maximize the factor $`(\kappa _{1,2}/G_{1,2})(G/\kappa )`$. Also, in the far-IR spectral range it is relatively easy to manipulate the refractive index of the IR mode and provide phase matching by only a modest doping of waveguide layers.

The authors gratefully acknowledge encouraging and helpful discussions with D. Depatie and R. Haden, A. Botero, Yu. Rostovtsev, and H. Taylor and thank the Office of Naval Research, the National Science Foundation, the Robert A. Welch Foundation and the Texas Advanced Technology Program for their support.


format: text

text:
# 1 Introduction

## 1 Introduction

Up to now a satisfactory understanding of baryon resonances in the $`N\omega `$ channel has not been achieved. A number of experiments on the photoproduction of $`\omega `$ mesons have displayed indications for resonant structures . However the data basis is not sufficiently broad to allow for an unique determination of the quantum numbers of the involved resonances. Thus – based mostly on quark models – a number of resonances in the mass range of approximately $`1.82`$ GeV have been proposed to show a sizeable coupling to the $`\omega `$ meson . On the other hand, no experimental evidence for the existence of baryon resonances in the $`N\omega `$ system has been reported from the reaction $`\pi N\omega N`$ .

An understanding of baryon resonances in the $`N\omega `$ channel is not only important for a description of scattering experiments, but might also be of interest as far as the properties of $`\omega `$ mesons in nuclear matter are concerned. This conjecture is guided from the experience with the $`\rho `$ meson, which has been shown to undergo significant modification in the nuclear medium due to the excitation of resonance-hole loops. In this scenario a central role is played by the subthreshold resonance $`D_{13}(1520)`$ , suggesting that especially information on the coupling strength of subthreshold resonances to the $`N\omega `$ system would be very valuable for a determination of the in-medium properties of the $`\omega `$ meson.

In this work we present estimates for the coupling strength of baryon resonances to the $`\omega N`$ channel $`f_{RN\omega }`$. Vector Meson Dominance (VMD) is utilized to relate $`f_{RN\omega }`$ to the isoscalar strength of the electromagnetic decay of the respective resonance, which is extracted from helicity amplitudes. The analysis covers all resonances for which helicity amplitudes are currently available , most of which are below the nominal $`\omega N`$ threshold. A coupling of subthreshold resonances to the $`N\omega `$ system has also been reported elsewhere . In  such a mechanism is found to be necessary for a satisfactory description of $`\pi N`$ scattering within a coupled-channel analysis, whereas the authors of  derive the coupling strength within a quark-model approach. We will compare our results with those of both works.

The paper is organized as follows. The details of the extraction of $`f_{RN\omega }`$ are presented in Sect. 2. In Sect. 3 we examine if our model gives reasonable results by comparing the VMD prediction for the isovector strength with fits to the hadronic decay width into the $`N\rho `$ channel. Sect. 4 is devoted to a presentation of the results for $`f_{RN\omega }`$ and a discussion of their compatibility with experimental information on the reactions $`\pi N\omega N`$ and $`\gamma N\omega N`$. Finally, in Sect. 5 the implications of excitations of resonance-hole loops for the in-medium properties of $`\omega `$ mesons are examined.

## 2 Determination of $`f_{RN\omega }`$

In this section we explain how we obtain an estimate for the magnitude of the coupling strength of a baryonic resonance to the $`N\omega `$ channel. Before going into the details of the calculation, we outline the basic idea of our approach.

Vector Meson Dominance (VMD) , a theory which describes photon-hadron interactions exclusively in terms of vector meson-hadron interactions, relates the hadronic coupling strength of resonances to vector mesons $`f_{RN\rho (\omega )}`$ and the isoscalar and isovector part of the photon-coupling, see Fig. 1:

$`f_{RN\omega }`$ $`=`$ $`g_sm_\omega {\displaystyle \frac{2g_\omega }{e}}`$
$`f_{RN\rho }`$ $`=`$ $`g_vm_\rho {\displaystyle \frac{2g_\rho }{e}}.`$
As values for $`g_\rho `$ and $`g_\omega `$ – the coupling strengths of $`\rho `$ and $`\omega `$ meson to the photon – we take $`g_\rho =2.5`$ and $`g_\omega =8.7`$ . The isoscalar and isovector coupling strength of the resonance to the $`N\gamma `$ system is given by $`g_s`$ and $`g_v`$, respectively, see Eq. 8.

Thus VMD gives access to both $`f_{RN\omega }`$ and $`f_{RN\rho }`$, if it is possible to obtain $`g_s`$ and $`g_v`$ from experimental data. In order to achieve this goal, it is clearly not sufficient to consider merely the total electromagnetic decay width of the resonance. Rather, the coupling has to be decomposed into an isoscalar and an isovector part, which is readily done by constructing suitable linear combinations of proton- and neutron-amplitudes. These amplitudes are known from experiment , thus allowing to deduce numerical values for $`g_s`$ and $`g_v`$.

In a nonrelativistic formulation, the coupling of baryon resonances to the photon-nucleon system can be formulated as :

$$_{RN\gamma }=\{\begin{array}{ccccc}\chi _R^{}ϵ_{ijk}\sigma ^iq^jϵ_\lambda ^k\chi _NI+h.c.\hfill & \text{for}\hfill & J^P=\frac{1}{2}^+\hfill & & \\ & & & & \\ \chi _R^{}ϵ_{ijk}S^iq^jϵ_\lambda ^k\chi _NI+h.c.\hfill & \text{for}\hfill & J^P=\frac{3}{2}^+\hfill & & \\ & & & & \\ \chi _R^{}R_{ij}F_\lambda ^{ij}\chi _NI+h.c.\hfill & \text{for}\hfill & J^P=\frac{5}{2}^+\hfill & & \\ & & & & \\ \chi _R^{}(\sigma _kϵ_\lambda ^k\omega \sigma _kq^kϵ_\lambda ^0)\chi _NI+h.c.\hfill & \text{for}\hfill & J^P=\frac{1}{2}^{}\hfill & & \\ & & & & \\ \chi _R^{}(S_kϵ_\lambda ^k\omega S_kq^kϵ_\lambda ^0)\chi _NI+h.c.\hfill & \text{for}\hfill & J^P=\frac{3}{2}^{}.\hfill & & \end{array}$$
(2)
Here $`\chi _R^{}`$ and $`\chi _N`$ denote resonance and nucleon spinors respectively. The $`\sigma ^i`$ are the Pauli matrices, $`S^i`$ denotes the spin-$`\frac{3}{2}`$ and $`R^{ij}`$ the spin-$`\frac{5}{2}`$ transition operator. $`q`$ is the c.m. momentum of the photon, $`\omega `$ its energy and $`ϵ_\lambda `$ its polarization vector. $`F_\lambda ^{ij}`$ denotes the spatial components of the electromagnetic field strength tensor $`F_\lambda ^{\mu \nu }`$ with $`F_\lambda ^{ij}=q^iϵ_\lambda ^jq^jϵ_\lambda ^i`$. By $`I`$ we denote the isospin part of the coupling which is given as:

$`I=\chi _R^I\left(g_s+g_v\left\{\begin{array}{cc}\sigma ^3\hfill & \\ S^3\hfill & \end{array}\right\}\right)\chi _N^I\text{for}\{\begin{array}{cc}I_R=\frac{1}{2}\hfill & \\ & \\ I_R=\frac{3}{2}\hfill & \end{array}.`$ (8)
The spinors $`\chi _R^I`$ and $`\chi _N^I`$ represent resonance and nucleon isospinors and $`I_R`$ denotes the isopin of the resonance. $`\sigma ^3`$ and $`S^3`$ are defined as above.

Note that the invariant amplitude – as derived from (2) and (8) – for the decay of a resonance into a proton, $`_p`$, is proportional to $`g_p=g_s+\alpha g_v`$, whereas for the neutron the amplitude $`_n`$ is proportional to $`g_n=g_s\alpha g_v`$. The factor $`\alpha `$ is defined as:

$$\alpha =\chi _R^I\left\{\begin{array}{cc}\sigma ^3\hfill & \\ S^3\hfill & \end{array}\right\}\chi _N^I.$$
(9)
For nucleon resonances it is $`1`$ and for $`\mathrm{\Delta }`$ resonances it is exactly the Clebsch-Gordan coefficient for the respective transition, in this case $`\alpha =\sqrt{\frac{2}{3}}`$ .

Thus, the linear combinations

$`_s`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left(_p+_n\right)`$
$`_v`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left(_p_n\right)\{\begin{array}{ccc}1& \text{for}& I_R=\frac{1}{2}\\ & & \\ \frac{1}{\alpha }& \text{for}& I_R=\frac{3}{2}\end{array}`$ (14)
are proportional to $`g_s`$ and $`g_v`$ respectively. $`_s`$ and $`_v`$ describe the $`RN\gamma `$ system in terms of the isospin of the photon rather than the isospin of the nucleon.

Experimental information on the decay amplitudes $`_p`$ and $`_n`$ exists in form the of measured helicity amplitudes $`A_{\frac{1}{2}}^{p/n}`$ and $`A_{\frac{3}{2}}^{p/n}`$. They describe the transition of the photon-proton and the photon-neutron system to a resonance state with $`j_z=\frac{1}{2}\text{or}\frac{3}{2}`$, where the quantization axis is defined along the direction of the photon momentum. From the helicity amplitudes the isoscalar and isovector part of the coupling are constructed in exactly the same way as shown in Eq. 2 for the Feynman amplitudes. The helicity amplitudes are known from experiment at the pole-mass of the resonance, allowing the determination of $`g_s`$ and $`g_v`$.

To this end we introduce the $`\gamma `$-width $`\mathrm{\Gamma }_{s/v}^\gamma `$, which – using the normalization from  – is defined in terms of the helicity amplitudes $`A_{s/v}`$ in the following way:

$$\mathrm{\Gamma }_{s/v}^\gamma (m_R)=\frac{q^2}{\pi }\frac{2m_N}{(2j_R+1)m_R}\left(|A_{\frac{1}{2}}^{s/v}|^2+|A_{\frac{3}{2}}^{s/v}|^2\right),$$
(15)
with $`j_R`$ and $`m_R`$ denoting spin and pole-mass of the resonance. Note that after the angular integration has been performed no interference terms between $`A^{\frac{1}{2}}`$ and $`A^{\frac{3}{2}}`$ appears in Eq. 15.

Clearly, $`\mathrm{\Gamma }_{s/v}^\gamma `$ can also be expressed using Feynman amplitudes:

$$\mathrm{\Gamma }_{s/v}^\gamma (k^2)=\frac{1}{(2j_R+1)}\frac{q}{8\pi k^2}|_{s/v}|^2,$$
(16)
where $`\sqrt{k^2}`$ is the invariant mass of the resonance and $`q`$ the momentum of the photon in the rest frame of the resonance. After summing over the photon polarizations, $`|_s|^2`$ assumes the following form:

$$|_s|^2=4m_Nm_R\kappa g_s^2q^2F(k^2).$$
(17)
The same relation holds for $`|_v|^2`$ with $`g_s`$ replaced by $`g_v`$. For the formfactor $`F(k^2)`$ at the $`RN\gamma `$ vertex we assume the following functional dependence :

$$F(k^2)=\frac{\mathrm{\Lambda }^4}{\mathrm{\Lambda }^4+(k^2m_R^2)^2}$$
(18)
with $`\mathrm{\Lambda }=1`$. If we consider on-resonance decays ($`\sqrt{k^2}=m_R`$) this formfactor is equal to $`1`$. The numerical factor $`\kappa `$ depends on the quantum numbers of the resonance. For isospin-$`\frac{1}{2}`$ resonances the values for $`\kappa `$ are listed in Table 2. In the case of $`\mathrm{\Delta }`$ resonances they need to be multiplied by a factor of $`F^2=\frac{2}{3}`$ due to isospin. The two expressions Eqs. 15 and 16 can now be equated allowing to solve for $`g_{s/v}`$:

$$g_{s/v}^2=\frac{4}{\kappa }\frac{|A_{\frac{1}{2}}^{s/v}|^2+|A_{\frac{3}{2}}^{s/v}|^2}{q}\{\begin{array}{ccc}1& \text{for}& I_R=\frac{1}{2}\\ & & \\ \frac{1}{F^2}& \text{for}& I_R=\frac{3}{2}\end{array}.$$
(19)
Again, in the case of $`\mathrm{\Delta }`$ resonances an additional factor $`\alpha ^2`$ needs to be introduced.

Thus it is possible to obtain $`g_{s/v}`$ from helicity amplitudes. The hadronic couplings $`f_{RN\omega (\rho )}`$ are then readily deduced from the VMD relation Eq. 2. The corresponding values are listed in Tables 1 and 2.

A modified version of VMD (called VMD1 in ) contains also a direct coupling between hadrons – in our case baryon resonances – and the photon. In particular, at the photon point $`q^2=0`$ the photon-resonance interaction proceeds without intermediate vector mesons, whereas the decay into a massive photon consists of both a direct coupling term and a term with an intermediate vector meson. In principle, it would be preferable to perform an analysis such as ours within the framework of the modified version of VMD by studying – currently unavailable – dilepton production data on the nucleon. This process involves massive photons and would, in combination with the data for real photons, allow one to disentangle the direct photon and vector meson contribution. As a further advantage, the coupling strength of the vector mesons then does not need to be extrapolated over the large mass range from the photon point to the on-shell mass of $`\rho `$ and $`\omega `$ meson.

## 3 Results for the $`\rho `$: How Reliable is VMD ?

In this section we study the applicability of VMD in the resonance region. This is done for the isovector part of the analysis, where the VMD predictions for $`f_{RN\rho }`$ can be compared with the measured resonance decay width into the $`N\rho `$ channel.

We consider all resonances with $`j_R<\frac{7}{2}`$, for which both the hadronic decay widths into $`N\rho `$ and the helicity amplitudes are known. This excludes a few light resonances like the $`S_{11}(1535)`$ and the $`P_{33}(1600)`$, for which only upper limits for the $`N\rho `$ branching ratio exist . For both the helicity amplitudes and for the partial $`N\rho `$ decay width we use different parameter sets in order to provide an estimate for the experimental uncertainties entering this analysis. The helicity amplitudes are taken from Arndt et al  and the PDG . Furthermore we use the values obtained in the work of , which is of particular interest as it represents the only available approach which describes simultaneously hadronic-and electromagnetic reactions over the full energy region of interest here. The $`\rho `$ decay widths are taken from the analysis of Manley et al  and the PDG . Furthermore, we study the additional model dependence introduced by different parameterizations of the spectral function of the $`\rho `$ meson.

The hadronic coupling constant $`f_{RN\rho }`$ is obtained via the following expression for the $`RN\rho `$ decay width (taken on the pole mass of the resonance):

$$\mathrm{\Gamma }_{RN\rho }=\frac{1}{8\pi m_R^2}\frac{1}{2j_R+1}\underset{2m_\pi }{\overset{m_Rm_N}{}}𝑑m\mathrm{\hspace{0.17em}2}mA_\rho (m)|_{RN\rho }|^2𝐪.$$
(20)
The decay amplitude $`_{RN\rho }`$ is obtained from the same Lagrangians as in the electromagnetic case and is proportional to $`f_{RN\rho }`$ . $`𝐪`$ is the three momentum of the $`\rho `$ meson in the rest frame of the resonance. Note that there is no formfactor in Eq. 20 because we consider the on-resonance decay.

The spectral function of the $`\rho `$ meson $`A_\rho (m)`$ is taken as:

$$A_\rho (m)=\frac{1}{\pi }\frac{m\mathrm{\Gamma }(m)}{(m^2m_\rho ^2)^2+m^2\mathrm{\Gamma }^2(m)}.$$
(21)
The decay into $`N\rho `$ of low-lying resonances, for example the $`D_{13}(1520)`$, proceeds mainly through the low-mass tail of $`A_\rho `$. The shape of the tail is quite sensitive to the parameterization of the $`\rho `$ decay width. To study the effect on $`f_{RN\rho }`$ we compare two different parameterizations of the $`\rho `$ decay width:

$`\mathrm{\Gamma }(m)`$ $`=`$ $`\left({\displaystyle \frac{m_\rho }{m}}\right)^2\mathrm{\Gamma }_0\left({\displaystyle \frac{𝐩_m}{𝐩_{m_\rho }}}\right)^3`$ (22)
and
$`\mathrm{\Gamma }(m)`$ $`=`$ $`\left({\displaystyle \frac{m_\rho }{m}}\right)\mathrm{\Gamma }_0\left({\displaystyle \frac{𝐩_m}{𝐩_{m_\rho }}}\right)^3{\displaystyle \frac{1+r^2𝐩_{m_\rho }^2}{1+r^2𝐩_m^2}}.`$ (23)
The first version follows from a one-loop approximation of the self energy of the $`\rho `$ meson , the second one is taken from the PDG . $`m`$ is the invariant mass of the $`\rho `$ meson, $`m_\rho =0.768`$ GeV its pole mass and $`\mathrm{\Gamma }_0=150`$ MeV its decay width on the pole mass. $`𝐩_{m_\rho }(𝐩_m)`$ denotes the 3-momentum of the pions measured in the rest frame of a decaying $`\rho `$ meson with mass $`m_\rho (m)`$. The range parameter $`r`$, appearing in the second expression has the numerical value $`r=5.3`$/GeV .

In Table 1 the results for the coupling constants and the corresponding error-bars - which are calculated from the error-bars assigned to the partial decay width and the helicity amplitude in the respective analysis \- are given. As a general tendency, our finding is that VMD works well within a factor of two. It tends to somewhat underestimate the coupling constants, leaving some room for an additional direct coupling of the resonance to the photon. This can be seen particularly well in the case of the $`D_{13}(1520)`$ and the $`F_{15}(1680)`$, which are the most prominent resonances in photon-nucleon reactions, and whose $`\rho `$ decay widths are also well under control. Note that in the case of the $`S_{31}(1620)`$ and the $`P_{13}(1720)`$ large discrepancies occur between the VMD predictions based on different sets of helicity amplitudes, thus reflecting the uncertainties inherent in this analysis. The helicity amplitudes from Arndt  and PDG  produce nearly identical results for the $`\rho `$ coupling of the $`P_{13}(1720)`$, which might seem suprising as the helicity amplitudes from both works are very different. However, most of these differences cancel out in the isovector channel after proton- and neutron- amplitudes have been subtracted from each other (see Eq. 2). In the isoscalar sector this cancellation does not occur and both sets lead to very different predictions for $`f_{RN\omega }`$, see Sect. 4.

For the $`P_{13}(1720)`$ and the $`F_{35}(1905)`$ VMD is off by an order of magnitude. We argue that this mismatch does not hint on a failure of VMD, but has to be attributed to the unsatisfactory experimental information on these two resonances. Comparing the helicity amplitudes of the $`P_{13}(1720)`$ obtained from different analyses  reveals that neither their sign nor their magnitude are under control at all. Also the partial decay width into the $`\rho N`$ channel is subject to large uncertainties, here the PDG and Manley differ by about a factor of three . Obviously, the experimental observation of this resonance is very complicated and its parameters might be sensitive to the details of the underlying theoretical model, such as the treatment of the non-resonant background. For a conclusive VMD analysis of these resonances it is therefore mandatory to enlarge the data base and to describe hadronic- and photoinduced reactions within one model. In  a VMD analysis for exactly these resonances was performed, leading to the general conclusion - in contradiction to this work - that VMD is not at all reliable in the resonance region.

Different parameterizations of the $`\rho `$ decay width mainly influence the coupling constants of low-lying resonances. Since the PDG parameterization of the $`\rho `$ decay (see Eq. 23) distributes more strength at invariant mass below the $`\rho `$ meson pole mass, it leads to smaller values for the coupling constant, thus leading to an improved agreement between VMD and the hadronic fits. This effect is most notable for the $`D_{13}(1520)`$ (see fourth and fifth column in Table 1).

Summarizing the results, we conclude that VMD can be applied in the resonance region on a phenomenological basis. We have shown that for the coupling constants $`f_{RN\rho }`$ it leads to an approximate agreement with values adjusted to the hadronic decay width $`\mathrm{\Gamma }_{RN\rho }`$ with a tendency to somewhat underestimate these values. Therefore, our approach should yield reasonable predictions – at least for the lower limits – of the unknown coupling constants $`f_{RN\omega }`$.

## 4 Results for the $`\omega `$: Strong Coupling of Subthreshold Resonances

In this section we present our results for the coupling constants $`f_{RN\omega }`$ following from the VMD analysis and discuss their compatibility with experimental information obtained from pion- and photon-induced $`\omega `$-production cross sections.

All nucleon resonances with $`j_R<\frac{7}{2}`$, for which helicity amplitudes have been extracted, are included. Thus we consider only one resonance above the $`N\omega `$ threshold in our analysis, namely the $`D_{13}(2080)`$. Three different parameter sets for the helicity amplitudes are used, the results from the analysis of Arndt et al  and Feuster et al  as well as the values listed in the PDG . The corresponding results for $`f_{RN\omega }`$ together with the error-bars are given in Table 2. We find a strong coupling to the $`N\omega `$ channel in the $`S_{11}`$, $`D_{13}`$ and $`F_{15}`$ partial waves; especially the $`S_{11}(1650)`$, the $`D_{13}(1520)`$ and the $`F_{15}(1680)`$ resonances show a sizeable coupling strength to this channel. In most cases the three parameter sets for the helicity amplitudes lead to similar results. As in the case of the $`\rho `$ meson, however, for a few resonaces ($`P_{11}(1440)`$, $`S_{11}(1650)`$ and $`P_{13}(1720)`$) the differences are more pronounced. In particular, for the $`P_{13}(1720)`$ the poor experimental situation does not allow a reliable extraction of the coupling strength, thus emphasizing the difficulties concerning this resonance which were already discussed in the previous section. The relatively large coupling constant for the $`D_{13}(1520)`$ might at first seem surprising as the helicity amplitudes suggest that its electromagnetic decay is mainly an isovector one, corresponding to a small value for $`g_s`$. However, $`f_{RN\omega }`$ is proportional to $`g_\omega `$ (see Eq. 2), which is about three times larger than $`g_\rho `$. As a result, the $`\rho `$ and $`\omega `$ coupling are of the same size, but with a much larger error-bar for the $`\omega `$ coupling.

It is noteworthy that the resonances with the largest coupling are well below the $`N\omega `$ threshold. Subthreshold resonances in the $`N\omega `$ channel are also reported elsewhere . In  it is shown that within a coupled-channel analysis the description of $`\pi N`$ scattering enforces resonant structures in the $`N\omega `$ channel. As in our work, strong contributions from the $`S_{11}`$ and $`D_{13}`$ partial waves are found. Quantitatively, both approaches yield quite different predictions for the coupling strength of the $`D_{13}(1520)`$, however. Whereas our VMD analysis predicts a coupling strength of about $`3`$ (see Table 2), in  a value of $`f_{RN\omega }6.5`$ is given. At the same time it is found in a quark model calculation  that the coupling constant of this resonance should be $`f_{RN\omega }2.6`$, which is surprisingly close to our result. However, from their quark model calculation the authors of  obtain initially coupling constants for hadronic Lagrangians which display very different energy dependences compared to ours. To be specific, in the case of the $`D_{13}(1520)`$ close to threshold $`\omega `$ and nucleon couple in a relative $`d`$-wave in their formalism, rather than forming an $`s`$-wave state as following from our Lagrangian, see Eq. 2 (which for example has also been employed in ). Thus a direct comparison of the results as done in  is possible only after simplifying assumptions, if at all.

A direct experimental observation of a resonant coupling of subthreshold resonances is very complex, since the resonances can add to the cross-sections only through the high-mass tails of their mass distribution and are therefore hard to disentangle from background effects. Sensitivity to the contribution of subthreshold resonances in $`\omega `$ production experiments can probably only be achieved by studying polarization observables. The current data do not permit such a project, however. As was pointed out in various previous works , an analysis of dilepton production on the nucleon in combination with VMD might provide further insight on this issue.

In spite of these difficulties, it is still rewarding to compare the VMD predictions with existing data on the reactions $`\pi ^{}p\omega n`$ and $`\gamma p\omega p`$. Thus it is possible to infer if the predicted coupling constants of the subthreshold resonances produce a satisfactory overall agreement with data above the threshold. Comparison with data allows also to discuss the results for the $`D_{13}(2080)`$, the only resonance in our analysis above threshold. Based on recent data from photoproduction experiments, which display a richer structure than expected within a simple meson-exchange model, the search for resonances in the $`N\omega `$ channel has so far concentrated on the mass range of approximately $`1.82`$ GeV. At the current stage the data do not allow to pin down the quantum numbers of the involved resonances, however, and a variety of predictions exist . The $`D_{13}(2080)`$ is not amongst them. In contrast to the photoproduction data, no experimental signature of a coupling of baryon resonances to the $`N\omega `$ channel above the threshold has been detected in pion induced reactions . We find for the $`D_{13}(2080)`$ an $`\omega `$ decay width of about $`70`$ MeV and argue that its contribution to both reactions is too small to be seen in experiment.

As a first approximation, we take the full production amplitude as an incoherent sum of Breit-Wigner type amplitudes, describing $`s`$-channel contributions:

$$\sigma _{ab}=\underset{R}{}\frac{2j_R+1}{(2j_N+1)\mathrm{\Omega }_a}\frac{\pi }{k_a^2}\frac{\mathrm{\Gamma }_a\mathrm{\Gamma }_b}{(E_aE_R)^2+\mathrm{\Gamma }_{tot}^2/4}.$$
(24)
Here $`a`$ stands for the incoming channel ($`\pi ^{}p`$ or $`\gamma p`$) and $`b`$ for the $`\omega N`$ channel. $`\mathrm{\Gamma }_{a(b)}`$ denote the respective partial decay width of the resonance and $`\mathrm{\Gamma }_{tot}`$ its total decay width, taking into account the $`\omega N`$ decay as well. The energy dependence of the photon and the $`\omega `$ decay amplitudes is obtained from the Lagrangians given in Eq. 2 and the formfactor of Eq. 18 with $`\mathrm{\Lambda }=1.0`$. Whenever possible, $`f_{RN\omega }`$ is obtained from the helicity amplitudes from Arndt et al, otherwise the PDG estimates are used. For the remaining channels we take the parameterization from Manley et al. In Eq. 24 $`\mathrm{\Omega }_a`$ is $`1`$ for the pion and $`2`$ for the photon. $`k_a`$ is the cm-momentum of the particles in channel $`a`$ and $`E_a`$ denotes their cm-energy. All resonances listed in Table 2 are included. We stress that the only free parameter is the cutoff parameter $`\mathrm{\Lambda }`$.

The results for both reactions are shown in Fig. 2 and Fig. 3. The data for the $`\pi `$-induced reaction are taken from  and we use the photoproduction data of . The cross-section for $`\pi ^{}p\omega n`$ is reproduced rather well, especially near threshold, whereas at high energies some strength is missing. This is in approximate agreement with the findings of the authors of , who are able to give a satisfactory description of this process around threshold in terms of subthreshold resonances. One can therefore conclude that the excitation of subthreshold resonances constitutes an essential ingredient to the production mechanism. This interpretation is confirmed by the fact that it is hard to understand this process in terms of non-resonant amplitudes only. Already the contribution from $`\rho _0`$ exchange overestimates the data; only by suppressing the amplitudes by the introduction of very restrictive form factors a rough agreement with experimental data can be achieved . Similar problems have been reported in  in a model that includes all Born terms. The contribution coming from the only resonance above threshold – the $`D_{13}(2080)`$ – is about $`0.1`$ mb, roughly $`10\%`$ of the total cross-section. This is certainly too small to produce a distinguishable resonant structure in the total cross-section.

On the other hand, the photoproduction data can not be saturated within the resonance model. In particular, it seems futile to look for the $`D_{13}(2080)`$ in this reaction. Since it is commonly assumed that the $`\pi ^0`$-exchange as invoked in  plays a key role in the photoproduction we – incoherently – added this contribution using the same parameters as given in . As shown in Fig. 3 the sum of both mechanisms yields a qualitative explaination of the data over the energy range under consideration.

It is not surprising that the resonance contribution is more likely to produce lower bounds for the cross-sections. We already discussed in Sect. 3, that the VMD analysis tends to underestimate the hadronic coupling constants. Also, at energies above the $`\omega N`$ threshold, the helicity amplitudes of only a few resonances are known and the experimentally observed cross-sections have to be explained by additional production mechanisms, such as the $`\pi `$-exchange in the photoproduction. Keeping this in mind, the predictions of the resonance model are in reasonable agreement with the data and can be viewed as a confirmation of the VMD analysis.

## 5 The $`\omega `$ Meson in Nuclear Matter

In the previous section we discussed the relevance of resonances in the $`\omega N`$ channel for an understanding of $`\omega `$ production data. Here we study up to which extent the existence of these resonant states affects the properties of $`\omega `$ mesons in nuclear matter. Of particular interest in this context is the fact that the VMD analysis predicts a strong coupling of subthreshold resonances, since in the case of the $`\rho `$ meson it is widely acknowledged that its in-medium properties are dominated by low lying resonances, especially the $`D_{13}(1520)`$.

The in-medium properties of the $`\omega `$ meson can be read off its spectral function, which is defined as:

$$A_\omega ^{med}(\omega ,𝐪)=\frac{1}{\pi }\frac{\text{Im}\mathrm{\Sigma }_{tot}(\omega ,𝐪)}{(\omega ^2𝐪^2m_\omega ^2\text{Re}\mathrm{\Sigma }_{med}(\omega ,𝐪))^2+\text{Im}\mathrm{\Sigma }_{tot}(\omega ,𝐪)^2}.$$
(25)
Energy and three-momentum of the $`\omega `$ meson are denoted by $`\omega `$ and $`𝐪`$. The total self energy of the $`\omega `$ meson $`\mathrm{\Sigma }_{tot}`$ is the sum of the vacuum and the in-medium self energies $`\mathrm{\Sigma }_{vac}`$ and $`\mathrm{\Sigma }_{med}`$. In this work we neglect $`\text{Re}\mathrm{\Sigma }_{vac}`$ and approximate $`\text{Im}\mathrm{\Sigma }_{vac}`$ with the dominating $`3\pi `$ decay width. We assume that this decay may only proceed via an intermediate $`\rho `$ meson. The coupling constant $`f_{\omega \rho \pi }`$ is adjusted to reproduce the vacuum decay width of $`8.41`$ MeV .

To lowest order in the nuclear density $`\rho _N`$, the in-medium self energy of the $`\omega `$ meson $`\mathrm{\Sigma }_{med}`$ in symmetric nuclear matter is given by:

$$\mathrm{\Sigma }_{med}=\rho _NT_{\omega N},$$
(26)
where $`T_{\omega N}`$ is the spin/isospin averaged $`\omega N`$ forward scattering amplitude. We approximate the scattering amplitude as a sum over all resonances which are discussed in this work and determine the resonance contributions in a non-relativistic approach.

For further details of the calculation the reader may consult our previous publication on the $`\rho `$ meson . Both calculations are carried out within the same framework, differences appearing only in some minor points: First, the explicit values for the coupling constants and isospin factors are different. More crucial, we now evaluate the self energy in the cm-frame rather than in the rest frame of the nucleon. We do so since it has been demonstrated in  that thus a much better approximation of a relativistic calculation is achieved. As explained in , this leads in general to a reduction of the results. Finally, at the $`\omega NR`$ vertex the same formfactor is taken as in Eq. 18, in contrast to the choice in .

In Fig. 4 we show the spectral function for both an $`\omega `$ meson at rest and moving with a momentum of $`0.4`$ GeV with respect to the nuclear medium. The most obvious observation is that our model predicts an in-medium $`\omega `$ meson which survives very well as a quasi-particle. At rest, its peak position is slightly shifted upwards by about $`20`$ MeV, due to level-repulsion. The in-medium width is roughly $`40`$ MeV as can be read off Fig. 5, where $`\mathrm{\Sigma }_{med}`$ is depicted for an $`\omega `$ at rest. The $`\omega N`$ forward scattering amplitude contains all elastic and inelastic channels. In the previous section we demonstrated that our model can give a good description of the important $`\pi N`$ channel close to threshold. This enhances our belief that the VMD approach yields reasonable results for the in-medium properties of the $`\omega `$ meson, especially for an $`\omega `$ at rest.

By construction, the in-medium self energy develops resonant structures in the subthreshold region. The most important contributions come from the $`D_{13}(1520)`$ and the $`S_{11}(1650)`$ resonances. In contrast to the case of the $`\rho `$ meson, these structures translate only into small bumps in the spectral function, reflecting \- as follows from the coupling constants - that the $`\omega `$ self energy is substantially smaller than that of the $`\rho `$ meson.

Our value for the in-medium broadening of the $`\omega `$ meson is in approximate agreement with that obtained in  and – in a dynamical simulation – in  and . Also in the coupled channel analysis of  an in-medium broadening of $`40`$ MeV for the $`\omega `$ meson is reported. However, in this work the authors find that the $`\omega `$ spectrum displays a richer structure in the subthreshold region, which is not surprising since their coupling constant for the $`D_{13}(1520)`$ is much larger.

At higher momenta the influence of resonance-hole excitations gets reduced, see Fig. 4. Also, the longitudinal and transverse channel display a very similar behaviour (we therefore only show the results for the transverse channel). This is different from the case of the $`\rho `$ meson, whose spectral function at large momenta receives a strong broadening in the transverse channel from the the coupling of $`p`$-wave resonances, for example the $`P_{13}(1720)`$, whilst the longitudinal mode is only slightly modified. In the case of the $`\omega `$ meson however, only the $`F_{15}(1680)`$ shows a sizeable coupling strength.

We have already mentioned that the $`\omega `$ decay in vacuum proceeds mainly via an intermediate $`\rho `$ meson. Due to kinematics, only the low-mass tail of the $`\rho `$ mass spectrum is involved in this process. Consequently, if a lot of $`\rho `$ strength is shifted to smaller invariant masses in the nuclear medium – as is generally believed, see for example  – the $`\omega `$ meson will receive additional broadening. We have estimated this effect by replacing the vacuum spectral function of the $`\rho `$ with the in-medium one as given in  for the calculation of the $`\omega `$ decay. To lowest order in the density, this corresponds to scattering processes like $`\omega N\pi R`$, where $`R`$ stands for any baryon resonance included in the calculation in ; such processes have so far not been included in any calculation. We find an additional in-medium width for the $`\omega `$ of about $`40`$ MeV at the $`\omega `$ mass from this process. It has to be pointed out, however, that we did not take into account the corresponding mass-shift of the $`\omega `$ meson. Since the $`\omega `$ decay width varies rapidly as a function of its mass, this is a strong simplification.

## 6 Summary

We have calculated the coupling strength of baryon resonances to vector mesons in a VMD approach. The isoscalar and isovector strength of the photon coupling are extracted by means of helicity amplitudes. For the isovector part we compare the VMD approach with fits to the hadronic decay width of the resonance into the $`N\rho `$ channel. This allows to test the reliability of VMD in the resonance region. We find that VMD works well within a factor of two, with a tendency to underestimate the coupling strength.

In the $`\omega `$ channel we find a strong coupling of subthreshold resonances, especially the $`D_{13}(1520)`$ and the $`S_{11}(1650)`$. A direct experimental observation of these couplings is very complex and not possible within the current data base, which provides only an upper bound for the resonance contribution. We have calculated the total cross-section for the reactions $`\gamma N\omega N`$ and $`\pi N\omega N`$ with a simple Breit-Wigner ansatz and have shown that our coupling constants do not overestimate the data. In particular, the pion-induced reaction is very nicely described near threshold. We have argued that these results can be interpreted as a confirmation of the VMD model.

Guided by the observation, that the $`\rho `$ meson properties in nuclear matter are strongly modified from the excitation of resonance-hole loops, we estimate these effects for the $`\omega `$ meson. Since on average the $`\omega `$ couples weaker to baryon resonances than the $`\rho `$, the effects are not as pronounced. In contrast to the $`\rho `$ meson, the $`\omega `$ survives very well as a quasi-particle. For an $`\omega `$ at rest, we find an upward mass-shift of about $`20`$ MeV and an in-medium broadening of about $`40`$ MeV. At higher momenta the medium modifications are further reduced. By taking into account also the possibility that in nuclear matter the $`\omega `$ meson might decay into a modified $`\rho `$ meson, we calculate an additional broadening of the $`\omega `$ of about $`40`$ MeV.

## 7 Acknowledgements

The authors gratefully acknowledge many stimulating discussions with V. Koch and the hospitality of the Nuclear Science Division of the LBNL, Berkeley, where parts of this work were completed. This stay was made possible through support by the NATO science program, NATO-collaborative research grant number 970102.


text:
# 1 A plot of the ratio 𝑅(𝑠) versus √𝑠 in the most sensitive √𝑠 interval. Up to √𝑠=2.125 GeV we show by continuous lines the upper and lower bounds of the sum of the exclusive channels. However in the “exclusive” analysis we only use these data up to √𝑠=1.9 GeV. The value quoted in Table 1 corresponds to the average of the bounds. Above √𝑠=1.46 GeV we show the inclusive measurements of 𝑅(𝑠), together with the interpolation (dashed curve) used in (). For 2.8<√𝑠<3.74 perturbative QCD is used to evaluate 𝑅; the central QCD prediction is shown by the continuous curve. In the region 1.46<√𝑠<1.9 we used, in turn, the inclusive and exclusive data to evaluate Δ𝛼_had, see Table 1.

## Acknowledgements

We thank Martin Grünewald, Andrei Kataev, Thomas Teubner and Zhengguo Zhao for informative discussions. One of us (MGR) thanks the Royal Society for support.


text:
# A Probabilistic Quantification of Galaxy Cluster Membership

## 1 Introduction

In the wake of the Hubble deep field campaigns (Williams et al. , 1996), and with the growth of large photometric surveys (e.g., DPOSS, SDSS, etc.), the techniques involved in calculating photometric redshifts have become increasingly more sophisticated and precise (see Weymann et al. , 1999, for a recent synopsis of the field). Accuracies of $`\sigma _z=0.05`$ can be routinely achieved (Brunner et al. , 1999). Therefore, it is now possible to use photometric redshifts to answer specific scientific questions (e.g.,  Connolly et al. , 1998; Brunner et al. , 2000). In light of this, we have begun a program to use accurate photometric redshifts to study the galaxy populations in high-redshift clusters of galaxies. For this study, we have examined the most well-studied sample of optically-selected, high-redshift clusters to date – that of Oke, Postman & Lubin (1998). Over the past few years, Oke et al. (1998) have undertaken a detailed photometric, spectroscopic, and morphological survey of nine candidate clusters at $`z0.6`$. The survey consists of deep $`BVRIK^{^{}}`$ photometry, over 130 low-resolution Keck spectra per cluster field, and high-spatial-resolution WFPC2 imagery from the Hubble Space Telescope (HST).

Previous work on similar galaxy cluster samples at various redshifts has relied on statistical corrections in order to remove contamination due to field galaxies when spectroscopic identifications are not available for all galaxies in the cluster field (e.g.,  Aragón-Salamanca et al. , 1993; Lubin, 1996; Smail et al. , 1997; Dressler et al. , 1997; Couch et al. , 1998; Stanford et al. , 1998; Lubin et al. , 1998). These types of corrections, however, are fraught with potential systematic effects, including cosmological variance, weak lensing, and sample selection criteria. While these effects will only serve to increase the noise of any measurement made over a statistically complete sample, they will systematically affect any measurements made on an individual cluster (e.g., the Morphological fraction or the Morphology-Density relationship), introducing a bias into the actual measurement process. Furthermore, these line-of-sight contamination corrections apply to the entire galaxy sample and do not provide any information for specific galaxies. As a result, we have decided to explore the role photometric redshifts could play in improving the reliability of identifying galaxy cluster members (Kodama et al. , 1999; Lubin & Brunner, 1999).

In this paper, we first describe in Section 2 the details of the data observations, reduction, and catalog generation, emphasizing relevant points which are important for the estimation of photometric redshifts. In Section 3 we describe the details of our redshift estimation technique. Next, we detail the cluster galaxy identification procedure in Section 4. We conclude with a discussion of our results and future applications.

## 2 Data

For the photometric redshift analysis, we have used the photometric and spectroscopic data from the original cluster survey of Oke et al. (1998). This survey consists of nine candidate clusters at redshifts of $`z0.6`$ : Cl 0023+0423, Cl 0231+0048, Cl 0943+4804, Cl 1324+3011, Cl 1325+3009, Cl 1604+4304, Cl 1604+4321, Cl 1607+4109, and Cl 2157+0347. Oke et al. (1998) have determined that six of the nine candidate clusters are indeed real density enhancements, consistent with that of a cluster of galaxies. The observational data relevant for the work described in this paper consists of $`B`$, $`V`$, $`R`$, $`I`$, and $`K^{^{}}`$ band imaging and spectroscopic observations of selected galaxies in each field.

### 2.1 Observations and Data Reduction

All of the ground-based optical observations, both broad band and spectroscopic, were taken with the Low Resolution Imaging Spectrometer (LRIS; Oke et al. 1995) on either the Keck I or the Keck II 10-m telescopes. The infrared imaging was performed with the Infrared Imaging Camera (IRIM) on the Mayall 4-m telescope at Kitt Peak National Observatory. In this section, we briefly describe the observations and the data reduction; however, the reader is referred to Oke et al. (1998) for a complete account of these observations.

#### 2.1.1 Broad-Band Optical

The photometric survey was conducted in four broad band filters, $`BVRI`$, which closely match the Cousins system. The response curves of these filters are shown in Figure 1 of Oke et al. (1998). In imaging mode, LRIS covers a field-of-view of $`6\times 8`$ arcmin. The total exposure times in each filter are 3600, 2000, 1200, and 900 seconds for $`B`$, $`V`$, $`R`$, and $`I`$, respectively, and were chosen to give fairly uniform errors in the photometry. The $`B`$, $`V`$, and $`R`$ observations consisted of two equal exposures to allow for accurate cosmic ray rejection. The $`I`$ exposure time was broken into 3 equal exposures in order to cope with the cosmic rays, as well as to avoid approaching the CCD saturation level.

The LRIS imaging data were reduced in the standard fashion. All frames were bias-subtracted, and pixel-to-pixel sensitivity variations were removed using dome flats. Large-scale gradients were removed by dividing each frame by a normalized two-dimensional spline fit to the sky values in a sky-flat. The sky-flat was created by generating a median image from a stack of frames for each night and passband. For the $`I`$ band data, a two-dimensional fringe map was also created from the median filtered image by removing large-scale gradients within the median image. Fringing was removed from each $`I`$ band frame by subtracting a suitably scaled version of the fringe map.

Image registration for a given cluster field was performed by identifying approximately 10 unsaturated stars (detectable in all 4 passbands) to be used as astrometric reference points. The mean $`X`$ and $`Y`$ offsets of these stars in every frame taken of the cluster were computed relative to their locations in a fiducial $`B`$ band image. All image data for the cluster were shifted to match the $`B`$ band coordinate system using a flux conserving Lagrangian interpolation scheme to achieve registration at the sub-pixel level. Once all frames of a given cluster were registered to a common coordinate system, the independent exposures in each passband were co-added to produce the final four $`BVRI`$ images.

The Keck imaging has been calibrated to the standard Cousins-Bessell-Landolt (Cape) system through exposures of a number of Landolt standard star fields (Landolt 1992). For a circular aperture with a radius of $`3^{^{\prime \prime }}`$, the approximate limiting magnitudes are $`B=25.1`$, $`V=24.1`$, $`R=23.5`$, and $`I=21.7`$ for a 5-$`\sigma `$ detection.

#### 2.1.2 Spectroscopic

Multi-slit observations of the cluster field were made with LRIS in spectroscopic mode using a $`300\mathrm{g}\mathrm{mm}^1`$ grating blazed at 5000 Å. The chosen grating provided a dispersion of 2.35 Å per pixel and a spectral coverage of 5100 Å. The grating angle was set in order to provide coverage from approximately 4400 Å to 9500 Å in the first order. A GG495 glass filter was used to eliminate the overlapping second order spectrum; there is, therefore, no second order contamination below 9700 Å. In order to obtain the full wavelength range along the dispersion axis, the field-of-view of the spectral observations was reduced from that of the imaging mode to approximately $`2\times 8`$ arcmin.

Spectroscopic candidates were chosen from preliminary $`R`$ band imaging of each cluster field. All objects brighter than approximately $`R=23.5\pm 0.1`$ within the spectroscopic field-of-view of LRIS were included as candidates for slit-mask spectra. For each cluster field, six different slitmasks were made with approximately 30 objects per mask (including guide stars and duplicate observations). The exposure time for each mask was 1 hour. Flat-fielding and wavelength calibration were performed using internal flat-field and arc lamp exposures which were taken after each science exposure.

In practice, about 130 spectra were acquired per cluster field, with $`100`$ yielding measurable redshifts. A redshift measurement procedure which relies, in part, on visual inspection was used. The quality of the redshift identification was ranked with a number from 1 to 4, which roughly corresponds to the number of features used to identify the redshift. A quality of 4 means that the redshift is certain; a quality of 1 means that only one emission line was observed, and the redshift is only possible. Because the fore/background contamination rate is so large at our cluster redshifts, approximately 50 – 85% of the galaxies turn out to be field galaxies, rather than cluster members (Oke et al. 1998; Postman et al. 1998, 2000).

#### 2.1.3 Infrared

Deep infrared imaging of the cluster sample was taken with IRIM which contains a $`256\times 256`$ NICMOS3 HgCdTe array and has a resolution of 0.60 arcsecond per pixel on the 4-m telescope. The field-of-view ($`154\times 154`$ arcseconds) covers the central region of each cluster. In order to make the deepest possible observations over this region the entire LRIS field-of-view was not mosaiced. Observations were made in the $`2.2\mu K^{^{}}`$ band. Each central field was observed using a $`4\times 4`$ dither pattern with a stepsize of $`10^{^{\prime \prime }}`$ and a total extent of $`30^{^{\prime \prime }}`$. Each exposure had an effective exposure time of 1 minute, with 4 co-additions of individual, background–limited 15 second integrations. The total integration time on an individual cluster varied between 3 and 4.4 hours. Because the fields were not excessively crowded, we were able to use in-field dithering to create a global sky flat.

The $`K^{^{}}`$ cluster data were reduced using the Deep Infrared Mosaicing Software (DIMSUM), a publicly available package of IRAF scripts. This software generates not only a final stacked image but also a corresponding exposure image, where each pixel encodes the total number of seconds in the corresponding stacked image. This exposure image was appropriately scaled and used as a weight image by SExtractor during the catalog generation (see §2.2).

The data were linearized, trimmed to exclude masked columns and rows on the edges of the arrays, and dark–subtracted using dark frames of the same exposure length as the observations. All images of a given night were flattened by a super flat made from a series of dome flats taken during the previous day. As part of the DIMSUM procedure, sky subtraction was done by subtracting a scaled median of nine temporally adjacent exposures for each frame. A first-pass reduction was used to create an object mask for each frame. This mask was created from a fully stacked mosaic image. It therefore excluded not only the bright objects, but also those objects too faint to be detected in an individual exposure. In the second pass, the object mask was used to avoid object contamination of the sky flat in the production of sky frames. Final mosaicing of the images of each cluster were made with a replication of each pixel by a factor of 4 in both dimensions. This procedure conserves flux while eliminating the need for interpolation when the individual frames are co-aligned. A bad pixel mask was used to exclude bad pixels from the final summed images.

Absolute photometric transformations were derived from the observations of Persson et al. (1998) HST standard stars. Each standard star was observed every night in five separate array positions. The observations of the standard stars were reduced in a similar manner to the cluster data (see Oke et al. 1998). An approximate limiting magnitude of $`K^{^{}}=20`$ for a $`5\sigma `$ detection was reached in a standard aperture of radius $`3\stackrel{}{\mathrm{.}}0`$.

### 2.2 Catalog Generation

We used SExtractor version 2.1.5 (Bertin & Arnouts, 1996) to perform the source detection and photometry because this program is able to detect objects in one image and analyze the corresponding pixels in a separate image. Applied uniformly to multi-band image data (i.e. use the same detection image for all measurement images), this method will produce a matched aperture photometric catalog. We generated an optimal detection image from the $`B,V,R,\&I`$ images using a $`\chi ^2`$ process (Szalay et al. , 1998). Briefly, this process involves convolving each input image with a Gaussian kernel matched to the seeing. The convolved images were squared and normalized so that they have zero mean and unit variance. The four processed images (corresponding to the original $`B,V,R,\&I`$ images) were coadded, forming the $`\chi ^2`$ detection image.

In order to determine the optimal threshold parameters for source detection, we compared a histogram of the pixel distribution in the $`\chi ^2`$ image with a $`\chi ^2`$ function with four degrees of freedom. (This distribution theoretically corresponds to the background pixel distribution for our coadded $`\chi ^2`$ image.) By taking the difference between the two histograms (pixel minus theory), we generated the histogram of pixel values which were due to the objects that we were trying to detect. We defined the Bayesian detection threshold as the intersection of the sky (or theoretical prediction) and object pixel distributions (i.e. where the object pixel flux becomes dominant). To convert this empirical threshold for use with SExtractor, we scaled the $`\chi ^2`$ threshold (which is a flux per pixel value) into a surface brightness threshold (which is in magnitudes per square arcsecond).

The infrared images were cataloged separately due to their radically different field-of-view and pixel scale. The catalog was generated using SExtractor as before but in single image mode (i.e. the detection image was also the measurement image). In order to maximize the number of detections, we utilized the exposure image generated by DIMSUM during the data reduction (see §2.1.3) to provide a weight map in order to properly determine object parameters in an image with varying signal-to-noise (due to the different total integration times in different pixels).

All of our resultant analyses use the total magnitudes calculated by SExtractor. The optical magnitudes were, by design, matched apertures which reduces the scatter in any photometric redshift technique. Since the $`K^{^{}}`$ band images were not co-registered with the optical images, the aperture can be slightly different. However, since the total apertures were rather large in angular extant, any “aperture effects” on our resulting analysis should be small ($`0.02^m`$) compared to the more dominant effects of photometric errors.

### 2.3 Cross-identifications

The last step in preparing these data for analysis was to combine the optical imaging, near-infrared imaging, and the spectroscopic observations into a single dataset for each cluster field. Since the original spectroscopic targets were determined from a separate catalog, we first cross-matched the spectroscopic catalog with the optical catalog using a growing annulus technique where the angular distance $`\psi `$ was determined using the formula

$$\psi =\mathrm{arccos}(\mathrm{sin}(\varphi _s)\mathrm{sin}(\varphi _p)+\mathrm{cos}(\varphi _s)\mathrm{cos}(\varphi _p)\mathrm{cos}(\theta _s\theta _p))$$

where the angles ($`\theta `$ Declination, $`\varphi `$ Right Ascension) have been properly converted into radians, and the subscript, $`s`$, refers to spectroscopic target and the subscript, $`p`$, refers to photometric object. In all cases, the match-up occurred for very small angles because spectroscopic targets were selected from the optical images.

The cross-identification between the infrared and optical catalogs, however, was slightly more complicated, although it was done in an identical manner. This was due to the different telescopes and detectors used during the observations, and the corresponding geometric distortions present between the two catalogs. Furthermore, the variation of the spectral energy distributions of the sources in our sample between the $`I`$ and $`K^{^{}}`$ bands can be quite large, especially for galaxies at relatively high redshift. As a result, some objects in each catalog remained unmatched below the angular limit which is approximately the confusion beam width (approximately five arcseconds).

## 3 Empirical Photometric Redshift

The technique that we have developed for determining cluster membership requires redshift estimates for all galaxies in the field. Although the technique is independent of the actual redshift estimation technique, we have used the empirical photometric redshift technique (Brunner et al. , 1999) to generate redshift estimates for all galaxies in our sample. This is the result of the simplicity of the empirical technique in determining redshifts for a photometric sample of galaxies for which an adequate training set exists. Furthermore, the empirical technique, since it is based only on the data under investigation, is more robust than alternative techniques to uncertainties in the spectral energy distributions of galaxies. These alternative techniques, such as template-based methods, could also be used, as long as a realistic estimate of the error in the estimated redshift can be determined.

The empirical photometric redshift technique requires a calibration sample to determine the coefficients of the empirical relation. For this work, we have used the spectroscopic data which was acquired as part of the original observational campaign to determine the cluster velocity dispersions (Postman et al. , 1998, 2000). The reliability (i.e. the inverse of the intrinsic error) of a spectroscopic redshift is generally quite high, with quoted errors often measured in the tenths of a percent (i.e. $`\delta <0.001`$). In reality, however, a redshift is accurate only when the spectral identification is also accurate. As a result, we restricted the spectroscopic calibrators to include only objects with redshifts that are determined by at least two spectral lines, i.e. a quality of 2 or higher (for details, see Postman et al. , 1998). The next step was to restrict the sample to those objects which were below the $`10\%`$ photometric error limit, which did not have bad detection flags (e.g., object near edge of frame, incomplete aperture data), and which were below our high redshift cut-off of $`z=1.5`$. We also removed one additional galaxy which, having an $`E+A`$ spectrum (Dressler & Gunn, 1983, 1992; Zabludoff et al. , 1996), skewed the results of our calibration for normal galaxy spectra (which are the dominant population in our sample). The final spectroscopic calibration data consisted of 130 galaxies from five different cluster fields.

With the final calibration sample, we determined a fourth order polynomial fit to the four magnitudes $`B,R,I,\&K^{^{}}`$, and the spectroscopic redshift. Our sample of 130 calibrators is nearly twice the number of coefficients required (a fourth order polynomial in four variables has 70 coefficients); we are, therefore, properly constraining the degrees of freedom when determining the coefficients of our polynomial fit. From Figure 1, there are no obvious systematic variations in our calculated relationship, which has an intrinsic dispersion (the standard deviation of the residual differences for the calibrating sample) of $`\sigma _z=0.072`$. In addition, the distribution of the residuals (see Figure 2) are approximately Gaussian in shape.

### 3.1 Estimating Photometric Redshift Errors

In order to reliably determine cluster membership, we also need to have an estimate for the error in our photometric redshift calculation. While the measurement of the intrinsic dispersion in our relationship provides one estimate (which is commonly used throughout the relevant literature, e.g.,  Kodama et al.  1999), we have also developed an alternative, extrinsic estimate. Briefly, the intrinsic error is calculated directly from any available spectroscopic redshifts; that is, it is simply the dispersion in the derived relationship between spectroscopic and photometric redshifts. Extrinsic errors, on the other hand, are calculated in a Monte-Carlo fashion, accounting for all known measurement uncertainties.

Following Brunner et al.  (1999), we derive extrinsic redshift error estimates for all galaxies in our dataset using a Monte-Carlo technique to simulate our uncertainty in both the observables (i.e. flux measurements), as well as our sampling of the galaxy distribution in the $`B,R,I,\&K^{^{}}`$ four-dimensional flux space. Specifically, we use a bootstrap with replacement algorithm to randomly draw a new set of calibrating galaxies from the original training set. This algorithm allows for duplicate galaxy entries and is designed to emphasize any incompleteness in the sampling of the true distribution of galaxies in the four dimensional space $`B,R,I,\&K^{^{}}`$ by the calibration redshifts. Furthermore, the magnitudes of the newly derived calibrating sample were drawn from a Gaussian probability distribution function with the mean given by the measured magnitude and sigma by the magnitude error. This accounts for any uncertainties in the flux measurements.

This technique was used to generate one hundred different realizations of the photometric redshift relationship. For each different realization (all of which form our ensemble), a photometric redshift was calculated for every galaxy in the photometric catalog. As a result, this produces one hundred different redshift estimates for each galaxy. By applying a quantile cut, we can approximate the redshift distribution for each individual galaxy by a Gaussian. Therefore, the error in the photometric redshift for each object was defined by $`\sigma _z=(Q5Q2)/2.0`$, where $`Q5`$ and $`Q2`$ refer to the fifth and second quantiles of the distribution of each galaxy’s estimated redshift, respectively (see Brunner et al. , 1999). Together, our original photometric redshifts and redshift error estimates completely define the photometric redshift distribution for each individual galaxy.

Of the two types of redshift error estimates, the intrinsic error is by far the most commonly used (e.g.,  Kodama et al. , 1999), as it is the easiest to calculate. The extrinsic error, however, does a better job of quantifying the actual error in the redshift estimate of an individual galaxy, which can be both smaller or considerably larger than the intrinsically estimated global error (see Figure 15 in Brunner et al. , 1999). Fundamentally, the extrinsic method provides an error estimate separately for each individual galaxy, while the intrinsic method provides only a global estimate.

## 4 Cluster Membership

While the optimal technique for unambiguously determining cluster membership is to obtain spectroscopic redshifts for all sources in the field, this process can be prohibitively expensive in telescope time allocation, especially at relatively high redshifts ($`z>0.5`$). As a result, many previous studies (e.g.,  Aragón-Salamanca et al. , 1993; Lubin, 1996; Smail et al. , 1997; Dressler et al. , 1997; Couch et al. , 1998; Stanford et al. , 1998; Lubin et al. , 1998) have utilized a background contamination estimator in order to statistically quantify both the number of galaxies which belong to a given cluster and the morphological fraction of the cluster members. In addition to possible systematic effects which are implicit in using techniques of this type, background corrections provide no information on the likelihood of an individual galaxy being a member of a cluster.

Photometric redshifts, however, provide a relatively inexpensive redshift estimate; therefore, they can be used to estimate cluster redshifts (e.g.,  Gal et al. , 2000) and to determine likely cluster members (Kodama et al. , 1999; Lubin & Brunner, 1999). Previous techniques to select cluster members based on photometric redshifts (e.g.,  Kodama et al. , 1999) are based on a simplistic technique for defining galaxy cluster members. Basically, any galaxy which lies within a redshift shell centered on the galaxy cluster’s redshift (i.e. $`z_c\pm \delta _z`$) are identified as members. Because photometric redshifts have a large uncertainty (generally a factor of ten to fifty higher than spectroscopic redshifts), a large redshift shell (e.g., $`\delta _z0.05`$) is used to identify likely cluster members. While certainly useful, this approach does not provide a reliable means for identifying line-of-sight contaminating galaxies.

As a result, we have developed an alternative technique which relies on the probabilistic interpretation of a photometric redshift to determine cluster membership. We define the probability density function, $`\mathrm{\Phi }(z)`$, for an individual galaxy’s redshift to be a Gaussian probability distribution function with mean ($`\mu `$) given by the estimated photometric redshift and standard deviation ($`\sigma `$) defined by the estimated error in the photometric redshift.

$$\mathrm{\Phi }(z)=\frac{1}{\sigma \sqrt{2\pi }}e^{\left(\frac{(z\mu )^2}{2\sigma ^2}\right)}$$

Using this interpretation, we can calculate the probability that a galaxy has an actual redshift within a given redshift interval.

$`P(Cluster|z_c,\mathrm{\Delta }z,\mu ,\sigma )`$ $`=`$ $`N{\displaystyle _{z_c\frac{1}{2}\mathrm{\Delta }z}^{z_c+\frac{1}{2}\mathrm{\Delta }z}}\mathrm{\Phi }(z)𝑑z`$
$`=`$ $`N{\displaystyle _{z_c\frac{1}{2}\mathrm{\Delta }z}^{z_c+\frac{1}{2}\mathrm{\Delta }z}}{\displaystyle \frac{1}{\sigma \sqrt{2\pi }}}e^{\left(\frac{(z\mu )^2}{2\sigma ^2}\right)}𝑑z`$
$`=`$ $`{\displaystyle \frac{N}{2}}(\gamma [{\displaystyle \frac{1}{2}},z_H]\gamma [{\displaystyle \frac{1}{2}},z_L])`$
where $`N`$ is a suitable normalization factor, $`z_c`$ is the cluster redshift, $`\mathrm{\Delta }z`$ is the width in redshift space which you are sampling, the limits of integration are

$$z_H=\frac{(z_c+\frac{1}{2}\mathrm{\Delta }z)\mu }{\sqrt{2}\sigma },z_L=\frac{(z_c\frac{1}{2}\mathrm{\Delta }z)\mu }{\sqrt{2}\sigma },$$

and $`\gamma `$ is the incomplete gamma function. (Of course, $`\gamma [\frac{1}{2},z]`$ is also known as the error function, $`erf(z)`$).

Within this formalism, the only undeclared quantity is $`\sigma `$ or, alternatively, the uncertainty in the estimated redshift for a given galaxy. In our framework, this value can either be determined from the intrinsic error in our photometric redshift relation or in a separately calculated extrinsic error. As a result, this approach is strongly dependent on the actual technique used to estimate the photometric redshift error and corresponding redshift error. In this paper, we have applied this technique using empirically derived photometric redshifts since we have a well-defined calibration dataset of spectroscopic redshifts. Empirically defined photometric redshifts are much less sensitive to uncertainties in the shape and the evolution with redshift of the spectral energy distribution of galaxies than competing techniques such as template photometric redshifts (Brunner, 1997).

On the other hand, the method we have introduced can easily be adapted to work with alternative redshift estimation procedures, including template photometric redshift algorithms, by utilizing the derived dispersion estimate. The empirical-based techniques do have shortcomings; most notably, the requirement of a high quality training set and the fundamental limitation that any resultant analysis is limited to the same region of flux space delineated by the training set of calibration redshifts. In our case, the spectroscopic survey of Oke et al.  (1998) provides a wonderful training set, which completely samples the region in flux space which encompasses the vast majority of cluster galaxies. When this is not the case, template based photometric redshift estimation methods, which can be applied to data where there are no calibrators, are required. In this case, however, extreme care must be used to minimize any systematic effects, such as morphological variation of the redshift error estimates or template incompleteness.

By using the calibration galaxies as a training set, we can empirically determine the optimal probability value which provides the threshold for cluster membership. This can be tuned to either maximize the completeness or minimize contamination. In this paper, our sample is small enough that a simple visual threshold determination is sufficient; however, there is no reason that a more powerful maximum likelihood or regression analysis technique could not be employed, especially for larger datasets (e.g., the SDSS). Furthermore, additional information (such as evolutionary color tracks) can be used to select likely cluster members which can be used to provide a bootstrap estimator for the probability threshold.

As a demonstration of this technique, we have applied this procedure to the data described earlier in this paper. We consider two clusters from our survey (see Section 2 for more information), the first is Cl 0023+0423 at a redshift of $`z=0.84`$ (Postman et al. 1998; Lubin, Postman & Oke 1998; Lubin et al. 1998), and the second is Cl 1604+4321 at a redshift of $`z=0.92`$ (Postman et al. , 2000; Lubin et al. , 2000). These two clusters provide a useful demonstration because the first cluster is well sampled by our dataset, while the second cluster pushes the photometric limits of our data. Since the utility of this technique is so strongly dependent on the exact method used for estimating the redshift and redshift error for every galaxy in the field, we have intentionally relegated a more complete discussion of the application of this technique to our own data to a later, more targeted paper (Lubin et al. , 2000).

In Figures 3 and 4, we display the results of applying this technique to the spectroscopic targets in the Cl 0023+0423 cluster field, as well as the full photometric sample for this field. In both cases, we used the intrinsic error of the photometric redshift relationship as our redshift error. From the first panel, we clearly do quite well (i.e. $`100\%`$ completeness with no contamination). Of course, translating this success to the full sample is not guaranteed; however, as long as our calibration sample adequately samples the four dimensional flux space occupied by the galaxies in our full sample (i.e. we have a fair and unbiased training set), there is no reason to expect significant variation from these results.

In Figures 5 and 6, we show the results of applying this technique to the spectroscopic targets in the Cl 1604+4321 cluster field, as well as the entire photometric sample for this cluster. This time, however, we have used the extrinsic error which is calculated independently for each galaxy. By judiciously selecting the probability threshold, we can maximize the completeness (which is $`88\%`$), while minimizing the contamination (which is $`21\%`$). While we would clearly prefer full completeness and no contamination, this cluster demonstrates one of the strengths of our technique, namely the ability to empirically quantify both our completeness and contamination.

In fact, our technique predicts 37 cluster members within the field of view of our combined image (see Figure 6). If we use the predictions of deep number counts (e.g.,  Driver et al. , 1995; Smail et al. , 1995; Abraham et al. , 1996) to estimate the contamination due to foreground and background galaxies, we estimate there should be $`33\pm 14`$ cluster members in the image (note the large uncertainty in this statistical calculation), which is the exact number that we find after correcting for contamination and incompleteness \[i.e. $`37\times (121\%)/88\%=33`$\]. We have made this comparison with the statistical background correction simply in order to show that we are finding a reasonable number of cluster galaxies; however, unlike the background correction technique, we have explicitly identified the likely cluster members.

## 5 Discussion

In this paper, we have presented a novel technique for the probabilistic determination of cluster membership based solely on photometric redshift and redshift error estimates. Our technique can be tuned to either maximize the completeness of cluster member identification (i.e. identify all actual cluster members), or, alternatively, to minimize the contamination due to foreground and background galaxies in the field. Furthermore, this technique provides a robust estimation for both the completeness and contamination fractions of identified cluster members as a function of likelihood. A comparison between our new technique with the traditional technique of statistical background correction shows remarkable agreement, with the caveat that our technique actually identifies the galaxies which are cluster members.

The results in this paper are based on an empirical photometric redshift estimation, which we have used because (a) we have a reliable training set from our extensive spectroscopic survey, and (b) these photometric redshifts are less sensitive to uncertainties in spectral evolution which can affect other redshift estimation techniques. Alternative redshift estimation techniques, however, can easily be incorporated into the algorithm, as all that is needed is a redshift estimate and a corresponding redshift error estimate. In the case of template photometric redshifts, the redshift error could be determined from either the intrinsic dispersion in the photometric redshift relation, an interpolation over the $`\chi ^2`$ goodness-of-fit parameter, or even more appropriately through a combination of the $`\chi ^2`$ parameter and a Monte-Carlo bootstrap approach to mimic photometric uncertainty and template incompleteness. In this case, however, caution must be used to minimize any possible systematic effects which might depend on various fundamental assumptions, such as the spectral type.

The only free parameter in the entire technique is the probability threshold used to determine cluster membership. Optimally, as we have done in this paper, a training set of spectroscopically confirmed cluster members can be used to empirically set the threshold value. When this is not the case, alternative techniques, such as color-selection to identify the early-type sequence at the appropriate redshift, can be used to identify likely cluster candidates. On the other hand, the expected number of cluster candidates can be calculated from the statistical background correction arguments. This expected number of cluster members can be used as an independent estimator as we can set the probability threshold value to reproduce the estimated number of cluster galaxies.

Based on our work, we believe that photometric redshifts are an ideal way to study the galaxy populations in high-redshift clusters. For a large number of the tests of common cluster properties (e.g., the Morphology-Density relationship), this method eliminates the need for extensive spectroscopic surveys and the uncertainty of estimating the background contamination. As an additional example of the utility of this approach, consider the magnitude-limited, spectroscopic survey used in this analysis. We can approximate the total number of spectroscopic targets observed by multiplying the number of clusters used in this analysis (five) by the approximate number of targets per cluster field (130), implying that a total of 650 objects were observed spectroscopically. The redshift identification success rate was nearly $`77\%`$. Of these redshifts, $`50\%`$$`85\%`$ belong to non-cluster galaxies. These numbers imply that the spectroscopic survey was only $`12\%`$$`38\%`$ efficient in identifying cluster candidates, depending on the richness of the cluster.

In addition, our technique provides the ability to tune the cluster candidate identification to either maximize the completeness or, alternatively, to minimize the contamination. As a practical example, maximizing completeness is essential for any observational program which is interested in completely sampling the cluster (e.g., spectroscopic targeting). On the other hand, minimizing contamination is important for selecting a characteristic set of cluster members for a detailed study of cluster member properties (e.g., the Morphological fraction or the Morphology-Density relationship; Lubin et al. , 2000). Furthermore, if a training set of spectroscopically confirmed galaxies is available, we can actually characterize the percentages of both the completeness and contamination that we expect in our final cluster candidate catalog.

Like any photometric redshift based technique, our technique is limited in its efficacy by the quality of the photometric data used in the analysis. However, our technique can be applied using a variety of photometric redshift estimation methods. Therefore, we expect that, in the near future, we will be able to fully utilize the capabilities of this technique for cluster research.

We are very grateful to Bev Oke and Marc Postman for the enormous time and effort that they put into obtaining and reducing all of the data used in this paper. We also wish to thank the anonymous referee for a thorough reading of the original paper and many interesting comments. This research has made use of NASA’s Astrophysics Data System Abstract Service. LML is supported by NASA through Hubble Fellowship grant HF-01095.01-97A from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555.