# Oxygen Metallicity Determinations from Optical Emission Lines in Early-type
Galaxies

Alex E. Athey  The University of Texas at Austin  Applied Research
Laboratories  10000 Burnet Rd  Austin, TX 78758<EMAIL_ADDRESS>Joel
N. Bregman  University of Michigan  Department of Astronomy  500 Church St.
Ann Arbor, MI 48109-1090<EMAIL_ADDRESS>
###### Abstract

We measured the oxygen abundances of the warm (T$\sim 10^{4}K$) phase of gas
in seven early-type galaxies through long-slit observations. A template
spectra was constructed from galaxies void of warm gas and subtracted from the
emission-line galaxies, allowing for a clean measurement of the nebular lines.
The ratios of the emission lines are consistent with photoionization, which
likely originates from the UV flux of post-asymototic giant branch (PAGB)
stars. We employ H II region photoionization models to determine a mean oxygen
metallicity of $1.01\pm 0.50$ solar for the warm interstellar medium (ISM) in
this sample. This warm ISM $0.5$ to $1.5$ solar metallicity is consistent with
modern determinations of the metallicity in the hot (T$\sim 10^{6}-10^{7}K$)
ISM and the upper range of this warm ISM metallicity is consistent with
stellar population metallicity determinations. A solar metallicity of the warm
ISM favors an internal origin for the warm ISM such as AGB mass loss within
the galaxy.

Galaxies:Elliptical and Lenticular, Galaxies:ISM, Galaxies:Abundances

##  1 Introduction

###  1.1 Warm ISM Discovery

Early-type galaxies were once thought to contain very little gas in their
interstellar medium (ISM) (e.g. Mathews & Baker (1971); Bregman (1978); White
& Chevalier (1983)), but are now known to host massive amounts of extremely
hot (T$\sim 10^{6}-10^{7}K$) gas that emits in the X-rays (Forman et al.,
1979). The reason early-type galaxies were originally thought to be void of
gas is that photographic-plate surveys did not reveal the warm (T$\sim
10^{4}K$) ISM that is abundant in spiral galaxies (e.g. Mayall (1936); Sandage
(1957); Mayall (1958)). Only with the advent of large telescopes and efficient
detectors was the warm ISM detected, although the masses are an order of
magnitude lower than that of spirals. Spectroscopic surveys in the mid-1980’s
determined warm ISM detection frequencies in the cores of early-type galaxies
of 30%-50% (Caldwell, 1984; Phillips et al., 1986; Veron-Cetty & Veron, 1986).
Subsequent narrow-band $H\alpha$ imaging studies revealed extended emission in
general agreement with the the stellar morphology and detection rates
consistent with the 50% rate determined from the early spectroscopic surveys
(Demoulin-Ulrich et al., 1984; Kim, 1989).

###  1.2 Key Warm ISM Studies

There are three large sample, multi-wavelength studies of early-type galaxies
that have provided a foundation for our understanding of the cold and warm
ISM. Trinchieri & di Serego Alighieri (1991, hereafter TdSA) conducted one of
the first quantitative early-type galaxy emission-line studies with $H\alpha$
imaging and long-slit optical spectroscopy of 13 early-type galaxies known to
have hot X-ray halos. (Goudfrooij et al., 1994a, b, c; Goudfrooij & de Jong,
1995, hereafter G94-5) investigated the origin and evolution of the ISM in an
extensive study of a complete sample of 56 elliptical galaxies using broad-
band imaging, narrow-band $H\alpha$ imaging, long slit spectroscopy, and
_Infrared Astronomical Satellite_ (IRAS) far-IR imaging. A study of the warm
and cold phases of the ISM in early-type galaxies was carried as part of the
ESO Key-Program 1-004-43K (Macchetto et al., 1996; Ferrari et al., 1999; Caon
et al., 2000; Ferrari et al., 2002, hereafter MFCF). Similar to G94-5, the
MFCF study probed the origin and evolution of the warm and cold ISM phases
through a combination of broad-band, narrow-band and far-IR imaging and long-
slit spectroscopy of 73 early-type galaxies. The two later studies had the
benefit of FIR imaging capabilities and added the cold ISM (i.e. dust) to the
discussion of the ISM in early-type galaxies. These studies sought to
determine the frequency of detection of warm and cold ISM in early-type
galaxies and then measure fundamental properties, such as mass, density, and
luminosity. By correlating these properties with stellar and hot ISM
properties, inferences are made about the origin and evolution of the ISM and
how it relates to the origin and evolution of the galaxy itself. Below we
briefly summarize the findings of these key studies and note the areas where
the studies agree and where questions remain.

The three key studies find similar fractions of their samples that host a warm
ISM, ranging from 60%-80%, despite the differences in sample selection
criteria and sizes of the samples. TdSA found 75% of the 13 early-type
galaxies surveyed had detectable $H\alpha$ emission which was extended over a
range of 5-10 $kpc$. The extent and elliptical shape of the emission confirmed
the galaxy-wide nature of this phenomena and excluded a nuclear confined
process, such as AGN activity. The TdSA survey was conducted on galaxies known
to have significant masses of hot ISM, and the authors concede that the high
detection rate is likely connected to sample selection. The G94-5 optically
complete sample (chosen without regard for X-ray luminosity) offers a less
biased measure of the frequency of the warm ISM; 60% of the early-type
galaxies surveyed were detected to have a warm ISM, while a cold (T$\sim 10K$)
dust component is detected in 40% of these galaxies. Through the FIR imaging,
the G94-5 study was able to show that unlike spiral galaxies, early-type
galaxies display a smooth dust component distributed throughout the galaxy.
The MFCF authors report an ionized gas detection frequency at 78% for the 73
galaxies in the optically complete sample. The cold dust component was
detected at a 75% detection rate with smooth, extended morphologies similar to
that of the warm ISM. Two key observations link the warm and the cold ISM
together; first, the similar spatial extent and morphology of the cold ISM
(dust) and warm ISM and second, the fact that the presence of dust is a strong
predictor of the presence of emission lines. It is concluded that these ISM
phases are closely related and likely of the same origin.

The mass of the warm ISM is orders of magnitude below the X-ray emitting hot
ISM and the stellar mass. TdSA calculated a warm ISM gas mass by assuming
case-B recombination, a filling factor of $10^{-5}-10^{-6}$, and densities
corresponding to the pressure implied by the X-ray observations. The resultant
masses were calculated to be $\sim 10^{5}-10^{6}M_{\sun}$ for a typical galaxy
in their sample. The $H\alpha$ imaging in G94-5 revealed warm ISM masses
similar to the FIR determined dust masses of $\sim 10^{4}-10^{5}M_{\sun}$. The
MFCF study determined masses for the warm ISM in the range of
$10^{3}-10^{5}M_{\sun}$. Masses of dust at $10^{5}M_{\sun}$ have direct
implications on the origin of the dust as the dust destruction rate is high in
the presence of the hot ISM and is further discussed below.

The correlations between the luminosities of $L_{FIR}$, $L_{H\alpha}$,
$L_{X}$, $L_{B}$ and $L_{X}/L_{B}$ explore the links between the difference
phases of the ISM and the stars. For numerous of reasons, the luminosity
correlations have been difficult to accurately measure and there are
contradictions between the key studies. TdSA observed that galaxies with more
massive X-ray halos had stronger line emission, indicating that the warm gas
is somehow linked with the hot gas. No other luminosity correlations were
observed. In a related work to G94-5, the authors find a strong anti-
correlation between the masses of dust and X-ray gas (Goudfrooij, 1994, 1995).
However in a later study with a larger sample (based in part on MFCF data),
Goudfrooij (1997) reveals a positive correlation between the warm ISM with
both the hot ISM and the stellar blue luminosity. MFCF determined the
luminosity of the ionized gas to be correlated with both blue and X-ray
luminosity. The scatter in the $L_{H\alpha}$ and $L_{B}$ relation is greatly
reduced when the blue luminosity is confined to the same region as the
emission lines. The scatter in all the $L_{H\alpha}$ to $L_{X}/L_{B}$
relations is large and plagued by upper limits on non-detections.

In general, the three studies agree that the ionization source for the warm
ISM is consistent with PAGB stars, but all comment that the data are unable to
exclude other excitation sources such as internal shocks, AGN heating or
electron conduction from the hot gas. TdSA show that the observed line
intensities are a factor of $\sim 10^{2}-10^{3}$ too large to result from
cooling flow luminosities (e.g. Fabian (1991)) and conclude that UV flux from
post asymptotic giant branch (PAGB) stars are the most likely ionizing
mechanism. The MFCF observed correlation between $H\alpha$ and blue luminosity
is an indication of photoionization of the warm ISM from PAGB stars. However,
the observed correlation with X-ray luminosity is consistent with warm ISM
excitation via electron conduction from the hot ISM (Sparks et al., 1989);
Neither ionization mechanism is preferred in MFCF’s large sample of 73
galaxies. See Section 4.1.1 for a continued discussion on the warm ISM
ionization source.

The key studies narrow down the origin of this dust and warm gas to either be
mass shed from AGB stars or the result of small galaxy accretion and mergers.
G94-5 reasons that cold dust originated in the circumstellar envelopes of AGB
stars are quickly sputtered away by the hot ISM with a typical dust grain
lifetime ($\sim 10^{7}yrs$) much shorter than the age of the galaxy (Draine &
Salpeter, 1979). And further, the balance between the dust injection and
destruction rate does not lead to the dust masses inferred by the far-IR
imaging (Faber & Gallagher, 1976; Knapp et al., 1992). The surprising result
from the confluence of the G94-5 studies is that the majority of the dust in
elliptical galaxies is external in origin. Since the distribution of the
ionized gas (i.e. the warm ISM) is observed to be similar to the dust
morphology, it is presumed that the warm ISM is likewise external to the
galaxy, and will only be seen in galaxies that either have recent merger
activity or significant inflows. The high detection frequency of the warm ISM,
implies mergers and smaller galaxy acquisitions are common in early-type
galaxies. Adding to the puzzle, in a subsample of the gas-rich galaxies, the
warm ISM was frequently observed to be kinematically distinct from the main
stellar population, also strongly suggesting an external origin for this gas
(Caon et al., 2000).

###  1.3 Warm ISM Open Issues

The warm ISM is important to study in early-type galaxies and despite the
ambitious efforts of the studies described above and subsequent works, several
fundamental questions remain. It is critical to resolve these open issues
since the intertwined relationships between the cold, warm, and hot gas phases
and the stellar population ultimately defines galaxy evolution. Even though
the warm ISM mass is orders of magnitude less than the stellar mass,
discriminating between internal versus external origin of this material
provides pivotal data in formation models. Also, even though minute in mass
compared to the cold and hot phases of the ISM, the optical line emissions of
the warm ISM can locally dominate the energy output for a region. Finally, the
optical emission lines of the warm ISM provides information on the heating and
cooling in these systems which demands further exploration since the observed
hot ISM properties along with the non-detection of cooling flows creates a
“cooling crisis.”

The metallicity of the warm ISM is unknown. Because the metallicities of the
hot ISM and the stellar metallicities (cf. Trager et al. (2000b, a)) are
known, determining the metallicity of the warm ISM can create a link between
the gas phases and the stellar population. Also, the metallicity of the warm
ISM contains important discriminating information concerning the outstanding
issues that the key warm ISM studies were not able to resolve. Specifically,
the line ratios of the ionization species provides information on the
excitation mechanism of the gas; currently, both electron conduction and
photoionization remain viable excitation mechanism.

Also the metallicity of the warm ISM contains clues to the internal versus
external origin debate. The gas injected into the ISM from stellar mass loss
will generally be more metal rich than gas accreted from dwarf galaxies (cf.
Lequeux et al. (1979) and Aaronson (1986)).

Because of the important information contained in the metallicity of the warm
ISM and the unresolved issues concerning its nature and origin, we have
conducted a program to determine the oxygen metallicity of the $T\sim 10^{4}K$
gas in a small sample of early-type galaxies. Section 2 describes the sample
and observations. The data processing techniques and spectral analysis is
discussed in Section 3. Results from the emission line analysis is presented
in Section 4. In Section 5, we conclude with a discussion of the implications
of this work on both the origin of the warm ISM and the relationships between
the stellar population and the other phases of the ISM.

