© ESO, 2010

1. Introduction

Starting with the discovery of the gamma-ray burst-supernova (GRB-SN) association GRB 980425/SN 1998bw (Galama et al. 1998; Iwamoto et al. 1998; Sollerman et al. 2000), there is by now convincing evidence that the progenitors of long GRBs are massive stars exploding as type Ic SNe (for a review see Woosley & Bloom 2006). Within this picture, the optical light observed after a long GRB is the superposition of the afterglow light, a supernova component, and light from the underlying host galaxy (plus potential additional radiation components at very early times, which we will not consider here). Phenomenologically, this immediately unveils two observing strategies to reveal a massive-star origin of a GRB: (i) via the detection of a late-time SN bump (method i) in the optical light curve (e.g., Reichart 1999; Galama et al. 2000; Dado et al. 2002; Zeh et al. 2004) and (ii) by the spectroscopic confirmation (method ii) of associated supernova light (the best case so far being GRB 030329: Hjorth et al. 2003; Kawabata et al. 2003; Matheson et al. 2003; Stanek et al. 2003). In addition to both observing strategies, there are two further methods by which a massive-star origin can be revealed. Some afterglow spectra showed blue-shifted absorption line systems (method iii), which can be understood as signatures from the expanding pre-explosion wind escaping from the GRB progenitor (e.g., Mirabal et al. 2003; Schaefer et al. 2003; Klose et al. 2004; Starling et al. 2005; Berger et al. 2006; Fox et al. 2008; Castro-Tirado et al. 2010). Some authors, however, notice that several properties of the putative blue-shifted absorption line systems disagree with the expectations from Wolf-Rayet (WR) winds, e.g., line widths, ionisation levels and metallicities (Chen et al. 2007; Prochaska et al. 2007; Fox et al. 2008). Finally (method iv), the circumburst medium determines the spectral and temporal evolution of the afterglow, allowing us to discern between a constant-density medium, n(r) = const., and a wind medium, n(r) ∝ r^-2 (Sari et al. 1998; Chevalier & Li 1999, 2000). This method was successfully applied in, e.g., Starling et al. (2009) and Curran et al. (2010).

Naturally, these various approaches have their observational advantages and disadvantages. While method (ii) can provide the strongest observational evidence for a massive-star origin of the GRB under consideration, it can only be applied to the nearest and hence brightest events up to a redshift of about 0.5. Even 13 years after the first discovery of an afterglow (Costa et al. 1997; van Paradijs et al. 1997), i.e., after more than 500 GRBs with detected (X-ray, optical, radio) afterglow light¹, secure evidence for a spectroscopically associated SN was only reported for roughly 1% of all events (GRBs 980425: Galama et al. 1998; 030329: Hjorth et al. 2003; Kawabata et al. 2003; Matheson et al. 2003; Stanek et al. 2003; 031203: Malesani et al. 2004; 060218: Ferrero et al. 2006; Mirabal et al. 2006; Modjaz et al. 2006; Pian et al. 2006; Sollerman et al. 2006; 081007: Della Valle et al. 2008; 100316D: Chornock et al. 2010; Starling et al. 2010). Contrary to this approach, method (i) can basically reveal a SN component up to a redshift of about 1 (Zeh et al. 2004), assuming a non-extinguished SN 1998bw as a template (the most distant SN bump was found for GRB 000911 at z = 1.06; Masetti et al. 2005). At notably higher redshifts a GRB-SN becomes too faint to be discovered even with 8m-class optical telescopes because of line blanketing (e.g., Filippenko 1997).

In principle, methods (iii) and (iv) do not have redshift constraints, because both rely on the observation of the afterglow and not on the (expected) SN component. In addition, method (iv) splits into different approaches and basically works for the optical and X-ray band in the same way. Moreover, it can even be applied to the most distant GRBs, which are already affected by Lyman dropout in the optical bands. Its main disadvantage is that it usually requires substantial observational efforts. In particular, it relies on the measurement of the light-curve evolution, including the determination of the spectral energy distribution (SED).

In this paper, applying method (iv), we use late-time data of afterglows detected by Swift to tackle the question of the preferred density profiles in a statistical sense. Our goal is to combine all publicly available optical data with Swift/XRT data in order to determine the corresponding circumburst density profile, n(r) ∝ r^−k. Qualitatively, this splits into either a constant-density medium (hereafter referred to as the interstellar medium or ISM) (k = 0) or a free wind profile (k = 2). While the complex wind history of an evolved massive star might produce density profiles different from the ideal case (k = 0,2; Crowther 2007), in general the data do not allow for a more accurate determination of k, but only to distinguish between these two cases.

It should be stressed that the intention of our study is not to provide ultimate conclusions on the density profiles found for individual bursts. Instead, it is meant as a statistical approach using bursts with the best available X-ray as well as optical data. The questions we want to address are: (1) What is, in a statistical sense, the preferred circumburst density profile? (2) What is the ratio between events with ISM and with wind profiles? (3) What does this tell us about the typical radius of the wind termination shock that is expected to exist in the ambient medium surrounding a long burst GRB progenitor?

Throughout the paper we use the convention ${\mathit{F}}_{\mathit{\nu }}{}^{\mathrm{\right(}}\mathit{t}^{\mathrm{\left)}}\mathrm{\propto }{\mathit{t}}^{\mathrm{-}\mathit{\alpha }}{\mathit{\nu }}^{\mathrm{-}\mathit{\beta }}$ for the flux density, where α is the temporal slope and β is the spectral slope. All errors are 1σ uncertainties unless noted otherwise.