##  2 Observations

The MDM 2.4 meter Hiltner Telescope coupled with a Boller and Chivens
spectrograph was used to observe 21 nearby, early-type galaxies in four
observation runs over a two year period from 2000 to 2002. The sample was
selected to span a range of hot ISM properties as determined from the Brown &
Bregman (1998) X-ray study as well as a range of warm ISM properties
determined in the narrow-band $H\alpha$ imaging surveys by TdSA, G94-5 and
MFCF. The galaxies with their relevant fundamental properties are listed in
Table 1. Our warm ISM isolation and detection method depends on observing
galaxies that contain no emission line gas, therefore, several galaxies were
chosen to have very weak detections of the hot and warm phases.

Table 1: Early-type Galaxy Properties: Warm Phase Sample  Galaxy | RAaaValues taken from NED (NASA/IPAC Extragalactic Database). | DecaaValues taken from NED (NASA/IPAC Extragalactic Database). |  $B_{T}^{0}$bbRadial velocity distance corrected for local flows from Faber et al. (1989) unless otherwise noted. | DbbRadial velocity distance corrected for local flows from Faber et al. (1989) unless otherwise noted. |  $log\,L_{B}/L_{\odot}$bbRadial velocity distance corrected for local flows from Faber et al. (1989) unless otherwise noted. |  $log\,L_{x}$ccX-ray luminosity from Brown & Bregman (1998) unless otherwise noted. |  $log\,L_{H_{\alpha}}$ddH$\alpha$ narrow-band imaging luminosity. |  ${H_{\alpha}}$ Reference eeReference for H$\alpha$ narrow-band imaging luminosity with the following abbreviations: M96 $=$ Macchetto et al. (1996), G94 $=$ Goudfrooij et al. (1994b), and T91 $=$ Trinchieri & di Serego Alighieri (1991).
---|---|---|---|---|---|---|---|---
| (J2000.0) | (J2000.0) |  | $(km/s)$ | $log~{}(erg~{}s^{-1})$ |
NGC 1407 | 03 40 11.8 | -18 34 48 | 10.57 | $1990\pm 187$ | 11.16 | 41.34 | 39.32 | M96
NGC 1600 | 04 31 39.8 | -05 05 10 | 11.79 | $4019\pm 489$ | 10.67 | 40.8411O’Sullivan et al. (2001). | 40.15 | M96
NGC 2768 | 09 11 37.5 | +60 02 15 | 10.93 | $1532\pm 325$ | 10.79 | 40.41 |  $\cdots$ |  $\cdots$
NGC 3115 | 10 05 13.9 | -07 43 07 | 9.95 | $1021\pm 215$ | 10.83 | 39.74 | $<37.6$ | M96
NGC 3377 | 10 47 42.3 | +13 59 08 | 11.13 | $857\pm 126$ | 10.21 | 39.42 | 38.95 | G94
NGC 3379 | 10 47 49.6 | +12 34 55 | 10.43 | $857\pm 126$ | 10.49 | 39.78 | 39.04 | M96
NGC 3489 | 11 00 18.3 | +13 54 05 | 11.1222RC3, de Vaucouleurs et al. (1991). |  $1039\pm 101$22RC3, de Vaucouleurs et al. (1991). | 10.3522RC3, de Vaucouleurs et al. (1991). |  $\cdots$ | 39.34 | M96
NGC 3607 | 11 16 54.3 | +18 03 10 | 10.53 | $1991\pm 242$ | 11.18 | 40.82 | 39.92 | M96
NGC 4125 | 12 08 05.8 | +65 10 27 | 10.58 | $1986\pm 295$ | 11.16 | 41.01 | 40.30 | T91
NGC 4261 | 12 42 02.4 | +11 38 48 | 11.32 | $2783\pm 590$ | 10.35 | 41.1822RC3, de Vaucouleurs et al. (1991). | 39.38 | G94
NGC 4374 | 12 25 03.7 | +12 53 13 | 10.13 | $1333\pm 71$ | 10.99 | 41.09 | 39.56 | G94
NGC 4406 | 12 26 11.7 | +12 56 46 | 9.87 | $1333\pm 71$ | 11.10 | 41.80 | 40.50 | T91
NGC 4472 | 12 29 46.8 | +08 00 02 | 9.32 | $1333\pm 71$ | 11.32 | 41.77 | 39.60 | T91
NGC 4494 | 12 31 24.1 | +25 46 28 | 10.69 | $695\pm 147$ | 10.20 | 39.28 |  $\cdots$ |  $\cdots$
NGC 4552 | 12 35 39.8 | +12 33 23 | 10.84 | $1333\pm 71$ | 10.71 | 40.92 | 39.26 | M96
NGC 4636 | 12 42 50.0 | +02 41 17 | 10.20 | $1333\pm 71$ | 10.96 | 41.81 | 39.69 | M96
NGC 4649 | 12 43 39.6 | +11 33 09 | 9.77 | $1333\pm 71$ | 10.96 | 41.48 | 39.83 | T91
NGC 4697 | 12 48 35.9 | -05 48 02 | 10.03 | $794\pm 168$ | 10.58 | 40.13 | 39.63 | G94
NGC 5044 | 13 15 23.9 | -16 23 08 | 11.25 | $2982\pm 314$ | 10.34 | 42.39 | 40.73 | M96
NGC 5322 | 13 49 15.2 | +60 11 26 | 11.09 | $1661\pm 352$ | 10.80 | 40.11 | 39.74 | G94
NGC 5846 | 15 06 29.2 | +01 36 21 | 10.67 | $2336\pm 284$ | 11.26 | 42.01 | 40.25 | M96

The configuration of the spectrograph was chosen to maximize galaxy light
input with a $2.1\arcsec$ wide by $5\arcmin$ long slit, while retaining the
ability to discriminate between the $H\alpha$ and [N II]$\lambda 6583$
emission lines. The available and appropriate grating for these requirements
was a 350 lines/mm grating blazed at $4026$Å, resulting in $7.1$Å/pixel and a
spectral range of $\sim 1600$Å over the 1200x800 pixel Loral CCD. The CCD was
characterized by relatively low read noise at seven electrons with the nominal
gain set at 2.1 electrons per ADU. The one drawback to the chosen spectrograph
configuration is that it was necessary to obtain separate blue and red spectra
for each galaxy in order to obtain all of the important emission lines from [O
II]$\lambda 3727$ to [Si II]$\lambda\lambda 6717,6731$. The grating is servo
controlled from the control computer and our tests indicated a $5$Å accuracy
in repositioning. Therefore, internal HgNe and Ne lamps were observed before
and after each grating reposition for wavelength calibration. Additional
calibrations included bias frames, evening and morning twilight flats when the
skies were clear, internal flats illuminated from an incandescent source, and
spectrophotometric standards.

The four observation runs were in March 2000, May 2000, February 2001, and
February 2002. Galaxy exposure times and program-type (E $=$ emission-line, T
$=$ Template, or L $=$ LINER) are listed in Table 2. The meaning of the
program types is discussed below. During the first three runs, the nights were
cloudy and over half of the full run was lost to weather and instrument
problems. Some of these observations are of marginal value and only reported
here for completeness. The final run was entirely photometric.

Table 2: MDM 2.4 Meter Observations Galaxy | Exp Blue (s) | Exp Red (s) | Emission Type aaEmission Type. Column display how the observations were ultimately used: T= template galaxy, L=LINER galaxy, E=emission galaxy, W=weak emission.
---|---|---|---
NGC 1407 | $8400$ | $6600$ | T
NGC 1600 | $8400$ | $5400$ | W
NGC 2768 | $12000$ | $11400$ | L
NGC 3115 | $6000$ | $4500$ | T
NGC 3377 | $7200$ | $6000$ | W
NGC 3379 | $6000$ | $6000$ | T
NGC 3489 | $6000$ | $4500$ | E
NGC 3607 | $9600$ | $7200$ | E
NGC 4125 | $12000$ | $10800$ | L
NGC 4261 | $6000$ | $4500$ | E
NGC 4374 | $6000$ | $4500$ | E
NGC 4406 | $9600$ | $6000$ | W
NGC 4472 | $9600$ | $6000$ | W
NGC 4494 | $1200$ | $0$ | W
NGC 4552 | $6000$ | $6600$ | W
NGC 4636 | $6000$ | $9900$ | E
NGC 4649 | $6000$ | $4500$ | W
NGC 4697 | $6000$ | $1800$ | W
NGC 5044 | $6000$ | $3600$ | E
NGC 5322 | $9600$ | $9900$ | W
NGC 5846 | $7200$ | $6600$ | E

The observing strategy was chosen to suppress spurious emission line
detections and obtain uniform observations of program and template galaxies.
For each galaxy, we attempted to obtain five integrations, ensuring that
reasonable statistics could be employed to eliminate cosmic rays from the
combined data. Because the emission lines are faint and spread over only a few
pixels in both spectral and spatial directions, a fortuitous single cosmic ray
can result in a spurious warm ISM detection. The total integration time was
chosen to result in a S/N of 5-10 in the template subtracted spectra (See
Section 3). We did not constrain the angle of the slit due to complications
with the instrument rotator during the first two runs; this has little
programatic impact because of the choice of a wide slit, and the non-
preferential orientation of the warm ISM observed in the narrow-band imaging
surveys (TdSA, G94-5, MFCF).

##  3 Reductions and Spectral Analysis

The data were reduced in the standard manner using tasks within IRAF (overscan
fitting and subtraction, bias frame construction and subtraction). Each of the
runs was calibrated separately but in a similar manner. Internal lamps were
used to create a response image and divided through the data, eliminating
differences in pixel-to-pixel sensitivity on the CCD. Twilight flats were used
to construct an illumination image, correcting for large scale structure and
slit illumination. For the first two runs, the incandescent internal lamp
produced too few photons in the far blue and only added noise to the data. For
these runs we used wavelength calibrated twilight flats to produce response
images. A HgNe lamp was used for wavelength calibration of blue-tuned
observations and a Ne lamp was used for red-tuned observations. Residuals in
the linear, spectral solution were much better (typically $<1$Å) than the
resolution element ($\sim 8$Å) over the entire chip.

The galaxies were summed over the spatial dimension based on the visual
inspection of the 2-D spectra, with background regions selected from the outer
regions of the chip. Typically the summed regions for the galaxies were about
$1.2\arcmin$ of the $5.2\arcmin$ unvignetted slit length. The spectra were
flux calibrated from observations of spectrophotometric standards (Massey et
al., 1988) and corrected for extinction with the standard correction for Kitt
Peak distributed with IRAF 2.11.

Figure 1: Composite Blue and Red Spectra of NGC3489. Important emission lines
that are visible in the un-subtracted data are marked aobve the continuum
level. In addition, calcium H & K lines are marked for reference.

Custom code was developed to implement the combining of the individual 1-D
spectra, template matching and subtraction, as well as the extraction of line
fluxes. The 1-D spectra combination procedure employed an iterative
$1.5\sigma$ rejection algorithm which eliminated cosmic rays that were
overlooked in the aperture summing rejection methods. The wavelength
calibration was tuned to night sky lines and the measured redshift was checked
and tuned to strong stellar population lines (e.g. Ca H&K in the blue and [N
II]$\lambda 6583$ in galaxies showing emission lines in the red). An example
summed spectra is shown in Figure 1.

Although we have upper-limits from the $H\alpha$ imaging surveys, it is
unknown which galaxies are void of warm ISM emission. Therefore, the
determination of template galaxies, or galaxies without any detected warm ISM,
is an involved and iterative process. First the $H\alpha$ and [N II]$\lambda
6583$ region of the spectra is examined for obvious emission. Using the local
continuum, upper limits on the line fluxes are determined for galaxies with no
obvious $H\alpha$ and [N II]$\lambda 6583$ emission. Ten of the original 21
galaxies match this criteria and act as potential template galaxies. We take
these ten galaxies and combine them into a master template galaxy. When
combined into the template galaxy, the strongest of these weak emitters will
not be included in the template creation because of the aggressive clipping
algorithm in the combination routine. To determine the the weak emitters the
template is subtracted from the individual galaxy observations and examined
for emission lines. Three of the ten template galaxies revealed oxygen lines
and were removed from the master template. This process was repeated with the
emitters removed and reveal four more weak emitters, leaving three galaxies
without emission line activity: NGC 1407, NGC 3115, NGC 3379. As an
independent check of our template definition methodology, all seven of the
potential template galaxies were subtracted from the galaxy with the strongest
emission lines (NGC 2768) and we determined that the three template galaxies
selected resulted in the high reported fluxes (i.e. lowest background).

The three template galaxies, NGC 1407, NGC 3115 and NGC 3379, lack of warm
emission is consistent with results from the narrow-band imaging detections
for these galaxies. NGC 3379 is reported to have a weak $H\alpha$ in MFCF but
an upper limit below this reported detection is determined in G94-5. NGC 3115
is detected only as an upper limit by MFCF and has been previously used as a
stellar template galaxy by the Palomar spectral survey of the nuclear regions
of 500 nearby galaxies (Ho et al., 1993). NGC 1407 is marginally detected in
both the G94-5 and MFCF studies, but the fluxes reported differ by a factor of
five. The inconsistencies of weak detections in the narrow-band surveys reveal
intrinsic limitations of the narrow-band technique. Broad-band red imaging is
used to subtract an underlying stellar population from the narrow-band data; a
process which is susceptible to scaling errors. The stellar template
constructed from our sample contains $<10^{-16}~{}erg\,s^{-1}\,cm^{-2}$
combined $H\alpha$ and [N II]$\lambda 6583$ flux; this is an order of
magnitude more sensitive as a non-detection than the narrow-band imaging
surveys.

Once the template galaxies are defined, the stellar population was subtracted
from each of the other galaxies, revealing emission lines and any differences
in stellar populations. An example subtraction for NGC 4374 is shown in Figure
2. In all galaxies, the residuals in the subtractions are flat with large-
scale variation at a 8-10% level in the blue and better than 5% in the red.
Once the subtractions are made, the relative line fluxes are measured by
defining a low order polynomial to describe the local background and using a
gaussian profile to extract the line flux. A gaussian decomposition is only
necessary in the crowded [N II]$\lambda 6548$, $H\alpha$, [N II]$\lambda 6583$
region and the fitting method gives equivalent results to raw counts above a
background for all other lines. Note that because of the coarse spectral
sampling and the relatively low velocity dispersions, the width of the
gaussian profiles contain no useful kinematic information and is thus not
reported. In eight of the non-template galaxies, the signal-to-noise in the
data are insufficient for emission line work, or only H$\alpha$ is observed.
For the rest of the sample, extracted relative line fluxes normalized to
$H\beta$ are reported in Table 3. The reported flux ratios have been
dereddened based on an average value for the reddening from Burstein & Heiles
(1984) and Schlegel et al. (1998). The observed relative fluxes reported here
do not correct for internal galaxy reddening, displaying the Balmer decrements
as observed; Below we apply the correction (See Section 4.2). The errors
reported are determined by calculating the range of acceptable background
levels within the noise of the local continuum for isolated lines and
combining this background uncertainty in quadrature with the range of
acceptable gaussian widths that reproduce the total observed line flux for the
crowded [N II]-$H\alpha$ complex.

Table 3: Observed Nebular Emission Line Ratios Relative to H$\beta$ Galaxy | [O II] | [Ne III] | $H\beta$ | [O III] | [O III] | [O I] | [N II] | $H\alpha$ | [N II] | [S II] |
---|---|---|---|---|---|---|---|---|---|---|---
| $\lambda 3727$ | $\lambda 3869$ | $\lambda 4862$ | $\lambda 4959$ | $\lambda 5007$ | $\lambda 6300$ | $\lambda 6548$ | $\lambda 6563$ | $\lambda 6583$ | $\lambda\lambda 6717,6731$ |
NGC 2768 | $70\pm 7$ | $14\pm 6$ | $10\pm 3$ | $1\pm 4$ | $24\pm 4$ | $8\pm 2$ | $20\pm 4$ | $29\pm 5$ | $45\pm 7$ | $50\pm 8$ |
NGC 3489 | $19\pm 5$ | $2\pm 3$ | $10\pm 2$ | $2\pm 5$ | $13\pm 5$ | $3\pm 2$ | $30\pm 8$ | $34\pm 5$ | $57\pm 8$ | $37\pm 7$ |
NGC 3607 | $32\pm 8$ | $6\pm 3$ | $10\pm 3$ | $7\pm 3$ | $20\pm 4$ | $2\pm 4$ | $25\pm 5$ | $48\pm 6$ | $59\pm 8$ | $22\pm 5$ |
NGC 4125 | $74\pm 8$ | $14\pm 4$ | $10\pm 2$ | $23\pm 5$ | $71\pm 10$ | $3\pm 3$ | $19\pm 4$ | $34\pm 8$ | $40\pm 5$ | $23\pm 7$ |
NGC 4261 | $16\pm 4$ | $6\pm 5$ | $10\pm 3$ | $6\pm 2$ | $12\pm 3$ | $3\pm 2$ | $15\pm 5$ | $22\pm 8$ | $39\pm 5$ | $15\pm 6$ |
NGC 4374 | $16\pm 4$ | $<1$ | $10\pm 2$ | $5\pm 2$ | $12\pm 4$ | $2\pm 4$ | $21\pm 9$ | $35\pm 6$ | $49\pm 8$ | $29\pm 5$ |
NGC 4636 | $24\pm 5$ | $2\pm 4$ | $10\pm 3$ | $4\pm 4$ | $13\pm 3$ | $2\pm 5$ | $10\pm 3$ | $33\pm 5$ | $49\pm 7$ | $21\pm 5$ |
NGC 5044 | $24\pm 3$ | $2\pm 3$ | $10\pm 2$ | $6\pm 2$ | $20\pm 3$ | $3\pm 2$ | $33\pm 8$ | $32\pm 6$ | $53\pm 8$ | $33\pm 5$ |
NGC 5846 | $26\pm 5$ | $2\pm 3$ | $10\pm 4$ | $6\pm 3$ | $14\pm 4$ | $<1$ | $8\pm 5$ | $29\pm 4$ | $36\pm 8$ | $13\pm 10$ |

Figure 2: Blue (left) and Red (right) spectra of NGC4374 with template
galaxies subtracted, removing the stellar continuum. Large residuals at 5500Å
and 6800Å are due different night sky lines positions in the rest wavelength
of NGC4374 and template galaxies. The large residual at 3900Å is due to the
imperfect matching between the Ca H&K lines of the template galaxies and NGC
4374.

##  4 Warm ISM Emission Line Results

###  4.1 Excitation Mechanism

####  4.1.1 Normal Galaxies

There are many possible ionization mechanisms for the warm ISM that have been
previously explored and ruled out by the key warm ISM studies (Section 1.2)
and other works. Photoionization by hot young stars (Kim, 1989; Shields, 1991)
is inconsistent with both the optical colors (G94-5; MFCF) and integrated
optical spectra (Sadler, 1987; Heckman et al., 1988). Ionization from
radiation associated with an AGN (Fosbury et al., 1982) is ruled out by the
large radial extent of the warm ISM observed in the narrow-band imaging (TdSA;
G94-5; MFCF). Shock excitation (Ford & Butcher, 1979; Heckman et al., 1989)
would produce filamentary ionization regions, where smooth emission is
observed (TdSA; G94-5; MFCF). Further, the the input energy from shocks is too
low by two orders of magnitude for typical densities and velocities as noted
by Sparks et al. (1989). Condensation out of the hot phase and into warm phase
filaments has been investigated for a number of X-ray clusters (Cowie et al.,
1980; Hu et al., 1985) but the large number of recombinations required per
hydrogen atom (e.g. Johnstone et al. (1987); TdSA and references therein) is
inconsistent with the observed recombination rates for early-type galaxies
(Sparks et al., 1989).

The three remaining viable excitation mechanisms for the warm ISM are
photoionization from PAGB stars, electron conduction from the hot phase (e.g.
Sparks et al. (1989)), and photoionization by extreme UV (EUV) photons from
the hot phase (Voit & Donahue, 1990; Donahue & Voit, 1991). The data in the
key warm ISM studies are equally consistent with all of these ionization
mechanisms. However, photoionization by PAGB stars is the preferred mechanism
since these old stars are known to populate the galaxies and given the
significant correlation of $L_{H\alpha}$ with $L_{B}$ within a given region
(MFCF) strongly indicates a stellar ionizing source. Further, the presence of
the far ultraviolet (UV) flux has been confirmed in a small sample of 30
early-type galaxies with _International Ultraviolet Explorer_ (IUE) (Burstein
et al. (1988); For a review of Far UV flux in early-type galaxies see
O’Connell (1999)).

However, to adopt PAGB stars as the ionizing source, excludes the observed
correlations between $H\alpha$ and X-ray luminosity (TdSA; (Goudfrooij, 1997),
MFCF). The correlations between the X-ray and warm ISM morphologies (Singh et
al., 1995; Trinchieri et al., 1997; Trinchieri & Goudfrooij, 2002) provide
additional support for photoionization of the warm ISM by the EUV flux from
the hot ISM, but it is possible to construct physical conditions that
reproduce these morphological links (e.g. cooling of the hot ISM). Although
there have been detected correlations between the $L_{H\alpha}$ and
$L_{X}/L_{B}$, the correlations contain much more scatter and are less
significant than the correlation between $L_{H\alpha}$ with $L_{B}$ within a
given region. The large amount of scatter in the X-ray relationship and the
known high $L_{X}/L_{B}$ galaxies which contain weak or no $H\alpha$ emission
and strong $H\alpha$ emitters with low $L_{X}/L_{B}$ provide insight that the
hot and warm ISM link is more tenuous than the link between the stars and the
warm ISM. Finally, the key studies in in Section 1.2 and subsequent work
(Sarzi et al., 2006), indicate a strong connection between the dust (cold ISM)
and the emission line regions, but there are not strong correlations between
$L_{FIR}$ and $L_{X}/L_{B}$ in any of the key studies or in recent Spitzer
Observations (Temi et al., 2007).

Emission line ratios of the different ionization species have the ability to
discriminate between collisionally excitation and photoionization. Comparing
our sample’s line ratios to line ratios predicted for photoionization models
of Allen et al. (1998), Allen et al. (1999), and Binette et al. (1996) and
shock ionization of Allen et al. (2008), we conclude the data are consistent
with photoionization. Unfortunately, these data are not discriminating due to
the large errors in the data and the region populated by these galaxies
happens to land in an overlap region for collisional and photoionization.

Additional support for the photoionization of this gas is seen with the Far
Ultraviolet Space Explorer (FUSE) (Bregman et al., 2005). The observed UV
radiation field produces sufficient photons shortward of 912Å to ionize the
warm ISM and is consistent with the number of predicted extreme horizontal-
branch stars, or equivalently PAGB stars.

When the data presented here is combined with the literature, it presents a
reasonable case for photoionization of the warm ISM. This is the preferred
ionization mechanism in the literature and in discriminating between the two
sources of photoionization, PAGB stars or the hot ISM, it would seem more
probable that PAGB stars provide the necessary ionizing photons.

####  4.1.2 LINER Galaxies

Four of the galaxies listed in Table 3 are reported as low ionization nuclear
emission-line regions (LINERs) galaxies in the literature: NGC 2768 (Heckman,
1980), NGC 3115 (Ho et al., 1993), and NGC 4125 and NGC 4261 (Ho et al.,
1995). The question of identical treatment of these LINER galaxies and the
normal galaxies arises. LINERs are strictly defined to be galaxies with line
ratios of [O II]$\lambda 3727$ greater than [O III]$\lambda 5007$ and [O
I]$\lambda 6100$ greater than one-third [O III]$\lambda 5007$ (Heckman, 1980).
While useful for separating classes of objects, the LINER definition as noted
by its creator, Heckman, has no physical basis and in reality only provides a
crude boundary between photoionization and a low-level AGN activity.