2. Data selection

2.1. Data gathering

The optical data were taken from a photometric database maintained and updated by one of the co-authors (D.A.K.). In addition, we added data for GRB 090726 from Šimon et al. (2010), Fatkhullinet al. (2009), Haislip et al. (2009), Kelemen (2009), Landsman & Page (2009), Sakamoto et al. (2009) and Volnova et al. (2009). The properties of the optical afterglow sample, the data gathering and the deduced SEDs of the afterglows are discussed in Kann et al. (2006, 2008, 2010). Furthermore, we compared their SED results with the work of Schady et al. (2007, 2010). The latter authors added X-ray data to model the spectral energy distributions from roughly 1 eV to 10 keV. Finally, to create a denser light-curve coverage, all non-R_C band data in each light curve were shifted to the R_C band. In doing this, we used the colours of the SED, assuming no spectral evolution and omitting data where clear colour evolution was evident.

The X-ray data were retrieved from the Swift data archive and the light curves from the Swift light-curve repository (version February 2010) updated and maintained by Evans et al. (2007, 2009). Following Nousek et al. (2006), we reduced the data with the software package HeaSoft 6.6.1² together with the calibration file version v011³. Furthermore, we applied the methods detailed in Moretti et al. (2005), Romano et al. (2006), and Vaughan et al. (2006) to reduce pile-up affected data. We extracted SEDs at different epochs with approximately 500 background subtracted counts, to check for spectral evolution. If the properties of the SED (spectral slope and absorption) did not evolve, we extracted a new SED from the maximum possible time interval. In addition, we used the Tübingen absorption model by Wilms et al. (2000) and their interstellar-medium metal abundance template. The Galactic absorption was fixed to the weighted mean based on Kalberla et al. (2005). We included the Chandra light curve of GRB 051221A from Burrows et al. (2006).

2.2. Sample definition and light-curve fitting

Among all Swift GRBs observed until September 2009, we selected from our photometric database, which contains long and short bursts, those 90 bursts that have an optical and an X-ray afterglow as well as a measured spectroscopic or a photometric redshift. From these we selected those bursts with the best-sampled optical and X-ray afterglow light curves as follows.

First, we required information on the spectral slope of the afterglow in the X-ray band. We derived this information from publicly available data. In the optical bands, however, multi-band data are usually not available. In these cases we used the results from Kann et al. (2008, 2010) and Schady et al. (2007, 2010), with the latter being based on optical-to-X-ray SED fits.

Second, we required that an afterglow light curve can be fitted with a multiply broken power-law and is not dominated by flares or bad data sampling. After excluding time intervals affected by flares, we fitted the light curves with a smoothly broken power-law of the order m (see Appendix A for its definition) with a Simplex and a Levenberg-Marquardt algorithm (Press et al. 2007). Furthermore, we transformed R_c-band light curves to flux densities with the zero-point definition of Bessell (1979). For the X-ray regime we followed Gehrels et al. (2008); the flux density, F_ν,x, in μJy at the frequency ν_x is then given by $\begin{array}{ccc}{\mathit{F}}_{\mathit{\nu ,}\mathrm{x}}\mathrm{=}\mathrm{4.13}\mathrm{\times }{\mathrm{10}}^{\mathrm{11}}\frac{\mathrm{(}\mathrm{1}\mathrm{-}{\mathit{\beta }}_{\mathrm{x}}\mathrm{)}{\mathit{F}}_{\mathrm{x}}}{\mathrm{(}\mathrm{10}\hspace{0.17em}\mathrm{keV}{\mathrm{)}}^{\mathrm{1}\mathrm{-}{\mathit{\beta }}_{\mathrm{x}}}\mathrm{-}\mathrm{(}\mathrm{0.3}\hspace{0.17em}\mathrm{keV}{\mathrm{)}}^{\mathrm{1}\mathrm{-}{\mathit{\beta }}_{\mathrm{x}}}}{\mathit{E}}_{\mathrm{x}}^{\mathrm{-}{\mathit{\beta }}_{\mathrm{x}}}\mathit{,}& & \end{array}$where β_x is the spectral slope in the X-ray band, F_x is the measured flux in the 0.3−10 keV range in units of erg cm^-2 s^-1 and the reference energy E_x is given in keV. For all bursts we chose the logarithmic mean between 0.3 keV and 10 keV as a reference, i.e., E_x = 1.73 keV (ν_x = 4.19 × 10¹⁷   Hz). The numerical constant converts the energy in units of keV to a frequency in units of Hz and the flux density from erg cm^-2 s^-1 Hz^-1 to μJy.

Third, once power-law segments had been defined in the afterglow light curves, we excluded from further studies those bursts where the difference in the late-time decay slopes between the optical and the X-ray band could not be explained within the framework of the fireball model (Table 1); i.e., the difference in decay slopes, α_x − α_opt, was larger than 1/4 within 3σ. Because of this, a strict criterion for the beginning of the late-time evolution of an afterglow cannot be given. Evidence for an observed canonical light curve does not exist in the optical (Kann et al. 2008, 2010; Panaitescu & Vestrand 2010), but may exist in the X-rays (e.g., Nousek et al. 2006; Zhang et al. 2006; Evans et al. 2009). Therefore, as an operational definition an afterglow is in its late-time phase when its temporal and spectral evolution can be explained by the fireball model from a particular time after the corresponding GRB.

After applying these selection criteria, about half of the 90 bursts had to be rejected owing to bad sampling, poor data quality, flares, and other peculiarities. Another quarter had to be rejected because of  |α_x − α_opt|  > 1/4 within 3σ. All data were then corrected for host extinction in the optical, ${\mathit{A}}_{\mathit{V}}^{\mathrm{host}}$, if needed, and for Galactic and host absorption in the X-ray band, ${\mathit{N}}_{\mathrm{H}}^{\mathrm{host}}$.

In total, 27 afterglows passed these selection criteria. The sample consists of 25 long and one short GRB (GRB 051221A), and one controversial event (GRB 060614) in terms of the short/long classification scheme. Their input data are summarised in Tables B.1 and B.2.

3. Results and discussion

3.1. The circumburst medium density profile

3.1.1. Identifying the circumburst medium

As a first step in the identification of the circumburst medium, we defined nine possible spectral and dynamical regimes (Table 2). The spectral regimes separate according to the position of the cooling frequency with respect to the observer frame, while the dynamical regimes distinguish between a spherical and a jetted evolution.

To derive the most probable density profile into which a GRB jet propagates, we proceeded in the following way. The main criterion was that the result agrees with the optical as well as with the X-ray data. In doing so, we first analysed the closure relations for the nine models (Table 2) in the optical and X-ray bands and selected only those relations (models) which were fulfilled within 3σ. Second, we computed the difference in the decay slopes, α_x − α_opt, and distinguished between the models according to the three possible cases −1/4,0, +1/4 (Table 1), again within 3σ. Third, if possible, we took into account the difference in the spectral slope, β_x − β_opt, which is either 0 or 1/2 (see Table 1). If this criterion could be applied, we required that it is fulfilled within 3σ. In addition we required that the 1σ uncertainty in the spectral or temporal decay slopes was less than 0.2.

The final results of the light curve fits of the 27 bursts considered are presented in Fig. C.1. In addition, we plot here the observed flux-density ratio ${}^{\mathrm{(}}\mathit{F}_{\mathit{\nu ,}\mathrm{opt}}\mathit{/}{\mathit{F}}_{\mathit{\nu ,}\mathrm{x}}{}^{\mathrm{)}}\mathrm{\left(}\mathit{t}\mathrm{\right)}$ (middle panels in Fig. C.1) as a function of time. According to Table 1, the expected flux-density ratio is only allowed to take a certain value depending on the spectral and dynamical regime. We checked if the observed flux-density ratio, ${}^{\mathrm{(}}\mathit{F}_{\mathit{\nu ,}\mathrm{opt}}\mathit{/}{\mathit{F}}_{\mathit{\nu ,}\mathrm{x}}{}^{\mathrm{)}}\mathrm{\left(}\mathit{t}\mathrm{\right)}$ (middle panels in Fig. C.1), agreed with the model(s) that successfully passed the previous three criteria. The allowed parameter space of the flux-density ratio of all considered models (Table 1) is shown as grey box in the middle panels in Fig. C.1. The upper and lower boundary always refer to ν_c ≤ ν_opt with F_ν,opt/F_ν,x = (ν_opt/ν_x)^−p/2 and to ν_c ≥ ν_x with F_ν,opt/F_ν,x = (ν_opt/ν_x)^{−(p − 1)/2}. This criterion came into play when we were unable to distinguish whether the cooling frequency was redward of the optical band or blueward of the X-ray band, in other words when the flux-density ratio agreed either with the lower or upper boundary in Fig. C.1.

3.1.2. Ensemble properties

Combining these different criteria, we found that for about 60% (16/27) of all studied cases the cooling break was between the optical and the X-ray band (Table B.3), i.e., the difference in the decay slope is either  +1/4 or −1/4 (for a constant density profile and a free wind medium, respectively). In total, we could identify the circumburst density profile (ISM or wind) for 25 of the 27 investigated afterglows (Table B.3). However, we identified a wind medium for only six events (GRBs, 050603, 070411, 080319B, 080514B, 080916C and 090323). The other 19 bursts were consistent with an ISM profile except for GRBs 051221A and 060904B.

Our procedure to find physical descriptions of light-curve segments, in other words identifying the spectral regime and the circumburst medium, allowed us to find descriptions of more optical and X-ray afterglows than are presented in the literature (Table B.3). Our results usually agree with the literature (for references see Table B.3); they only differ for GRBs 070802, 080721, 090323, and 090328. The contradiction in the latter three bursts is due to the size of the optical data set. We used the maximum publicly available data set in contrast to Starling et al. (2009) (GRB 080721) and Cenko et al. (2010) (GRBs 090323 and 090328).

Four bursts are of particular interest in our sample:

• (i)

  GRB 051221A (z = 0.546; Soderberg et al. 2006) is a short burst (Burrows et al. 2006), i.e., most likely it originated from the merger of two compact objects (Blinnikov et al. 1984; Paczyński 1986; Goodman 1986; Eichler et al. 1989). Unfortunately, we cannot discern between a wind or an ISM profile because ν_c < ν_x. The cooling frequency was below the optical bands so that neither the optical nor the X-ray data can be used to reveal the circumburst medium except during a post-jet break phase without lateral spreading.