For our purposes, the LINER identification acts as a flag for further
investigation, and not as an automatic exclusion filter. The line ratios for
the four LINER galaxies are consistent within the photoionization, but the
excitation mechanism and thus the valid regions in diagnostic diagrams are
still uncertain and being debated for LINERs (for a review of the physical
mechanisms behind LINERs see Filippenko (1996) and Barth (2002)). Although it
is not entirely clear what to expect for LINER galaxies, the presence of
strong collisionally excited lines casts a cloud of doubt of pure
photoionization. The [Ne III]$\lambda 3869$ line is an indicator of the
importance of collisional excitation because the collisional cross-section is
over an order of magnitude larger than the photoionization transition
probability for this line. NGC 2768 and NGC 4125 have significant [Ne
III]$\lambda 3869$ fluxes and thus we conclude that these systems have a more
complex physical structure and excitation mechanism than photoionization.
These two galaxies are excluded in the subsequent sections where pure
photoionization is assumed. Note that NGC 4125 and NGC 4261 show no
distinguishing line ratios from the normal galaxies and, the data presented
here do not support the LINER classification. These two galaxies were
classified as LINERs by (Ho et al., 1995) which had a similar observational
aperture as this study, but only had wavelength coverage down to 4200Å,
excluding the useful [O II]$\lambda 3727$ information.

###  4.2 Oxygen Metallicity

For those galaxies where photoionization is the dominant ionization mechanism,
the line ratios can be used to determine the metallicity of the warm phase
gas. There are two critical assumptions to be addressed before using these
line ratios and photoionization models to determine global metallicities. The
first assumption is that UV flux from PAGB stars ionizing a diffuse ISM on a
galactic scale is well represented by photoionization models that were
developed for Str$\ddot{{o}}$mgren-sphere-type H II regions. First, we note
that the models are determining the balance between photoionization by UV
photons and recombinations to the ions; this is fundamentally the same
equation in a galaxy-wide ISM as an HII region, given a source of ionizing
photons internal to a volume of gas. The UV flux and electron densities of the
two environments are similar, although, PAGB stars have different spectral
shapes than high mass main sequence stars, but McGaugh (1991) shows that the
photoionization models are largely insensitive to the spectral shape of the
ionizing star. The models assume the gas is ejected from the ionizing star
(and thus of the same metallicity), this is a feature used to calibrate the
observed line ratios, but a direct link between the ionizing source and the
ejecta is not a necessary condition to employ the models. We argue that the
photoionization models will be applicable to the warm ISM, but we readily
acknowledge that there are differences in these two environments which will
lead to will lead to greater uncertainties in its application here.

The second key assumption is that the analysis of the linear combination of
many individual diffuse ionized regions observed in an integrated galaxy
spectra produces a result that is indicative of the mean properties of the
observed regions. Kobulnicky et al. (1999) conducted a study addressing these
issues in nearby galaxies, observing H II regions individually through small
apertures and globally with drift scanning with a long slit. The authors
conclude that global metallicity determinations, accurate to $0.2$ dex, are
possible with long-slit, integrated spectra.

Table 4: R23 parameter and Oxygen Metallicity Galaxy | log R23 | [O/H] | $Z_{Oxygen}/Z_{\sun}$
---|---|---|---
NGC 3489 | $0.63^{+0.11}_{-0.14}$ | $-3.26\pm 0.2$ | $1.13^{+0.49}_{-0.30}$
NGC 3607 | $0.93^{+0.08}_{-0.11}$ | $-3.66\pm 0.2$ | $0.43^{+0.29}_{-0.18}$
NGC 4261 | $0.51^{+0.06}_{-0.07}$ | $-3.14\pm 0.2$ | $1.44^{+0.47}_{-0.29}$
NGC 4374 | $0.58^{+0.09}_{-0.10}$ | $-3.20\pm 0.2$ | $1.26^{+0.50}_{-0.31}$
NGC 4636 | $0.68^{+0.09}_{-0.11}$ | $-3.31\pm 0.2$ | $0.98^{+0.41}_{-0.26}$
NGC 5044 | $0.75^{+0.05}_{-0.06}$ | $-3.37\pm 0.2$ | $0.86^{+0.36}_{-0.24}$
NGC 5846 | $0.68^{+0.08}_{-0.09}$ | $-3.30\pm 0.2$ | $1.00^{+0.37}_{-0.26}$

Metallicity determinations of ionized gas in individual H II regions are built
upon photoionization models, determining ratios of ionized species relative to
hydrogen. The most common diagnostic is the $R23$ parameterization,
$R23\equiv(([O\,II]\,\lambda 3727+[O\,III]\,\lambda\lambda
4959,\,5007)/H\beta)$, defined by Pagel et al. (1979) and refined by McGaugh
(1991). This parameterization is convenient since it employs the readily
available optical oxygen lines, however, it requires the [O II]$\lambda
3727$/[N II]$\lambda 6583$ information to resolve a degeneracy present in the
$R23$ models.

We take the line ratios reported in Table 3 and correct for reddening internal
to the galaxy using the interstellar extinction curve of Savage & Mathis
(1979) and assuming case B ratio for the $H\alpha$ to $H\beta$ at $T=10^{4}$
(Osterbrock, 1989). For NGC 4261 no reddening correction is applied, as the
Balmer decrement indicates a non-physical negative reddening.

These dereddened lines are used to generate $R23$ values for the galaxies.
Because of the faint nature of the $[O\,III]\,\lambda 4959$ line, we used the
conical, equilibrium value of 1/3 the line flux of $[O\,III]\,\lambda 5007$
when this line was buried in the noise of the subtracted spectra. All of the
galaxies end up on the upper branch of the $R23$ curve as indicated by the [O
II]$\lambda 3727$/[NII]$\lambda 6583$ diagnostic. To determine metallicities,
we employ McGaugh (1991) calibration as parameterized by Kobulnicky et al.
(1999) and reported the metallicity results in Table 4. For the zero point of
the solar metallicity we employ a solar photospheric oxygen abundance of
[O/H]$=-3.3$ (Asplund et al., 2005; Scott et al., 2008), but note the
uncertainty and current debate about the solar determinations. The errors
reported on the oxygen metallicities in these galaxies are at least $0.2$ dex
as these are the dominating and intrinsic uncertainties in the models used to
create the $R23$ grid and its application to global galaxy spectra (cf.
Kobulnicky et al. (1999)). The median metallicity for the seven early-type
galaxies is O/O☉=1.01. Within these $0.2$ dex absolute errors, the none of the
galaxies are significantly different from one another.

##  5 Implications and Discussion

The determined mean oxygen metallicity of solar for the warm phase of the ISM
has implications for the origin of this gas and its relationship to other
phases of the ISM and the stellar population. The limitations of this study
are apparent with a small sample of seven galaxies, large $0.2$ dex errors,
and a less than ideal application of H II region models to the warm ISM.
However, this work does present new information as it appears that there are
no determination of the metallicity of the warm ISM in early-type galaxies in
the literature. Sarzi et al. (2006) conducted the most recent large survey of
ionized gas in early-type galaxies and conclude at the end of the study that
the origin of the ionized gas is still unsolved and state that the metallicity
provides an insightful clue. Below we consider the implications of a near
solar metallicity for the warm ISM.

###  5.1 Origin of the Warm ISM

A warm ISM with a solar metallicity argues against an external origin model
for the emission line gas. The accretion galaxies involved in minor mergers
are presumably dwarf irregular (dIrr) galaxies, since dwarf spheroidal (dSph)
galaxies are most frequently devoid of gas. Gaseous and dusty galaxies with
near solar metallicity are luminous ($>10^{10}L_{\sun}$) and massive(Garnett,
2002). Mergers of these types of objects produce profound disturbances, which
are simply not observed in these early-type galaxies with $H\alpha$ emission.
It is possible that an intense gas enrichment by supernovae type-II occurs
once low metallicity gas is accreted, however, the broad-band colors and
global galaxy spectra (G94-5, MFCF) show little evidence for major star
formation. Further, the observed supernovae rates are too low to account for
such a drastic enrichment (Cappellaro et al., 1999).

We can compare a solar metallicity for the warm ISM with metallicities of the
hot ISM and the stars. The stellar population in early-type galaxies is
generally observed to have a solar or super-solar metallicity (Trager et al.,
2000b, a). A number of the seven galaxies observed in this work have stellar
Fe and $\alpha$ metallicity measurements in the literature, although, there
are some inconsistencies between authors (NGC 4261 is measured to have
$[Fe/H]$ of -0.03, 0.275, 0.29, and 0.188 by Trager et al. (2000b), Thomas et
al. (2005), Howell (2005), and Sánchez-Blázquez et al. (2006), respectively.)
Not all the galaxies have measured stellar metallicities but a five galaxy
average from all of the determinations in the 2005 and 2006 papers listed
above is $[Fe/H]=0.22$, ranging from $[Fe/H]=0.099$ to $[Fe/H]=0.44$
($[\alpha/H]$ is similar). It is important to remember that Trager et al.
(2000a) estimates the _absolute_ error on stellar metallicity determinations
close to a factor of two, while the relative errors are much smaller. The
_Chandra_ X-ray Observatory has the spatial resolution to resolve and remove
the stellar component from the gas which is leading to better data and a
revision of our (mis)understanding of the metallicity in the hot ISM.
Recently, Humphrey & Buote (2006) analyzed a sample of 28 early-type galaxies
observed with _Chandra_ and determine that the hot ISM abundances are
consistent and perhaps higher than the stellar abundances. In summary, a
$0.50$ to $1.50$ solar metallicity in the warm ISM is consistent with both hot
ISM and stellar metallicity determinations.

It has been thought that early-type galaxies would have a problem generating
and retaining dust and warm gas, due to the hostile conditions provided by the
hot ISM and the low injection rate of new material into the ISM from stars.
Parriott & Bregman (2008) have shown that is it possible for dust to condense
in the circumstellar envelopes of AGB stars and then subsequently be injected
into the hot ISM. Although this dust has a short lifetime in the hot ISM, the
injection rate has been observed to be greater than the destruction rate and
could account for the total observed mass given a simple evolution of this
injection rate given that more massive stars evolved off the main sequence in
the past (Athey et al., 2002). An internal origin for the dust is supported by
the G94-5 and MFCF studies indicating that the majority ($\sim 90\%$) of the
dust is distributed smoothly throughout the galaxy. The internal origin for
the cold ISM is further confirmed by a 2.29$\mu$m CO absorption feature study
that finds a lower velocity dispersion for the dust than that of the stars in
a sample of 25 nearby early-type galaxies (Silge & Gebhardt, 2003). Finally
from the analysis of _Spitzer_ observations of 46 galaxies in the FIR
($24\mu$m, $70\mu$m, and $160\mu$m), Temi et al. (2007) propose an AGB mass
loss, internal origin for the observed dust, which sinks to the core of the
galaxy and is periodically driven out by AGN activity. This outward transport
of dust may cool the hot gas down to the warm phase.

The major stumbling block to a galaxy internal origin for dust and gas is the
observations of gas kinematics which are distinct with respect to the stars.
Recent observational evidence has suggested that there may be different
evolutionary tracks for early-type galaxies. The classification of early-type
galaxies into fast and slow rotation categories by Emsellem et al. (2007) and
the key ingredient of gas for the formation and evolution of fast rotators
provides an important dividing line for future discussions on the origin of
the ISM in early-type galaxies.

The authors would like to acknowledge useful discussions with P. Goudfrooij
and E. D. Miller. We would also like to thank an anonymous referee for
suggestions that significantly improved this paper.