• (ii)

  GRB 060614 (z = 0.125; Della Valle et al. 2006b) is a quite controversial event that does not easily fit into the classical short/long classification scheme (Della Valle et al. 2006a; Fynbo et al. 2006; Gal-Yam et al. 2006; Gehrels et al. 2006; Mangano et al. 2007; Zhang et al. 2007; Kann et al. 2008; Zhang et al. 2009). The difference in the spectral slope, β_x − β_opt = 0.00 ± 0.11 (Table B.3), rules out that the cooling frequency lies between the optical and X-ray bands. Furthermore, we did not find evidence for a wind medium. The optical and X-ray data exclude a wind medium with high confidence. The deviation between the observed and predicted temporal decay slope is  <1.3σ and  >3.5σ for an ISM and a wind medium, respectively. In our sample this burst belongs to a small number of cases where the flux density ratio could be used to distinguish between ν_c > (ν_opt,ν_x) and ν_c < (ν_opt,ν_x). The small observed flux density ratio agrees with ν_c > (ν_opt,ν_x).

• (iii)

  GRB 060904B has a well-defined light curve (Fig. C.1, Table B.2) and SED (Table B.1) in the optical and X-ray bands, respectively. The difference in the decay slopes, α_x − α_opt = 0.20 ± 0.04 (Table B.3), favours an ISM profile with the cooling frequency lying between the optical and the X-ray bands. The difference in the spectral slopes, however, does not support this scenario, β_x − β_opt = 0.00 ± 0.16 (Table B.3). Rather, both afterglow components seem to be in the same spectral regime. Therefore, we could not find a consistent description of the optical and X-ray afterglow. The main reason could be that the optical data are mainly based on preliminary data (see Kann et al. 2010 for details on the data gathering).

• (iv)

  GRB 080319B (z = 0.937; D’Elia et al. 2009) is the only burst in our sample with a photometrically detected supernova component (Bloom et al. 2009; Tanvir et al. 2008). An ISM profile can be ruled out with very high confidence. The difference in the spectral slopes, β_x − β_opt = 0.48 ± 0.12 (Table B.3), favours the cooling break to be in between the optical and the X-ray bands. The deviation between Δα_obs = α_x − α_opt and Δα_predicted is 1σ for a wind medium but 9σ for an ISM profile. The closure relations support this finding. The deviation between the observed and predicted optical decay slope is 0.1σ for a wind medium and 4.6σ for an ISM profile (see also Racusin et al. 2008).

3.1.3. The electron index

Finally, the identification of the light-curve segments allowed us to derive the electron index, p. It is shown in the top panels in Fig. C.1 for every burst. In agreement with other studies (e.g., Panaitescu & Kumar 2001; Shen et al. 2006; Zeh et al. 2006; Starling et al. 2008; Curran et al. 2009; Ghisellini et al. 2009; Curran et al. 2010) we find that the distribution extends from p ~ 2 to p ~ 3, indicating that p is no universal value.

3.2. Constraining the wind – ISM transition zone

Once we had identified the circumburst density profiles, we could investigate if there is observational evidence for the position of the wind termination-shock radius, where the density profile changes from that of a free wind (k = 2) to a density profile with k = 0 (independent of whether it is the shocked wind or the ISM; see Pe’er & Wijers 2006; van Marle et al. 2007). While the theory of a blastwave crossing such a density discontinuity has been worked out (Pe’er & Wijers 2006), finding the corresponding observational signature is difficult. Even though data on several hundred afterglows exist, they do not provide this kind of information in a convincing way (Starling et al. 2008; Curran et al. 2009). This leaves open the question of observationally determining the wind termination-shock radius.

Here we cannot determine the radius of the wind termination shock for any (long) GRB progenitor either. However, we can characterise its position in a statistical sense. The light curves in Fig. C.1 show F_opt/F_x as a function of time. Its first logarithmic derivative, ${}^{\mathrm{(}}\mathit{\alpha }_{\mathrm{x}}\mathrm{-}{\mathit{\alpha }}_{\mathrm{opt}}{}^{\mathrm{)}}{}^{\mathrm{\right(}}\mathit{t}^{\mathrm{\left)}}$, reveals either the maximum time up to which a wind profile is identified or it reveals the minimum time after which an ISM profile agrees with the data. The function ${}^{\mathrm{(}}\mathit{\alpha }_{\mathrm{x}}\mathrm{-}{\mathit{\alpha }}_{\mathrm{opt}}{}^{\mathrm{)}}\mathrm{\left(}\mathit{t}\mathrm{\right)}$ is a smooth function in time in contrast to the fit values, because we used a smoothly broken power law to describe the light-curve evolution. Table 3 summarises the time intervals for which the asymptotic values for α_x − α_opt were reached, i.e. −1/4, 0, or 1/4 within 1σ.

In the observer frame, the radius of the fireball propagating into a free wind medium is (Chevalier & Li 2000) $\mathit{R}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}\mathrm{1.1}\hspace{0.17em}\mathrm{\times }\hspace{0.17em}{\mathrm{10}}^{\mathrm{17}}{\left(\frac{\mathrm{2}\mathit{t}\hspace{0.17em}{\mathit{E}}_{\mathrm{iso}\mathit{,}\mathrm{52}}}{\mathrm{(}\mathrm{1}\mathrm{+}\mathit{z}\mathrm{)}\hspace{0.17em}{\mathit{A}}_{\mathit{\star }}}\right)}^{\mathrm{1}\mathit{/}\mathrm{2}}\mathrm{cm}\mathit{,}$(1)where t is measured in units of days, E_iso is the isotropic equivalent energy in units of 10⁵²   erg and A _⋆  is defined via A = Ṁ_w/4πv_w = 5    ×    10¹¹   A _⋆    g   cm^-1, with Ṁ_w being the mass-loss rate, and v_w the wind velocity. The quantity A _⋆  refers to a mass-loss rate of Ṁ_w = 10^-5   M_⊙   yr^-1 and a velocity of the stellar wind of v_w = 10⁸   cm   s^-1.

Using the aforementioned approach and fixing for simplicity A _⋆  = 1, Table 3 summarises the deduced upper and lower limits on the wind termination-shock radius. The distribution extends over three orders of magnitude (from  ≈ 10^-3 pc to 1 pc). This spread is partly due to the strong correlation with E_iso ($\mathit{R}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{\propto }{\mathit{E}}_{\mathrm{iso}}^{\mathrm{1}\mathit{/}\mathrm{2}}$) and is thus related to the spread in the energy released during the prompt emission in gamma-rays (width  ≈ 4.3 dex). On the other hand, since R(t) scales with ${\mathit{A}}_{\mathit{\star }}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{2}}$, A _⋆  would have to vary by a factor of 100 in the right way to reduce the width of this distribution by only a factor of 10.

Fig. 1 Shown here are lower (to-the-right pointing triangles) and upper (to-the-left pointing triangles) limits on the position of the wind termination shock based on Eq. (1), assuming A ⋆  = 1 in all cases (Table 3). Note that GRB 060614 is a much debated burst (see Sect. 3.1.1). The step curves are the cumulative distributions to the lower limits up to which a wind profile is identified and upper limits after which a constant density medium (ISM) agrees with the data.

Is the width of the distribution we have found for the upper limits on the wind termination shock reasonable? Based on numerical wind models of Wolf-Rayet stars, Fryer et al. (2006) found that for A _⋆  = 1 the radius of the termination shock is approximately given by ${\mathit{R}}_{\mathrm{T}}\mathrm{=}\mathrm{(}{\mathit{n}}_{\mathrm{ISM}}\mathit{/}\mathrm{100}\mathrm{cm}-3{\mathrm{)}}^{\mathrm{-}\mathrm{1}\mathit{/}\mathrm{2}}\mathrm{pc}\mathit{.}$(2)ISM densities of the order of 10⁶ cm^-3 are then required to reduce R_T to 0.01 pc. These large-scale gas densities are rather unique and only typical for dense cores of molecular clouds (with the Rho Ophiuchi Cloud as an example, e.g., Klose 1986).

Observationally, gas densities could only be derived for a few bursts because of the lack of radio data. So far, the highest values are about 600 cm^-3 (Frail et al. 2006; Thöne et al. 2010), while the nominal value is about a factor of 100 smaller (Frail et al. 2006). These measurements do not necessarily rule out the model of Fryer et al. (2006). Successful radio observations might have picked out a certain class of GRBs. Furthermore, radio observations are very challenging at early times because the brightness of the afterglow is increasing in the radio bands while the afterglow is already decaying in the optical and X-ray bands (Zhang & Mészáros 2004). Thus, it is difficult to extract information on the direct vicinity of the progenitor. On the other hand, if the average particle density of the circumburst medium is about 1−10 cm^-3, then additional mechanisms are required to bring the wind termination-shock radius closer to the star (e.g., van Marle et al. 2006).

The other possibility to reduce R_T is to decrease A _⋆ , since ${\mathit{R}}_{\mathrm{T}}\mathrm{\propto }{\mathit{A}}_{\mathit{\star }}^{\mathrm{1}\mathit{/}\mathrm{2}}$ (Chevalier et al. 2004). For example, for A _⋆  = 0.01 all data points in Fig. 1 would shift along the R-coordinate to higher values by a factor of 10 (Eq. (1)), while R_T would decrease by the same factor. In this case, lower circumburst gas densities would be required. This touches upon the question on how small A _⋆  can be. Studies of WR stars do not favour values of less than 0.01 in polar directions (Eldridge 2007, his Table 1). Moreover, observations of nearby WR stars do not show evidence for these low values either (Nugis & Lamers 2000). However, the majority of nearby WR stars are surely not seen pole-on, in contrast to GRBs. Therefore, it is difficult to decide if these observational constraints on WR stars can be applied to GRB progenitors.

On the other hand, the six GRBs that favour a free wind medium (Table B.3) have large lower limits on the wind termination-shock radius (Table 3). This matches theoretical models by van Marle et al. (2007, 2008), which allow R_T to extend up to several parsecs.

The separation between the lower and upper limits for wind and ISM-profiles, respectively, on the wind termination-shock radius could be even larger. Refining the lower and upper limits is difficult, however. The lower limits, which are deduced from $\begin{array}{}{\mathrm{˜}}_{}\\ {\mathit{t}}_{\mathrm{end}}\end{array}$ (Fig. C.1, Table 3), depend on the observing strategies due to the brightness of the afterglow and the brightness of the underlying host galaxy in the optical bands. On the other hand, the upper limits, which are deduced from $\begin{array}{}{\mathrm{˜}}_{}\\ {\mathit{t}}_{\mathrm{start}}\end{array}$ (Fig. C.1, Table 3), can be affected by additional radiation components at early times.