## References

  * Aaronson (1986) Aaronson, M. 1986, in Star Forming Dwarf Galaxies and Related Objects, 125
  * Allen et al. (1998) Allen, M. G., Dopita, M. A., & Tsvetanov, Z. I. 1998, ApJ, 493, 571
  * Allen et al. (1999) Allen, M. G., Dopita, M. A., Tsvetanov, Z. I., & Sutherland, R. S. 1999, ApJ, 511, 686
  * Allen et al. (2008) Allen, M. G., Groves, B. A., Dopita, M. A., Sutherland, R. S., & Kewley, L. J. 2008, ApJS, 178, 20
  * Asplund et al. (2005) Asplund, M., Grevesse, N., & Sauval, A. J. 2005, in Astronomical Society of the Pacific Conference Series, Vol. 336, Cosmic Abundances as Records of Stellar Evolution and Nucleosynthesis, ed. T. G. Barnes, III & F. N. Bash, 25–+
  * Athey et al. (2002) Athey, A., Bregman, J. N., Bregman, J. D., Temi, P., & Sauvage, M. 2002, ApJ, 571, 272
  * Barth (2002) Barth, A. J. 2002, in ASP Conf. Ser. 258: Issues in Unification of Active Galactic Nuclei, 147
  * Binette et al. (1996) Binette, L., Wilson, A. S., & Storchi-Bergmann, T. 1996, A&A, 312, 365
  * Bregman (1978) Bregman, J. N. 1978, ApJ, 224, 768
  * Bregman et al. (2005) Bregman, J. N., Miller, E. D., Athey, A. E., & Irwin, J. A. 2005, ApJ, 635, 1031
  * Brown & Bregman (1998) Brown, B. A. & Bregman, J. N. 1998, ApJ, 495, L75
  * Burstein et al. (1988) Burstein, D., Bertola, F., Buson, L. M., Faber, S. M., & Lauer, T. R. 1988, ApJ, 328, 440
  * Burstein & Heiles (1984) Burstein, D. & Heiles, C. 1984, ApJS, 54, 33
  * Caldwell (1984) Caldwell, N. 1984, PASP, 96, 287
  * Caon et al. (2000) Caon, N., Macchetto, D., & Pastoriza, M. 2000, ApJS, 127, 39
  * Cappellaro et al. (1999) Cappellaro, E., Evans, R., & Turatto, M. 1999, A&A, 351, 459
  * Cowie et al. (1980) Cowie, L. L., Fabian, A. C., & Nulsen, P. E. J. 1980, MNRAS, 191, 399
  * de Vaucouleurs et al. (1991) de Vaucouleurs, G., de Vaucouleurs, A., Corwin, H. G., Buta, R. J., Paturel, G., & Fouque, P. 1991, Third Reference Catalogue of Bright Galaxies (Volume 1-3, XII, 2069 pp. 7 figs.. Springer-Verlag Berlin Heidelberg New York)
  * Demoulin-Ulrich et al. (1984) Demoulin-Ulrich, M.-H., Butcher, H. R., & Boksenberg, A. 1984, ApJ, 285, 527
  * Donahue & Voit (1991) Donahue, M. & Voit, G. M. 1991, ApJ, 381, 361
  * Draine & Salpeter (1979) Draine, B. T. & Salpeter, E. E. 1979, ApJ, 231, 438
  * Emsellem et al. (2007) Emsellem, E., Cappellari, M., Krajnović, D., van de Ven, G., Bacon, R., Bureau, M., Davies, R. L., de Zeeuw, P. T., Falcón-Barroso, J., Kuntschner, H., McDermid, R., Peletier, R. F., & Sarzi, M. 2007, MNRAS, 379, 401
  * Faber & Gallagher (1976) Faber, S. M. & Gallagher, J. S. 1976, ApJ, 204, 365
  * Faber et al. (1989) Faber, S. M., Wegner, G., Burstein, D., Davies, R. L., Dressler, A., Lynden-Bell, D., & Terlevich, R. J. 1989, ApJS, 69, 763
  * Fabian (1991) Fabian, A. C. 1991, in IAU Symp. 144: The Interstellar Disk-Halo Connection in Galaxies, 237–244
  * Ferrari et al. (1999) Ferrari, F., Pastoriza, M. G., Macchetto, F., & Caon, N. 1999, A&AS, 136, 269
  * Ferrari et al. (2002) Ferrari, F., Pastoriza, M. G., Macchetto, F. D., Bonatto, C., Panagia, N., & Sparks, W. B. 2002, A&A, 389, 355
  * Filippenko (1996) Filippenko, A. V. 1996, in ASP Conf. Ser. 103: The Physics of Liners in View of Recent Observations, 17
  * Ford & Butcher (1979) Ford, H. C. & Butcher, H. 1979, ApJS, 41, 147
  * Forman et al. (1979) Forman, W., Schwarz, J., Jones, C., Liller, W., & Fabian, A. C. 1979, ApJ, 234, L27
  * Fosbury et al. (1982) Fosbury, R. A. E., Boksenberg, A., Snijders, M. A. J., Danziger, I. J., Disney, M. J., Goss, W. M., Penston, M. V., Wamsteker, W., Wellington, K. J., & Wilson, A. S. 1982, MNRAS, 201, 991
  * Garnett (2002) Garnett, D. R. 2002, ApJ, 581, 1019
  * Goudfrooij (1994) Goudfrooij, P. 1994, PhD thesis, Ph. D. thesis, University of Amsterdam, The Netherlands, (1994)
  * Goudfrooij (1995) —. 1995, PASP, 107, 502
  * Goudfrooij (1997) Goudfrooij, P. 1997, in Astronomical Society of the Pacific Conference Series, Vol. 116, The Nature of Elliptical Galaxies; 2nd Stromlo Symposium, ed. M. Arnaboldi, G. S. Da Costa, & P. Saha, 338–+
  * Goudfrooij & de Jong (1995) Goudfrooij, P. & de Jong, T. 1995, A&A, 298, 784
  * Goudfrooij et al. (1994a) Goudfrooij, P., de Jong, T., Hansen, L., & Norgaard-Nielsen, H. U. 1994a, MNRAS, 271, 833
  * Goudfrooij et al. (1994b) Goudfrooij, P., Hansen, L., Jorgensen, H. E., & Norgaard-Nielsen, H. U. 1994b, A&AS, 105, 341
  * Goudfrooij et al. (1994c) Goudfrooij, P., Hansen, L., Jorgensen, H. E., Norgaard-Nielsen, H. U., de Jong, T., & van den Hoek, L. B. 1994c, A&AS, 104, 179
  * Habing (1996) Habing, H. J. 1996, A&A Rev., 7, 97
  * Heckman (1980) Heckman, T. M. 1980, A&A, 87, 152
  * Heckman et al. (1988) Heckman, T. M., Baum, S. A., van Breugel, W., & McCarthy, P. J. 1988, in NATO ASIC Proc. 229: Cooling Flows in Clusters and Galaxies, 245–250
  * Heckman et al. (1989) Heckman, T. M., Baum, S. A., van Breugel, W. J. M., & McCarthy, P. 1989, ApJ, 338, 48
  * Ho et al. (1995) Ho, L. C., Filippenko, A. V., & Sargent, W. L. 1995, ApJS, 98, 477
  * Ho et al. (1993) Ho, L. C., Filippenko, A. V., & Sargent, W. L. W. 1993, ApJ, 417, 63
  * Howell (2005) Howell, J. H. 2005, AJ, 130, 2065
  * Hu et al. (1985) Hu, E. M., Cowie, L. L., & Wang, Z. 1985, ApJS, 59, 447
  * Humphrey & Buote (2006) Humphrey, P. J. & Buote, D. A. 2006, ApJ, 639, 136
  * Johnstone et al. (1987) Johnstone, R. M., Fabian, A. C., & Nulsen, P. E. J. 1987, MNRAS, 224, 75
  * Kim (1989) Kim, D. 1989, ApJ, 346, 653
  * Knapp et al. (1992) Knapp, G. R., Gunn, J. E., & Wynn-Williams, C. G. 1992, ApJ, 399, 76
  * Kobulnicky et al. (1999) Kobulnicky, H. A., Kennicutt, R. C., & Pizagno, J. L. 1999, ApJ, 514, 544
  * Lequeux et al. (1979) Lequeux, J., Rayo, J. F., Serrano, A., Peimbert, M., & Torres-Peimbert, S. 1979, A&A, 80, 155
  * Macchetto et al. (1996) Macchetto, F., Pastoriza, M., Caon, N., Sparks, W. B., Giavalisco, M., Bender, R., & Capaccioli, M. 1996, A&AS, 120, 463
  * Massey et al. (1988) Massey, P., Strobel, K., Barnes, J. V., & Anderson, E. 1988, ApJ, 328, 315
  * Mathews & Baker (1971) Mathews, W. G. & Baker, J. C. 1971, ApJ, 170, 241
  * Mayall (1936) Mayall, N. U. 1936, PASP, 48, 14
  * Mayall (1958) Mayall, N. U. 1958, in IAU Symp. 5: Comparison of the Large-Scale Structure of the Galactic System with that of Other Stellar Systems, 23
  * McGaugh (1991) McGaugh, S. S. 1991, ApJ, 380, 140
  * O’Connell (1999) O’Connell, R. W. 1999, ARA&A, 37, 603
  * Osterbrock (1989) Osterbrock, D. E. 1989, Astrophysics of gaseous nebulae and active galactic nuclei (Research supported by the University of California, John Simon Guggenheim Memorial Foundation, University of Minnesota, et al. Mill Valley, CA, University Science Books, 1989, 422 p.)
  * O’Sullivan et al. (2001) O’Sullivan, E., Forbes, D. A., & Ponman, T. J. 2001, MNRAS, 328, 461
  * Pagel et al. (1979) Pagel, B. E. J., Edmunds, M. G., Blackwell, D. E., Chun, M. S., & Smith, G. 1979, MNRAS, 189, 95
  * Parriott & Bregman (2008) Parriott, J. R. & Bregman, J. N. 2008, ApJ, 681, 1215
  * Phillips et al. (1986) Phillips, M. M., Jenkins, C. R., Dopita, M. A., Sadler, E. M., & Binette, L. 1986, AJ, 91, 1062
  * Sadler (1987) Sadler, E. M. 1987, in IAU Symp. 127: Structure and Dynamics of Elliptical Galaxies, 125–132
  * Sánchez-Blázquez et al. (2006) Sánchez-Blázquez, P., Gorgas, J., Cardiel, N., & González, J. J. 2006, A&A, 457, 809
  * Sandage (1957) Sandage, A. 1957, ApJ, 125, 422
  * Sarzi et al. (2006) Sarzi, M., Falcón-Barroso, J., Davies, R. L., Bacon, R., Bureau, M., Cappellari, M., de Zeeuw, P. T., Emsellem, E., Fathi, K., Krajnović, D., Kuntschner, H., McDermid, R. M., & Peletier, R. F. 2006, MNRAS, 366, 1151
  * Savage & Mathis (1979) Savage, B. D. & Mathis, J. S. 1979, ARA&A, 17, 73
  * Schlegel et al. (1998) Schlegel, D. J., Finkbeiner, D. P., & Davis, M. 1998, ApJ, 500, 525
  * Scott et al. (2008) Scott, P., Asplund, M., Grevesse, N., & Sauval, A. J. 2008, ArXiv e-prints
  * Shields (1991) Shields, J. C. 1991, AJ, 102, 1314
  * Silge & Gebhardt (2003) Silge, J. D. & Gebhardt, K. 2003, AJ, 125, 2809
  * Singh et al. (1995) Singh, K. P., Bhat, P. N., Prabhu, T. P., & Kembhavi, A. K. 1995, A&A, 302, 658
  * Sparks et al. (1989) Sparks, W. B., Macchetto, F., & Golombek, D. 1989, ApJ, 345, 153
  * Temi et al. (2007) Temi, P., Brighenti, F., & Mathews, W. G. 2007, ApJ, 660, 1215
  * Thomas et al. (2005) Thomas, D., Maraston, C., Bender, R., & Mendes de Oliveira, C. 2005, ApJ, 621, 673
  * Trager et al. (2000a) Trager, S. C., Faber, S. M., Worthey, G., & González, J. J. 2000a, AJ, 120, 165
  * Trager et al. (2000b) —. 2000b, AJ, 119, 1645
  * Trinchieri & di Serego Alighieri (1991) Trinchieri, G. & di Serego Alighieri, S. 1991, AJ, 101, 1647
  * Trinchieri & Goudfrooij (2002) Trinchieri, G. & Goudfrooij, P. 2002, A&A, 386, 472
  * Trinchieri et al. (1997) Trinchieri, G., Noris, L., & di Serego Alighieri, S. 1997, A&A, 326, 565
  * Veron-Cetty & Veron (1986) Veron-Cetty, M.-P. & Veron, P. 1986, A&AS, 66, 335
  * Voit & Donahue (1990) Voit, G. M. & Donahue, M. 1990, ApJ, 360, L15
  * White & Chevalier (1983) White, R. E. & Chevalier, R. A. 1983, ApJ, 275, 69