4. Summary and conclusion

After applying several selection criteria (closure relations, the differences in the spectral and temporal slopes, and the flux density ratio, F_opt/F_x) we selected the best-sampled Swift GRBs with well-observed optical as well as X-ray afterglow data from June 2005 to September 2009. Altogether 27 bursts entered our sample, which was used to investigate the density profile of the circumburst medium (constant density medium or free wind), including one short burst (GRB 051221A) and one controversial event in terms of classification (GRB 060614), which successfully passed our selection criteria among all bursts. The other 25 events are classified as long bursts without doubt.

Combining optical with X-ray data is advantageous because optical data usually allow for a more precise determination of the temporal decay slope of an afterglow, while X-ray data can in general be used to extract the SED. Combining both emission components substantially improves our capability to distinguish between an ISM and a wind medium. Thereby, we concentrated on the late-time evolution, i.e., times when the proper afterglow is not affected anymore by flares and additional radiation components (e.g., the reverse shock, central engine activity).

Our study shows that only six of the 25 long bursts (24%) investigated here (GRBs 050603, 070411, 080319B, 080514B, 080916C and 090323) showed evidence for a free wind medium at late times. In the other cases (76%), except for the short burst GRB 051221A and 060904B, the blastwaves were propagating into a constant density-medium. In particular, the controversial burst GRB 060614 favours an ISM profile. This is not in disagreement with a massive-star origin as our result for long bursts indicates.

In addition, we were able to set limits on the wind termination-shock radii of the corresponding GRB progenitors. Only 24 of 27 bursts (Table 3) had good enough data to perform this analysis. Fixing the relative mass-loss rate to A _⋆  = 1, the distribution we deduced covers three orders of magnitude. We find a tentative grouping into (long) GRB progenitors with comparably small and comparably large termination-shock radii. Whether this points to two distinct populations of (long) GRB progenitors or if this is a selection effect, remains an open issue. At least theoretically it is well possible that the long burst population splits into single star progenitors and those belonging to a binary system (Georgy et al. 2009). Further observational data are required to reveal a potential binary nature of the long burst progenitors.

Online material

Appendix A: Smoothly broken power law of the order m

The equation for a smoothly broken power-law of the order m, F_ν,m(t), was derived by recursion in the following way. Let us assume the function F_ν,m(t) consists of m power-law segments connected by (m − 1) breaks. To add an additional power-law segment $\begin{array}{}{}_{\mathrm{˜}}\\ {\mathit{F}}_{\mathit{\nu ,m}\mathrm{+}\mathrm{1}}\end{array}\mathrm{\left(}\mathit{t}\mathrm{\right)}$, we first normalised the new power-law segment to the previous one, F_ν,m(t), at the break time t_b,m: $\begin{array}{}\mathrm{˜}\\ {\mathit{F}}_{\mathit{\nu ,m}\mathrm{+}\mathrm{1}}\end{array}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}{\mathit{F}}_{\mathit{\nu ,m}}\mathrm{\left(}{\mathit{t}}_{\mathit{b,m}}\mathrm{\right)}{\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,m}}}\right)}^{\mathrm{-}{\mathit{\alpha }}_{\mathit{m}\mathrm{+}\mathrm{1}}}\mathrm{·}$(A.1)Here α_m + 1 is the slope of segment (m + 1). Second, we followed Beuermann et al. (1999) and introduced a smoothness parameter n_m so that the smoothly broken power law of order the (m + 1) takes the form ${\mathit{F}}_{\mathit{\nu ,m}\mathrm{+}\mathrm{1}}\mathrm{\left(}\mathit{t}\mathrm{\right)}\mathrm{=}{\mathrm{\left(}{\mathit{F}}_{\mathit{\nu ,m}}^{\mathrm{-}{\mathit{n}}_{\mathit{m}}}\mathrm{\right(}\mathit{t}\mathrm{)}\mathrm{+}\begin{array}{}\mathrm{˜}\\ {\mathit{F}}_{\mathit{\nu ,m}\mathrm{+}\mathrm{1}}^{\mathrm{-}{\mathit{n}}_{\mathit{m}}}\end{array}\mathrm{(}\mathit{t}{\mathrm{\left)}}^{\mathrm{\right)}}}^{\mathrm{-}\mathrm{1}\mathit{/}{\mathit{n}}_{\mathit{m}}}\mathit{.}$(A.2)If the light-curve consists of m segments, both steps (adding and smoothing) have to be performed (m − 1)-times.

For example, let us derive the equation for a smoothly broken power-law (Beuermann et al. 1999). In this case m = 2, thus the function consists of two power-law segments connected by one break at the time t_b,1. The initial function is a simple power law F_ν,1(t) = C   t^−α₁. First, the second power-law segment, $\begin{array}{}{\mathrm{˜}}_{}\\ {\mathit{F}}_{\mathit{\nu ,}\mathrm{2}}\end{array}$, has to be connected to the first one at the time t_b,1 (step A.1) $\begin{array}{ccc}\begin{array}{}\mathrm{˜}\\ {\mathit{F}}_{\mathit{\nu ,}\mathrm{2}}\end{array}& \mathrm{=}& {\mathit{F}}_{\mathit{\nu ,}\mathrm{1}}\left({\mathit{t}}_{\mathit{b,}\mathrm{1}}\right){\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,}\mathrm{1}}}\right)}^{\mathrm{-}{\mathit{\alpha }}_{\mathrm{2}}}\\ & \mathrm{=}& \mathit{C}{\mathit{t}}_{\mathit{b,}\mathrm{1}}^{\mathrm{-}{\mathit{\alpha }}_{\mathrm{1}}}{\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,}\mathrm{1}}}\right)}^{\mathrm{-}{\mathit{\alpha }}_{\mathrm{2}}}\mathrm{·}\end{array}$Second, the transition has to be smoothed by weighting both functions at the point of intersection (step A.2) $\begin{array}{ccc}{\mathit{F}}_{\mathit{\nu ,}\mathrm{2}}& \mathrm{=}& {\mathrm{(}{\mathit{F}}_{\mathit{\nu ,}\mathrm{1}}^{\mathrm{-}{\mathit{n}}_{\mathrm{1}}}\mathrm{(}{\mathit{t}}^{\mathrm{)}}\mathrm{+}\begin{array}{}\mathrm{˜}\\ {\mathit{F}}_{\mathit{\nu ,}\mathrm{2}}^{\mathrm{-}{\mathit{n}}_{\mathrm{1}}}\end{array}{\mathrm{(}{\mathit{t}}^{\mathrm{)}}}^{\mathrm{)}}}^{\mathrm{-}\mathrm{1}\mathit{/}{\mathit{n}}_{\mathrm{1}}}\\ & \mathrm{=}& {\left({\mathit{C}}^{\mathrm{-}{\mathit{n}}_{\mathrm{1}}}{\mathit{t}}^{{\mathit{\alpha }}_{\mathrm{1}}\hspace{0.17em}{\mathit{n}}_{\mathrm{1}}}\mathrm{+}{\mathit{C}}^{\mathrm{-}{\mathit{n}}_{\mathrm{1}}}{\mathit{t}}_{\mathit{b,}\mathrm{1}}^{{\mathit{\alpha }}_{\mathrm{1}}\hspace{0.17em}{\mathit{n}}_{\mathrm{1}}}{\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,}\mathrm{1}}}\right)}^{{\mathit{\alpha }}_{\mathrm{2}}\hspace{0.17em}{\mathit{n}}_{\mathrm{1}}}\right)}^{\mathrm{-}\mathrm{1}\mathit{/}{\mathit{n}}_{\mathrm{1}}}\\ & \mathrm{=}& \mathit{C}\hspace{0.17em}{\mathit{t}}_{\mathit{b,}\mathrm{1}}^{\mathrm{-}{\mathit{\alpha }}_{\mathrm{1}}}{\left({\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,}\mathrm{1}}}\right)}^{{\mathit{\alpha }}_{\mathrm{1}}\hspace{0.17em}{\mathit{n}}_{\mathrm{1}}}\mathrm{+}{\left(\frac{\mathit{t}}{{\mathit{t}}_{\mathit{b,}\mathrm{1}}}\right)}^{{\mathit{\alpha }}_{\mathrm{2}}\hspace{0.17em}{\mathit{n}}_{\mathrm{1}}}\right)}^{\mathrm{-}\mathrm{1}\mathit{/}{\mathit{n}}_{\mathrm{1}}}\mathrm{·}\end{array}$This leads to the equation found by Beuermann et al. (1999) for a smoothly broken power-law. Repeating both steps leads to a smoothly broken power-law of the order 3 (double smoothly broken power-law; Liang et al. 2008). Thus, looping (m − 1)-times over both steps results in a smoothly broken power law of the order m.

Appendix B: Tables

Fig. C.1 Optical and X-ray afterglow light curves of the 27 bursts that entered our sample. Upper panel: the optical data in the Rc band are shown as dots and the X-ray data at 1.73 keV as bigger dots with an error bar in time. The light-curve fits are over-plotted. Upper limits are shown as downwards-pointing triangles. The grey box is the overlapping time interval of the late-time evolution. Vertical dotted and dashed lines indicate breaks in the optical and X-ray band. Information on the SEDs are shown in the bottom left (see also Table B.1). The given extinction, ${\mathit{A}}_{\mathrm{opt}}^{\mathrm{host}}$, is the observed host-extinction in the Rc band based on the deduced host extinction in the V-band, ${\mathit{A}}_{\mathit{V}}^{\mathrm{host}}$. Additionally, we deduced the electron index, p, from βx. The electron index is either p = 2β if νc < νx or p = 2β + 1 if νc > νx (e.g., Zhang & Mészáros 2004). Its error was computed by propagating the uncertainty in βx. Middle panel: the flux density ratio between the optical and X-ray afterglow is shown as a solid line and its error as a dashed line for the shared time interval of the late-time evolution. The grey box represents the allowed parameter space of the flux density ratio (Table 1). The upper boundary is the expected flux density ratio for νc ≤ νopt, while the lower one shows the expected ratio for νc ≥ νx. If the cooling break is in between the optical and the X-ray bands, the expected flux-density ratio lies be in between these boundaries. The expected flux density ratio depends on the electron index. Not all bursts could be corrected for host extinction. The error on the electron index was neither propagated into the error of the expected nor of the observed flux-density ratio. Lower panel: the first logarithmic derivative of the flux-density ratio, ${}^{\mathrm{(}}\mathit{\alpha }_{\mathrm{x}}\mathrm{-}{\mathit{\alpha }}_{\mathrm{opt}}{}^{\mathrm{)}}\mathrm{\left(}\mathit{t}\mathrm{\right)}$, is shown as a solid curve and its error is plotted as a dashed line. For t/tbreak ≇ 1, the first logarithmic derivative is identical to the difference in the decay slopes obtained from the light-curve fit (asymptotic values). Usually breaks in the light curves tend to be smooth instead of sharp. Because of this, the first logarithmic derivative deviates from the asymptotic value close to a break depending on the smoothness of the break. Two solid lines are plotted to highlight the time interval when the asymptotic decay slopes were reached within 1σ. The precise values are shown on the left and in Table 3. Within 3σ, the asymptotic difference in the decay slopes agrees either with  +1/4, 0, −1/4 depending on the spectral and dynamical regime and the circumburst density profile. Furthermore, an envelope is drawn around expected values,  +1/4, 0, −1/4, with a width of 0.1 to guide the eye.