∎

11institutetext: Dominik Bauer
11email<EMAIL_ADDRESS>
ORCiD: 0000-0002-8781-2882
Tel.: +49-621-3836551

1 Computer Assisted Clinical Medicine, Mannheim Institute for Intelligent
Systems in Medicine, Medical Faculty Mannheim, Heidelberg University, Germany
2 Carl E. Ravin Advanced Imaging Labs, Department of Radiology, Duke
University, United States
* These two authors share senior authorship.

# Generation of annotated multimodal ground truth datasets for abdominal
medical image registration

Dominik F. Bauer1 Tom Russ1 Barbara I. Waldkirch1 Christian Tönnes1 William P.
Segars2 Lothar R. Schad1 Frank G. Zöllner1∗ Alena-Kathrin Golla1∗

(Received: date / Accepted: date)

###### Abstract

Purpose Sparsity of annotated data is a major limitation in medical image
processing tasks such as registration. Registered multimodal image data are
essential for the diagnosis of medical conditions and the success of
interventional medical procedures. To overcome the shortage of data, we
present a method that allows the generation of annotated multimodal 4D
datasets.
Methods We use a CycleGAN network architecture to generate multimodal
synthetic data from the 4D extended cardiac‐torso (XCAT) phantom and real
patient data. Organ masks are provided by the XCAT phantom, therefore the
generated dataset can serve as ground truth for image segmentation and
registration. Realistic simulation of respiration and heartbeat is possible
within the XCAT framework. To underline the usability as a registration ground
truth, a proof of principle registration is performed.
Results Compared to real patient data, the synthetic data showed good
agreement regarding the image voxel intensity distribution and the noise
characteristics. The generated T1-weighted magnetic resonance imaging (MRI),
computed tomography (CT), and cone beam CT (CBCT) images are inherently co-
registered. Thus, the synthetic dataset allowed us to optimize registration
parameters of a multimodal non-rigid registration, utilizing liver organ masks
for evaluation.
Conclusion Our proposed framework provides not only annotated but also
multimodal synthetic data which can serve as a ground truth for various tasks
in medical imaging processing. We demonstrated the applicability of synthetic
data for the development of multimodal medical image registration algorithms.

###### Keywords:

CycleGAN Image Registration Image Synthesis Liver Multimodal Imaging

††journal: International Journal of Computer Assisted Radiology and Surgery

##  1 Introduction

###  1.1 Background

Multimodal imaging plays an important part in the diagnosis of cancers, such
as liver cancer tempany2015multimodal . A variety of treatment options are
available for hepatocellular carcinoma (HCC), the sixth most common malignancy
worldwide and the third leading cause of cancer-related deaths
memon2011radioembolization . These include interventional procedures such as
transarterial chemoembolizations (TACE) waldkirch2019multimodal  or
radioembolization spahr2019multimodal . The treatment planning benefits from
using multimodal registration to combine pre- and intrainterventional data.
Each imaging modality has strengths and weaknesses. Image registration enables
the fusion of complementary information of each modality.

The lack of convenient ground truth data is a major limitation in the field of
medical image segmentation and registration zollner2020image . The generation
of organ masks for segmentation requires labor-intensive manual annotation.
For the development of image registration algorithms (especially for non-rigid
image registration methods) and validation of registration accuracy, the
ground truth is generally not available zollner2020image . This is because the
patient positioning in-between scans usually cannot be reproduced,
particularly in the case of multimodal imaging. In the abdomen, the variable
content of the bladder and bowel and additional patient motion like
respiration and heartbeat further exacerbate the problem.

###  1.2 Related Work

To evaluate registration results or to train deep learning registration
approaches, either anatomical multi-label organ masks or landmarks are
required. However, generating these labels is labor-intensive, subjective or
even impractical for large datasets. Pluim et al. introduced a semi-automatic
framework that allows the creation of a large number of high-quality landmarks
with minimal user effort pluim2016truth . Established ground truth datasets
for the validation of image registration are usually only available for brain
imaging. For instance, the simulated brain database (BrainWeb) provides
simulated MRI imaging sequences (T1-weighted, T2-weighted, and proton density)
cocosco97brainweb . The images are perfectly aligned, since they are
calculated from the same model.

Image synthesis is able to reduce multimodal registration problems to
monomodal problems by first converting one modality into the other. Modality
reduction has shown improvements in registration accuracy for the brain
roy2014mr  and the pelvis cao2017dual .

For MRI-only radiotherapy planning Wolterink et al. demonstrated feasible
results using a CycleGAN approach for MRI-to-CT translation and showed that
training with unpaired images is superior to training with paired images
Wolterink2017 .

Using the digital 4D extended cardiac‐torso (XCAT) phantom Segars2010  instead
of patient images for image synthesis is beneficial, because organ masks and
motion displacement fields are provided by the phantom. Tmenova et al.
presented a CycleGAN to synthesize X-ray angiograms from the XCAT phantom,
which proved to be useful as a data augmentation strategy tmenova2019cyclegan
. Russ et al. synthesized abdominal CT images using a CycleGAN and the XCAT
phantom russ2019synthesis . They showed that a vessel segmentation network
trained on a combination of real CT and synthetic CT images achieved a
superior performance compared to a network trained only on real data.
Analytical models paganelli2017tool ; wissmann2014mrxcat  and a GAN approach
abbasi20204d  to transform the CT XCAT phantom into cardiac or abdominal MRI
images have been developed. To our knowledge, no multimodal registration
ground truth dataset of the abdomen created from the same XCAT or digital
phantom has been reported in the literature.

###  1.3 Contribution

Our approach to bypass the lack of ground truth data in image registration and
segmentation is the generation of a mulitmodal synthetic dataset from the XCAT
phantom. In this work, we focus on multimodal image registration. The
synthesis is performed via CycleGAN networks, which were seperately trained
for each modality. To improve the preservation of high-contrast structures, we
extend the CycleGAN generator loss with an intensity loss and a gradient
difference loss.

Interventions are often monitored via cone beam CT (CBCT), whereas CT and MRI
images are taken for diagnosis beforehand to assist the navigation during
intervention waldkirch2019multimodal . Thus, our multimodal dataset consists
of T1-weighted MRI, CT, and CBCT images. We use XCAT data in the inhaled and
exhaled motion state to synthesize images. Since the same XCAT phantom is used
as the starting point for all modalities, the resulting multimodal synthetic
data is perfectly co-registered. Displacement fields for respiratory movements
and segmentation masks for all organs are provided by the XCAT phantom.
Therefore, it serves as a ground truth dataset for registration.

To demonstrate the utility of the multimodal dataset for the optimization of
registration algorithms, we evaluate a multimodal non-rigid registration for
varying parameter settings. We focus on the registration of the liver,
however, the registration quality can be assessed for any other organ.

##  2 Materials and Methods

###  2.1 Image Synthesis Framework

A schematic of our simulation framework is shown in Fig. 1. Starting from the
CT XCAT phantom, CBCT and MRI XCAT versions are generated by applying a FOV
Mask or by simulating the volume interpolated breathold exam (VIBE) signal
equation paganelli2017tool , respectively. Organ masks for each modality are
extracted from the phantoms. Images are synthesized via CycleGAN networks
using the XCAT phantom as input. CycleGANs learn the mapping between two
domains $X$ and $Y$ given unpaired training samples $x\in X$ and $y\in Y$
zhu2017unpaired . The mapping functions $G:X\rightarrow Y$ and $F:Y\rightarrow
X$ are called generators. Two discriminators $D_{\mathrm{X}}$ and
$D_{\mathrm{Y}}$ aim to distinguish between real images and generated images.
Fig. 2 shows the complete CycleGAN network architecture for the XCAT and CT
image domain. Synthetic CT, CBCT, and MRI images are created via separately
trained CycleGAN networks. The cycle consistency loss $L_{\mathrm{cyc}}(G,F)$
enforces forward and backward consistency for the generators, i.e.
$F(G(x))\approx x$ and $G(F(y))\approx y$. With a least squares generative
adversarial loss $L_{\mathrm{adv}}(G,F,D_{\mathrm{X}},D_{\mathrm{Y}})$, the
generators are trained to generate images that cannot be distinguished from
real images by the discriminator. The discriminators are 70 x 70 PatchGANs,
which are trained with a least squares generative adversarial loss function.
For the generators we use a Res-Net architecture with an encoding stage, 9
residual blocks and a decoding stage.

Figure 1: Schematic of the simulation framework. Starting point is the CT
XCAT, from which CBCT and MRI versions are derived. Synthetic CT, CBCT, and
MRI images are created via separately trained CycleGAN networks. Organ masks
can be obtained from the XCAT phantoms. Patient images that are used to train
the CycleGANs are shown on the right hand side. Figure 2: CycleGAN network
architecture: The generators $G_{\mathrm{XCAT\rightarrow CT}}$ and
$F_{\mathrm{CT\rightarrow XCAT}}$ map images from the XCAT domain to the CT
domain and vice versa. CycleGAN networks for MRI and CBCT images are trained
analogously.

###  2.2 Training and Loss Functions

Training is performed using the Adam optimizer with a learning rate of 0.0002.
The network is trained with 256x256 pixel image patches and a batch size of 4.
For each axial slice one random patch is extracted. We train the networks for
150.000 steps.

To enhance the preservation of high-contrast structures, we extend the
generator loss with an intensity loss and a gradient difference loss:

$\displaystyle L_{\mathrm{int}}(G,F)=||(G(x)-x)||_{1}+||(F(y)-y)||_{1},$ (1)

$\displaystyle\begin{aligned}
L_{\mathrm{gdl}}(G,x)=\sum_{i,j}\;&||x_{i,j}-x_{i-1,j}|-\;|G(x)_{i,j}-G(x)_{i-1,j}||^{2}\\\
+\;&||x_{i,j}-x_{i,j-1}|-\;|G(x)_{i,j}-G(x)_{i,j-1}||^{2}.\end{aligned}$ (2)

The intensity loss preserves the signal intensity of the organs provided by
the XCAT phantom. As shown by Nie et al., the gradient difference loss
prevents blurring and therefore sharpens the synthesized images nie2018medical
. The total generator loss is a combination of the previously defined losses
with different weights:

$\displaystyle\begin{aligned}
L_{\mathrm{gen}}(G,F,D_{\mathrm{X}},D_{\mathrm{Y}})&=L_{\mathrm{adv}}(G,F,D_{\mathrm{X}},D_{\mathrm{Y}})+\;\lambda_{\mathrm{cyc}}\;L_{\mathrm{cyc}}(G,F)\\\
&+\;\lambda_{\mathrm{int}}\;L_{\mathrm{int}}(G,F)+\;\lambda_{\mathrm{gdl}}\;(L_{\mathrm{gdl}}(G,x)+L_{\mathrm{gdl}}(F,y)).\end{aligned}$
(3)

We train three CycleGAN networks for CBCT, CT and MRI with empirically chosen
combinations of weights
$\lambda_{\mathrm{cyc}}/\lambda_{\mathrm{int}}/\lambda_{\mathrm{gdl}}$, of
10/10/5, 10/10/5 and 10/0.4/0.4, respectively. As shown in our previous study,
a combination of the gradient loss and the intensity loss yields the best
results bauer2019synthesis . Further increasing the weighting factors leads to
excessive regularization and thus the networks learns an identity mapping. For
the MRI networks lower over-regularization thresholds are found.

###  2.3 Data

We train our CycleGAN network to map between XCAT phantom data and real
patient data. The goal is to obtain networks that generate realistic looking
synthetic data using the XCAT phantoms as input. In the following paragraphs
we will address the real patient and the XCAT training data separately.

####  2.3.1 Patient Data