Fig. C.1 continued.

Fig. C.1 continued.

Fig. C.1 continued.

Fig. C.1 continued.

Fig. C.1 continued.

Fig. C.1 continued.

Acknowledgments

We thank the referee for a very careful reading of the manuscript and a rapid reply. S.S. acknowledges support by a Grant of Excellence from the Icelandic Research Fund and Thüringer Landessternwarte Tautenburg, Germany, where part of this study was performed. D.A.K. acknowledges support from grant DFG Kl 766/16-1. A.R. acknowledges support from the BLANCEFLOR Boncompagni-Ludovisi, née Bildt foundation. T.K. acknowledges support by the DFG cluster of excellence “Origin and Structure of the Universe”. S.S. acknowledges Robert Chapman (U Iceland), Elisabetta Maiorano (CNR Bologna), Andrea Mehner (U Minnesota), Kim Page (U Leicester), Eliana Palazzi (CNR Bologna) and Gunnar Stefansson (U Iceland) for helpful discussions. This work made use of data supplied by the UK Swift Science Data Centre at the University of Leicester.

References

All Tables

All Figures

Fig. 1 Shown here are lower (to-the-right pointing triangles) and upper (to-the-left pointing triangles) limits on the position of the wind termination shock based on Eq. (1), assuming A ⋆  = 1 in all cases (Table 3). Note that GRB 060614 is a much debated burst (see Sect. 3.1.1). The step curves are the cumulative distributions to the lower limits up to which a wind profile is identified and upper limits after which a constant density medium (ISM) agrees with the data.

In the text

Fig. C.1 Optical and X-ray afterglow light curves of the 27 bursts that entered our sample. Upper panel: the optical data in the Rc band are shown as dots and the X-ray data at 1.73 keV as bigger dots with an error bar in time. The light-curve fits are over-plotted. Upper limits are shown as downwards-pointing triangles. The grey box is the overlapping time interval of the late-time evolution. Vertical dotted and dashed lines indicate breaks in the optical and X-ray band. Information on the SEDs are shown in the bottom left (see also Table B.1). The given extinction, ${\mathit{A}}_{\mathrm{opt}}^{\mathrm{host}}$, is the observed host-extinction in the Rc band based on the deduced host extinction in the V-band, ${\mathit{A}}_{\mathit{V}}^{\mathrm{host}}$. Additionally, we deduced the electron index, p, from βx. The electron index is either p = 2β if νc < νx or p = 2β + 1 if νc > νx (e.g., Zhang & Mészáros 2004). Its error was computed by propagating the uncertainty in βx. Middle panel: the flux density ratio between the optical and X-ray afterglow is shown as a solid line and its error as a dashed line for the shared time interval of the late-time evolution. The grey box represents the allowed parameter space of the flux density ratio (Table 1). The upper boundary is the expected flux density ratio for νc ≤ νopt, while the lower one shows the expected ratio for νc ≥ νx. If the cooling break is in between the optical and the X-ray bands, the expected flux-density ratio lies be in between these boundaries. The expected flux density ratio depends on the electron index. Not all bursts could be corrected for host extinction. The error on the electron index was neither propagated into the error of the expected nor of the observed flux-density ratio. Lower panel: the first logarithmic derivative of the flux-density ratio, ${}^{\mathrm{(}}\mathit{\alpha }_{\mathrm{x}}\mathrm{-}{\mathit{\alpha }}_{\mathrm{opt}}{}^{\mathrm{)}}\mathrm{\left(}\mathit{t}\mathrm{\right)}$, is shown as a solid curve and its error is plotted as a dashed line. For t/tbreak ≇ 1, the first logarithmic derivative is identical to the difference in the decay slopes obtained from the light-curve fit (asymptotic values). Usually breaks in the light curves tend to be smooth instead of sharp. Because of this, the first logarithmic derivative deviates from the asymptotic value close to a break depending on the smoothness of the break. Two solid lines are plotted to highlight the time interval when the asymptotic decay slopes were reached within 1σ. The precise values are shown on the left and in Table 3. Within 3σ, the asymptotic difference in the decay slopes agrees either with  +1/4, 0, −1/4 depending on the spectral and dynamical regime and the circumburst density profile. Furthermore, an envelope is drawn around expected values,  +1/4, 0, −1/4, with a width of 0.1 to guide the eye.

In the text

	Fig. C.1 continued.
In the text

	Fig. C.1 continued.
In the text

	Fig. C.1 continued.
In the text

	Fig. C.1 continued.
In the text

	Fig. C.1 continued.
In the text

	Fig. C.1 continued.
In the text