The patient training data is retrospectively extracted from our Picture
Archiving and Communication System (PACS). CT and T1-weighted MRI scans are
acquired as part of routine clinical practice before TACE patients undergo a
CBCT-guided TACE intervention. All scans are acquired on whole body clinical
devices (Siemens Healthineers, Forchheim, Germany; CT: Somatom Emotion 16;
CBCT: Artis Zeego; MRI: Magnetom Tim Trio). The MRI images are acquired at 3
Tesla with the VIBE sequence. For each modality the patient images are
resampled to a unified voxel spacing given in Table 1. All scans include the
whole liver and the narrow field of view of the CBCT scan is focused on the
liver. The MRI images include arms, the CT and CBCT images do not. As the XCAT
phantom does not include the patient couch, we removed the patient couch from
the CT patient volumes. For CBCT and MRI no couch is visible in the patient
images.

The image intensities are windowed to the ranges given in Table 1. For CT and
CBCT a fixed window was used. Since MRI intensities vary widely from image to
image, the 10th and 90th percentile of each volume (whole 3D matrix) was used
for windowing. For training, a linear intensity transformation is applied to
transform the intensities from the windowing interval to [-1,1]. Normalization
of training data is a crucial step in improving training performance,
regardless of the normalization method used jacobsen2019analysis .

Table 1: Training data statistics. For CT and MRI the number of slices per
image vary in the given interval.

| Parameter | Resolution (x/y/z) [mm] | Windowing | Volumes | Volume size | Arms | Age
---|---|---|---|---|---|---|---
CT Patient | 100-130 kVp | 1/1/2 | [-1024, 1500] HU | 22 | $512\times 512\times[52,151]$ | no | 66 $\pm$ 9
CT Phantom | 90-120 kVp | 1/1/2 | [-1024, 1500] HU | 56 | $512\times 512\times[80,124]$ | no | 51 $\pm$ 14
CBCT Patient | 93-124.7 kVp | 0.486/0.486/0.486 | [-1024, 2000] HU | 24 | $512\times 512\times 386$ | no | 67 $\pm$ 10
CBCT Phantom | 90-120 kVp | 0.486/0.486/0.486 | [-1024, 2000] HU | 56 | $512\times 512\times 386$ | no | 51 $\pm$ 14
MRI Patient | 3 T | 1/1/3 | [10th , 90th] percentile | 24 | $330\times 450\times[48,93]$ | yes | 67 $\pm$ 10
MRI Phantom | 3 T | 1/1/3 | [10th , 90th] percentile | 56 | $330\times 450\times[59,88]$ | yes | 51 $\pm$ 14

####  2.3.2 XCAT Phantom Data

The XCAT model provides highly detailed whole-body anatomies. Organ masks can
be easily obtained within the XCAT framework. Since CycleGANs maintain the
geometry provided by the XCAT, the organ masks can be used as segmentation
masks in the synthesized images. The phantom includes female and male models
for varying ages. The heart beat and respiratory motions can be simulated and
displacement fields of these motions can be generated. The anatomy and motion
can be adapted by various parameters. This allows the creation of highly
individual patient geometries. For the XCAT training data we generate one XCAT
volume per XCAT model for each modality with 56 different models of varying
ages. The XCATs include the whole liver and are generated with the same voxel
spacing, windowing and normalization as the resampled patient data. Arms are
included only in the MRI XCATs.

The XCAT phantom provides attenuation coefficients for all organs. We vary the
simulated tube energy of the CBCT and CT phantoms from 90-120 keV in steps of
5 keV. This leads to a variation of attenuation coefficients in the phantoms.
Afterwards, those are transformed into Hounsfield Units. To obtain CBCT and
MRI XCAT data, we need to convert the CT XCAT. For the CBCT XCAT we apply a
field of view mask obtained from the patient CBCTs, which is centered on the
liver. For the MRI phantoms we replace the attenuation coefficients for each
organ with simulated MRI values using the signal equation for the VIBE
sequence. It ensures that the MRI signal is initialized with realistic values
matching the MRI training data. This enables us to use the aforementioned
intensity and gradient loss for the generation of synthetic MRI images, since
the transformation with the CycleGAN is now monomodal. The signal intensities
($SI$) for the VIBE sequence in terms of acquisition parameters repetition
time TR, echo time TE, and flip angle $\alpha$ and tissue-specific $T1$, $T2$
relaxation times, and proton density $\rho$ is given by:

$\displaystyle
SI=\frac{\rho\sin{\alpha}(1-\exp{-\frac{TR}{T1}})}{(1-\cos{\alpha}\exp{\frac{-TR}{T1}})}\exp{\frac{-TE}{T2}}.$
(4)

We calculate the MRI intensity for all 44 abdominal organs present in the
XCAT. The imaging parameters $TE=4.54$ ms, $TR=7.25$ ms, and
$\alpha=10^{\circ}$ are obtained from the patient VIBE scans. The values for
the proton density $\rho$ are taken from wissmann2014mrxcat . T1 and T2
relaxation times for 3 T for blood and the spinal cord are obtained from
stanisz2005t1  and the rest from de2004mr . For organs with no available T1,
T2 or $\rho$ we use values of similar organs. To simulate some organ
variability, we randomly vary T1, T2, and $\rho$ by $\pm 5\,\%$ using a
uniform distribution.

###  2.4 Evaluation Metrics

Quantification of the synthetic image quality is difficult, since there are no
corresponding real images for comparison tmenova2019cyclegan . Therefore,
metrics that require a one-to-one correspondence like the mean absolute error
(MAE) cannot be calculated between synthetic and real images. Instead, we
calculate one-to-one corresponding metrics between the synthetic images and
the XCATs, to investigate the magnitude of change from the XCAT phantoms. Real
patient images and synthetic images are then compared by assessing their noise
characteristics and voxel intensity distributions.

####  2.4.1 Synthetic vs. XCAT

The axial slices of the synthetic CT volumes are compared to the corresponding
axial slices of the XCAT volumes with respect to anatomical accuracy. The MAE
is calculated to assess the change of the intensity values. We exclude the
background for the calculation of the MAE. The similarity of structure and
features is evaluated using structural similarity index measure (SSIM) and
feature similarity index measure (FSIM) Wang2004 ; Zhang2011 . Additionally,
we calculate the edge preservation ratio (EPR) and edge generation ratio (EGR)
Chen2016 ; russ2019synthesis .

####  2.4.2 Synthetic vs. Real Patient

Regarding realistic noise characteristics and intensity distribution, the 3D
synthetic volumes are compared to the 3D patient volumes. For the noise
characteristics, only liver voxels are considered. Limiting the noise
considerations to the liver is reasonable, since the liver is a large and
mostly homogeneous organ. We manually segmented the liver in 4 patients for
each modality. The liver segmentations for the 56 synthetic images are
provided by the XCAT phantom. The noise texture is evaluated using an
estimation of the radial noise power spectrum (NPS). The radial NPS of the
synthetic and patient images is compared by calculating the Pearson
correlation coefficient, further called the NPS correlation coefficient (NCC)
russ2019synthesis . In addition to noise texture, we calculate the noise
magnitude (NM), i.e. the standard deviation of the liver voxel intensities.

Furthermore, intensity distribution histograms of patient and synthetic images
are calculated. To quantify their similarity, the Pearson correlation
coefficient between them is calculated (HistCC).

###  2.5 Proof of Principle Registration Evaluation

We perform a proof of principle image registration to demonstrate the
feasibility of the multimodal dataset for evaluation and thus development of
registration algorithms. Our goal is to investigate different parameter
settings to optimize the registration result. We implement the registration in
Python 3.5 with SimpleITK 1.2.4. A non-rigid B-spline transform with a
gradient descent optimizer, a learning rate of 1 and a maximum of 300
iterations is used. Three different registration metrics are considered,
namely Mattes Mutual Information (MMI), Normalized Correlation (NC), and Mean
Squares (MS). For the MMI, 50 histogram bins are used. The MS metric is only
used for the monomodal CT to CT registration, since it is not suited for
multimodal images. Additionally, we vary the spacing of the B-spline control
points from 50 mm to 150 mm in steps of 20 mm. For the multimodal and
monomodal registrations this results in 12 and 18 different parameter
settings, respectively.

The registration is performed on the synthetic data from all 56 XCAT models.
We registered the CT, MRI and CBCT images in the inhaled state to the CT image
in the exhaled state. To evaluate the registration, we take advantage of the
liver organ masks obtained from the XCAT phantoms. The veins and arteries
inside the liver are included in the liver mask by using a morphological
closing operation. We apply the registration transform to the liver masks in
the inhaled state and compare the result to the CT liver mask in the exhaled
state. The overlap of the two masks is assessed by calculating the Dice
similarity coefficient (DSC).

##  3 Results

###  3.1 Synthetic Images

First, we consider the metrics that compare the synthetic images with the XCAT
phantoms shown in the upper half of Table 2. The FSIM and SSIM indicate that
image structures and features were well preserved in the CT and CBCT images,
whereas the synthetic MRIs showed little structural and feature similarity to
the XCATs. Regarding edges, the EPR is similar for all modalities, whereas the
EGR is largest for the CBCT images. The MAE is slightly larger than the NM
(synthetic) for every modality. The MAE for CBCT is more than twice as high as
the MAE for CT.

Table 2: Image quality metrics for the evaluation of the synthetic images. | CBCT | CT | MRI
---|---|---|---
SSIM | 0.85 $\pm$ 0.05 | 0.94 $\pm$ 0.02 | 0.59 $\pm$ 0.04
FSIM | 0.82 $\pm$ 0.03 | 0.82 $\pm$ 0.02 | 0.51 $\pm$ 0.02
EPR | 0.47 $\pm$ 0.06 | 0.43 $\pm$ 0.04 | 0.40 $\pm$ 0.03
EGR | 3.0 $\pm$ 0.7 | 1.9 $\pm$ 0.3 | 1.7 $\pm$ 0.2
MAE | 109 $\pm$ 14 | 51 $\pm$ 16 | 37 $\pm$ 6
NCC | 0.997 $\pm$ 0.001 | 0.980 $\pm$ 0.010 | 0.86 $\pm$ 0.04
NM (Synthetic) | 52 $\pm$ 13 | 39 $\pm$ 5 | 25 $\pm$ 3
NM (Patient) | 60 $\pm$ 16 | 39 $\pm$ 19 | 22 $\pm$ 5
HistCC | 0.994 $\pm$ 0.003 | 0.999 $\pm$ 0.002 | 0.94 $\pm$ 0.03

Secondly, we compare the synthetic images to the patient images. The two right
columns of Fig. 1 show axial synthetic and patient slices of each modality.
Qualitatively, the style of the synthesized images is in good agreement with
the real patient images. To quantify this observation we compared the noise
characteristics and voxel intensity distribution of the synthetic images to
the patient images, the results are listed in the lower half of Table 2. A
high NCC for all modalities indicates that the noise texture was emulated
realistically, albeit the NCC is slightly smaller for the synthetic MRI
images. For all modalities, the NM (synthetic) is in excellent agreement with
the NM (patient). In Fig. 3 the intensity histograms are shown. In general,
the synthetic intensity distributions match the patient intensity
distributions nicely. This is underlined by the overall high HistCC values in
Table 2. However, for CT and CBCT the soft tissue peaks are modeled a bit too
narrowly. The lung tissue peak is shifted towards higher CT numbers for the
CT. In the MRI, the soft tissues is slightly underrepresented.

###  3.2 Proof of Principle Registration

The DSC for the evaluation of the proof of principle registration is shown in
Fig. 4. The monomodal CT to CT registration yielded good results for all
registration metrics and grid point spacings, with the best result for MMI
with 50 mm grid point spacing. For CBCT, the MMI again worked well, whereas
the registrations using the NC mostly failed. The best results were again
obtained with MMI and a grid point spacing of 50 mm. For MRI, the
registrations with MMI and NC yielded similar results with the best result
obtained for NC with a grid point spacing of 150 mm. Overall, the monomodal CT
to CT registration achieved the best results.

|  |
---|---|---
Figure 3: Intensity histograms of the patient and synthetic images averaged over all volumes. Note that the background peaks are cropped. |  |
---|---|---
Figure 4: DSC for the proof of principle registrations with 56 data points
each. The mean is marked as a ”+” and the whiskers indicate the 10th and 90th
percentile. The dashed horizontal line shows the mean pre-registration DSC.

Coronal views of the registration results for the best settings of each
modality are visualized in Fig. 5. The registered images in the middle row
show a large similarity to the ground truth. This observation is further
supported by the overlaid liver contours. The post-registration liver contour
(yellow) is in high agreement with the ground truth liver contour (red).

|  Pre-Registration |  Post-Registration |  Ground Truth Image
---|---|---|---
CBCT |  |  |
CT |  |  |
MRI |  |  |
Figure 5: Pre- and Post-registered images and the corresponding ground truth.
Red contours indicate the ground-truth boundaries of the liver (target). Blue
and yellow contours represent the boundaries of the liver before and after
deformation respectively.

##  4 Discussion

The image quality metrics in Table 2 demonstrate that our framework provides
reaslistic multimodal image data. Low SSIM and FSIM for the MRI images
indicate, that image values of the homogeneous organs in the MRI XCAT phantoms
needed to be altered more strongly by the networks in comparison to the CT
XCATs. The lower SSIM, FSIM, and EPR for MRI are likely a result of lower
weighting of the gradient loss and intensity loss.

The ratio of MAE to NM (synthetic) is 2.1, 1.3, and 1.5 for CBCT, CT, and MR,
respectively. Assuming normally distributed noise, the ratio of MAE to NM is
only about 0.8 geary1935ratio . This means that the MAE cannot be attributed
to noise alone. The large MAE for the CBCT images compared to the CT images
could be due to the introduction of metal artifacts, as the patient CBCTs
showed metal artifacts in the liver caused by medical instruments. This is
supported by a large EGR for CBCT.

For all modalities realistic noise texture and magnitude was achieved.
Additionally, the voxel intensity distribution was modeled adequately. Most of
the discrepancies between the patient and synthetic histograms in Fig. 3 can
be explained by inspecting the XCAT phantoms. The deviation of the CT lung
peaks (synthetic -780 HU, patient -835 HU) can be explained by an
overestimated initial lung value of -760 HU given by the XCAT. The narrow soft
tissue peaks for CT and CBCT could be due to insufficient variation in organ
attenuation coefficients. The under representation of soft tissue in the
synthetic MRI is due to the body size of the patients and XCATs. We found that
in the MRI patient dataset 66.5 % of the image voxels show the body, whereas
for the MRI XCAT dataset, it is only 46.5 %. A rather large HistCC of 0.94
$\pm$ 0.03 was still achieved, since this under representation has only a
minor effect on the correlation between the histograms. We prepared the XCAT
data such that it matches the patient dataset as good as possible, see Table
1. In the future we will consider the patient body size beforehand and adjust
the XCAT body size accordingly.

An important requirement for using synthetic data to evaluate or train
registration algorithms is a realistic respiratory deformation model. The XCAT
framework allows the respiratory rate, the amount of diaphragmatic motion, the
amount of chest expansion, and the amount of cardiac motion due to respiration
to be varied. In addition, curves controlling the diaphragm motion and chest
expansion can be individualized. However, as the breathing patterns of
patients are complex and diverse, the XCAT model may still lack
generalizability.

The results of the proof of principle registration demonstrate that the
synthetic dataset can be used to evaluate different registration algorithms.
We were able to evaluate the performance of different registration algorithms
and to fine tune parameter settings. It is noticeable that the registration
for CT and CBCT works better for small B-spline grid point spacings, while it
is the opposite for MRI. Further studies are needed to assess whether this is
due to the smaller image size of the MRI data and thus the total number of
B-spline control points or to other imaging characteristics of MRI and CT.
Computation time can be measured and taken into consideration. For example,
registrations with smaller grid point spacings take much longer. Thus,
choosing MMI with 150 mm for CT and CBCT registrations might be a reasonable
trade-off, as the registration quality is only slightly lower, while the
registration time is substantially reduced. The availability of organ masks
enabled a rather simple registration evaluation.

##  5 Conclusion

The presented simulation framework can be used to extend small datasets by
transferring the style of the dataset onto the geometry given by the XCAT
phantom. The obtained datasets can serve as a ground truth for image
registration. A multimodal dataset consisting of T1-weighted MRI, CT and CBCT
images was created and used to demonstrate the refinement and evaluation of
multimodal image registration algorithms. In the future, the framework will be
extended to other modalities, such as T2-weighted MRI or PET, which can
further boost the performance of multimodal methods. An extension to other
body regions, such as the thorax or pelvis, is also possible. Synthetic images
over larger body regions are especially interesting for whole body
segmentation. Expansion of datasets using this method provides a promising
tool to overcome the dearth of medical training data.

###### Acknowledgements.

We gratefully acknowledge the support of NVIDIA Corporation with the donation
of the NVIDIA Titan Xp GPU used for this research.
Funding This research project is part of the Research Campus M2OLIE and funded
by the German Federal Ministry of Education and Research (BMBF) within the
Framework ’Forschungscampus - Public-Private Partnership for Innovation’ under
the funding code 13GW0388A.
Conflict of Interest The authors declare that they have no conflict of
interest.
Human and animal rights We are in compliance with ethical standards.
Ethics approval Ethical approval for this study was obtained from the ethics
committee of the medical faculty Mannheim (2016-863R-MA).

## References

  * (1) Abbasi-Sureshjani, S., Amirrajab, S., Lorenz, C., Weese, J., Pluim, J., Breeuwer, M.: 4d semantic cardiac magnetic resonance image synthesis on xcat anatomical model.  arXiv preprint arXiv:2002.07089 (2020)
  * (2) Bauer, D.F., Schnurr, A.K., Russ, T., Goerttler, S., Schad, L.R., Zoellner, F.G., Chung, K.: Synthesis of ct images using cyclegans: Enhancement of anatomical accuracy.  In: Proc. Int. Conf. Med. Imag. Deep Learning (2019)
  * (3) Cao, X., Yang, J., Gao, Y., Guo, Y., Wu, G., Shen, D.: Dual-core steered non-rigid registration for multi-modal images via bi-directional image synthesis.  Med. Image Anal. 41, 18–31 (2017)
  * (4) Chen, L., Jiang, F., Zhang, H., Wu, S., Yu, S., Xie, Y.: Edge preservation ratio for image sharpness assessment.  In: 2016 12th World Congress on Intelligent Control and Automation (WCICA), pp. 1377–1381. IEEE (2016)
  * (5) Cocosco, C.A., Kollokian, V., Kwan, R.K.S., Pike, G.B., Evans, A.C.: Brainweb: Online interface to a 3d mri simulated brain database.  NeuroImage 5, 425 (1997)
  * (6) De Bazelaire, C.M., Duhamel, G.D., Rofsky, N.M., Alsop, D.C.: Mr imaging relaxation times of abdominal and pelvic tissues measured in vivo at 3.0 t: preliminary results.  Radiology 230(3), 652–659 (2004)
  * (7) Geary, R.C.: The ratio of the mean deviation to the standard deviation as a test of normality.  Biometrika 27(3/4), 310–332 (1935)
  * (8) Jacobsen, N., Deistung, A., Timmann, D., Goericke, S.L., Reichenbach, J.R., Güllmar, D.: Analysis of intensity normalization for optimal segmentation performance of a fully convolutional neural network.  Z. Med. Phys. 29(2), 128–138 (2019)
  * (9) Lin Zhang, Lei Zhang, Xuanqin Mou, Zhang, D.: FSIM: A Feature Similarity Index for Image Quality Assessment.  IEEE Trans. Image Processing 20(8), 2378–2386 (2011)
  * (10) Memon, K., Lewandowski, R.J., Kulik, L., Riaz, A., Mulcahy, M.F., Salem, R.: Radioembolization for primary and metastatic liver cancer.  Semin. Radiat. Oncol. 21(4), 294–302 (2011)
  * (11) Nie, D., Trullo, R., Lian, J., Wang, L., Petitjean, C., Ruan, S., Wang, Q., Shen, D.: Medical image synthesis with deep convolutional adversarial networks.  IEEE Trans. Biomed. Eng. 65(12), 2720–2730 (2018)
  * (12) Paganelli, C., Summers, P., Gianoli, C., Bellomi, M., Baroni, G., Riboldi, M.: A tool for validating mri-guided strategies: a digital breathing ct/mri phantom of the abdominal site.  MBEC 55(11), 2001–2014 (2017)
  * (13) Pluim, J.P., Muenzing, S.E., Eppenhof, K.A., Murphy, K.: The truth is hard to make: Validation of medical image registration.  In: 2016 23rd International Conference on Pattern Recognition (ICPR), pp. 2294–2300. IEEE (2016)
  * (14) Roy, S., Carass, A., Jog, A., Prince, J.L., Lee, J.: Mr to ct registration of brains using image synthesis.  In: Proc. SPIE Med. Imag., vol. 9034, p. 903419. SPIE (2014)
  * (15) Russ, T., Goerttler, S., Schnurr, A.K., Bauer, D.F., Hatamikia, S., Schad, L.R., Zöllner, F.G., Chung, K.: Synthesis of ct images from digital body phantoms using cyclegan.  Int. J. Comput. Assist. Radiol. Surg. 14(10), 1741–1750 (2019)
  * (16) Segars, W.P., Sturgeon, G., Mendonca, S., Grimes, J., Tsui, B.M.W.: 4D XCAT phantom for multimodality imaging research.  Med. Phys. 37(9), 4902–4915 (2010)
  * (17) Spahr, N., Thoduka, S., Abolmaali, N., Kikinis, R., Schenk, A.: Multimodal image registration for liver radioembolization planning and patient assessment.  Int. J. Comput. Assist. Radiol. Surg. 14(2), 215–225 (2019)
  * (18) Stanisz, G.J., Odrobina, E.E., Pun, J., Escaravage, M., Graham, S.J., Bronskill, M.J., Henkelman, R.M.: T1, t2 relaxation and magnetization transfer in tissue at 3t.  Magn. Reson. Med. 54(3), 507–512 (2005)
  * (19) Tempany, C.M., Jayender, J., Kapur, T., Bueno, R., Golby, A., Agar, N., Jolesz, F.A.: Multimodal imaging for improved diagnosis and treatment of cancers.  Cancer 121(6), 817–827 (2015)
  * (20) Tmenova, O., Martin, R., Duong, L.: Cyclegan for style transfer in x-ray angiography.  Int. J. Comput. Assist. Radiol. Surg. 14(10), 1785–1794 (2019)
  * (21) Waldkirch, B., Engelhardt, S., Zöllner, F.G., Schad, L.R., Wolf, I.: Multimodal image registration of pre-and intra-interventional data for surgical planning of transarterial chemoembolisation.  In: Proc. SPIE Med. Imag., vol. 10951, p. 109512U. SPIE (2019)
  * (22) Wang, Z., Bovik, A.C., Sheikh, H.R.: Image quality assessment: From error measurement to structural similarity.  IEEE Trans. Image Processing 13(4), 600 – 612 (2004)
  * (23) Wissmann, L., Santelli, C., Segars, W.P., Kozerke, S.: Mrxcat: Realistic numerical phantoms for cardiovascular magnetic resonance.  J. Cardiovasc. Magn. Reson. 16(1), 63 (2014)
  * (24) Wolterink, J.M., Dinkla, A.M., Savenije, M.H.F., Seevinck, P.R., van den Berg, C.A.T., Išgum, I.: Deep mr to ct synthesis using unpaired data.  In: Proc. Int. Workshop Simul. Synth. Med. Imag., pp. 14–23 (2017)
  * (25) Zhu, J.Y., Park, T., Isola, P., Efros, A.A.: Unpaired image-to-image translation using cycle-consistent adversarial networks.  In: Proc. IEEE Int. Conf. Comput. Vis., pp. 2223–2232 (2017)
  * (26) Zöllner, F.G., Šerifović-Trbalić, A., Kabelitz, G., Kociński, M., Materka, A., Rogelj, P.: Image registration in dynamic renal mri—current status and prospects.  Magn. Reson. Mater. Phy. 33, 33–48 (2020)