Resonant and phonon-assisted ultrafast coherent control of a single hBN color center

Single-photon emitters in solid-state systems are important building blocks for scalable quantum technologies. Recently, quantum light emitters have been discovered in the wide-gap van der Waals insulator hBN. These color centers have attracted considerable attention due to their quantum performance at elevated temperatures and wide range of transition energies. Here, we demonstrate coherent state manipulation of a single hBN color center with ultrafast laser pulses and investigate in our joint experiment-theory study the coupling between the electronic system and phonons. We demonstrate that coherent control can not only be performed resonantly on the optical transition giving access to the decoherence but also phonon-assisted, which reveals the internal phonon quantum dynamics. In the case of optical phonons we measure their decoherence, stemming in part from their anharmonic decay. Dephasing induced by the creation of acoustic phonons manifests as a rapid decrease of the coherent control signal when traveling phonon wave packets are emitted. Furthermore, we demonstrate that the quantum superposition between a phonon-assisted process and the resonant excitation causes ultrafast oscillations of the coherent control signal. Our results pave the way for ultrafast phonon quantum state control on the nanoscale and open up a new promising perspective for hybrid quantum technologies.


I. INTRODUCTION
The ultrafast optical coherent manipulation of individual quantum systems is a key challenge in the development of quantum technologies [1]. Over the last decades several different solid-state systems, such as color centers in diamond [2], perovskite [3], colloidal [4], and self-assembled quantum dots [5], have been realized, improved, and utilized [6]. Recently, the family of quantum emitters has been joined by single-photon sources in layered van der Waals materials [7,8]. Compared to bulk materials these anisotropic crystals with their twodimensional character are particularly promising due to their structural flexibility [9]. Emitters have not only been found in semiconducting transition metal dichalcogenides [10][11][12][13][14] but also in insulating hexagonal boron nitride (hBN) in the form of color centers [15]. The latter have attracted massive attention because of their efficient single-photon emission at room temperature, emission over a wide spectral range spanning the nearinfrared [16], visible [17] and ultraviolet [18], and spectral tunability via external electric fields [19] or strain [20]. Furthermore, their characteristic coupling to phonons is remarkable. Due to the polar character, the color cen-ters couple to both, acoustic and optical phonon modes, giving rise to strong phonon sidebands at room temperature [21,22]. The bright emitters are not only found in deliberately grown samples [23] but also in commercially available nanocrystals, which makes their processability particularly easy [24]. Although some progress has been made regarding resonant optical excitation by demonstrating Rabi oscillations in g (2) -measurements [25] or the detection of lifetime-limited linewidths [26], the microscopic impact of phonons still requires deeper understanding. Studies indicate that there might be a connection between the phonon coupling and the emission linewidth of the emitter [27] but a consistent study in the spectral and the temporal domain, taking also phononassisted processes into account, has been pending. Here, we present the ultrafast coherent optical control of a single hBN color center. We demonstrate the role of different spectral jitter environments and reveal the impact of different phononic excitations. By combining experiment and theory, we explicitly quantify the emitter's pure and phonon-assisted coherence and at the same time access the impact of spectral noise. By performing a systematic study and thorough simulations, we reach a consistent picture of the emitter's coherence properties both in the spectral and the temporal domain and unravel the role of phonons in this context. To achieve an excellent agreement between the measurements and the simulation, we need to go beyond the calculation of correlation functions in the independent boson model [22] and also take into account the exact optical excitation and detection conditions, as well as phonon decay dynamics.

A. Optical emission
For our experiments, we use commercially available hBN nanopowder crystals, which host bright and stable single-photon emitters [22,24,28]. Details about sample preparation, the optical setup, and data acquisition are given in the Appendix and the Supplementary Material (SM), Sec. S1. Figure 1(a) displays the light emission from an individual color center. Single-photon emission is confirmed by a photon correlation measurement using a Hanbury Brown and Twiss setup presented in the SM, Sec. S1.C. We measure (dark lines) and calculate (bright lines) photoluminescence (PL, red) and PL excitation spectra (PLE, blue). The PL spectrum (dark red) is measured at a temperature of T = 80 K under pulsed excitation at 2.20 eV with a pulse duration of 230 fs. We find the characteristic sharp zero-phonon line (ZPL) at E ZPL = 2.036 eV, i.e., 160 meV below the excitation energy. The ZPL is accompanied by the typical peaked longitudinal optical (LO) phonon sidebands (PSBs) between 150 meV and 200 meV below the ZPL and a broad, continuous PSB stemming from longitudinal acoustic (LA) modes, exceeding E − E ZPL = −100 meV. We also detect multi-phonon sidebands for the sharp LO PSBs in the range of twice the LO phonon energies below the ZPL, i.e., between −300 meV and −400 meV. To measure the PLE spectrum, we scan the excitation laser energy above the ZPL and detect the PL integrated over ZPL and PSBs. The respective data in dark blue shows prominent maxima between E − E ZPL = 150 meV and 200 meV, which excellently agrees with the LO phonon energies. This demonstrates that the emitter is efficiently excited via a LO-phonon-assisted process [22,29,30]. Such a phonon-assisted excitation was for example combined with stimulated emission depletion experiments which, when compared to PLE measurements, revealed a red-shift of one LO mode in the excited state [31]. While at room temperature two prominent optical PSBs LO 1 at 165 meV and LO 2 at 200 meV are typically observed [21,22], at cryogenic temperatures the LO 1 emission exhibits two additional lines (see e.g. Fig. 6 in Ref. [22]), marked as LO 1 and LO 1 in Fig. 1(a). We assume that the resonances LO 1 and LO 2 stem from the phonon dispersion maxima at the Γ point, and the additional resonances LO 1 and LO 1 from Van Hove singularities at non-vanishing lattice momenta [32,33]. The 2nd LO-PSBs between E − E ZPL = −300 and −400 meV can be derived from these four lines LO 1 , LO 1 , LO 1 , and LO 2 by summing all ten possible combinations of them ( Fig. 1(a)).
The model, schematically displayed in Fig. 1(b), includes an optically driven two-level system (TLS) with the ground |G and the excited state |X , the spectrum of the laser pulse, the detection process, and the coupling to phonons giving rise to the characteristic sidebands (see Eqs. (5) and (9) for details). The coupling to phonons is considered via the independent boson model [34] taking into account a continuous distribution of LA phonons via a standard super-ohmic spectral density with Lorentz-Drude cutoff and four discrete LO modes via the Einstein model [? ]. Anharmonic decay of the LO modes into lower energy phonons and pure dephasing of the emitter, both leading to a finite width of the spectral lines, are included via Lindblad dissipators [35] (see Eq. (4) and SM, Sec. S2.A). The phonon parameters remain unchanged throughout this work and are summarized together with all other parameters used in the model in the SM, Sec. S3.E. The simulated spectrum of the 2nd LO-PSBs agrees well with the measured one, which indicates that all PSBs marked by green areas in Fig. 1(a) belong to the same ZPL. The narrow peaks at E − E ZPL ≈ −40 meV and −340 meV (marked by * ) do not belong to the investigated emitter (see SM, Sec. S1.E) and are therefore not considered in the following.
Overall, we find an excellent agreement between the measured spectroscopy data and the theoretical model considering one single bright state coupled to phonons ( Fig. 1(b)). All aspects of the measured time-integrated optical spectra of the emitter, namely all peak positions, their intensities and widths in both PL and PLE are reproduced by our model. The applicability of the model is further confirmed by measurements and simulations at a temperature of T = 8 K, which are shown in the SM, Sec. S3.C. In the following, we demonstrate that the ultrafast coherent control of the emitter is consistently reproduced and explained in the considered model.

B. Coherent control experiment
The time-integrated PL spectrum of a quantum emitter provides, however, only a very limited insight into its coherence properties and the impact of the specific environments such as spectral noise or lattice vibrations. To unravel these aspects, we perform ultrafast coherent control experiments on a single quantum emitter. We apply a double-pulse sequence, where the two laser pulses have equal intensities and a fixed phase relation, as schematically shown in Fig. 2(a) [36][37][38][39][40][41][42]. To gain a general understanding of the coherent control signal, it is instructive to first consider the TLS without any phonon coupling. The first laser pulse creates a coherent superposition |ψ = α G |G + α X |X of ground and excited state of the emitter, characterized by its occupation f = |X X| = |α X | 2 and coherence p = |G X| = |α * G α X |e −iφ0 , as depicted as blue and red curves in Fig. 2(a), respectively. The occupation typically decays within a few nanoseconds [22] and can therefore be assumed constant during the considered time scales in the picosecond range. The initial phase φ 0 of the quantum coherence φ is determined by the phase of the first laser pulse. During the following dynamics it rotates with the optical transition frequency of the emitter φ(t) = φ 0 + ω ZPL t, with ω ZPL = E ZPL / = 2π/T ZPL ≈ 3 fs −1 and its amplitude |p| is damped with the pure dephasing rate γ pd . After a delay time ∆t, the second laser pulse interacts with the quantum emitter. The final occupation of the quantum system depends on the relative phase between the second pulse and the quantum coherence. Because the applied pulse pair has a fixed phase relation, the phase difference between pulse and quantum state can be varied by tuning the time delay between the pulses ∆t. This is visualized in Fig. 2(a) by the two distinct timings for the second pulse (green Gaussians). When the second pulse arrives with a delay ∆t constr. such that the emitter's occupation gets efficiently enhanced ( Fig. 2(a), top), we call this "constructive". In contrast, as shown in the bottom of Fig. 2(a) for a delay of half a coherence period later ∆t constr. + T ZPL /2 = ∆t destr. the impact of the second pulse is "destructive" and the occupation is decreased. The occupation after the second pulse will therefore oscillate with the frequency ω ZPL as a function of the pulse delay ∆t. If the coherence has decayed entirely, the impact of the second pulse does not depend on the delay anymore, because the system is in a statistical mixture of |G and |X with |p| = 0.
After the second laser pulse, the occupation of the emitter decays fully into the ground state within its nanosecond radiative lifetime by emitting a photon. The next pulse pair hits the emitter at the repetition time of the laser (ns), which is long enough to guarantee the full decay of the emitter. Note, that a potential decay into a long-lived intermediate state with a lifetime longer than the repetition time does not disturb the coherent control experiment. Such a state would not be affected by the next pulse pair and therefore not contribute to the detected signal. The measured PL intensity is proportional to the occupation of the emitter, and we can compare the PL between constructive and destructive action of the second pulse within one period T ZPL . With this, we define the visibility of the coherent control experiment via which is a measure for the quantum coherence of the color center. In the experiment the visibility is retrieved from the data by a more involved procedure to reduce the impact of noise, as discussed in the SM, Sec. S1.D. Details on the simulation of the visibility are given in the SM, Sec. S2.D.3.
In the experiment, we use laser pulses (duration 1/α l = 230 fs) provided by a femtosecond Er:fiber laser system with nonlinear conversion stages (see Appendix). The visibility of the coherent control signal is detected by varying the delay on a femtosecond time scale (rotation of the phase), while the typical decay of the signal due to decoherence is on a picosecond time scale. Therefore, we use a Mach-Zehnder interferometer to create two pulses with a coarse delay varying on a picosecond time scale and a fine delay scanning over a few femtoseconds (see Fig. 2(a) and SM, Sec. S1.B). While changing the fine delay ∆t f , we detect the PL of the emitter via the PSBs and observe the oscillating signal depicted in Fig. 2(b). The dynamics cover several oscillations of the quantum coherence as marked by the relative phase between the emitter and the second pulse on the upper x-axis. For small coarse delay times (∆t = 1 ps), the oscillation amplitude is close to unity ( Fig. 2(b), left panel) due to the almost perfect constructive and destructive interaction for short delays, implying v PL ≈ 1 as discussed in the SM, Sec. S1.D. As mentioned earlier, for longer delays (∆t = 33 ps in the right panel) the decoherence of the system before the second pulse results in a reduced fine delay (fs) mean PL intens.  oscillation amplitude of the detected PL intensity and consequently a smaller coherent control visibility. Figure 2(c) displays measurements (dark dots) and corresponding simulations (bright line) of the coherent control PL visibility for delays between ∆t = 0 and 35 ps at a temperature of T = 8 K. The optical excitation is in resonance with the ZPL and results in a steady decay of the visibility over time. The first crucial observation is that this decay cannot be fitted by a single mono-exponential or Gaussian function (see SM, Sec. S3.D). Instead, the measured decoherence dynamics are excellently reproduced by the theory when considering a combination of an exponential (homogeneous/pure) and a Gaussian (inhomogeneous) dephasing contribution via (see Eq. (8) and discussion after Eq. (10))
As illustrated in Fig. 3, these two dephasing contributions can be identified with spectral jitter components (also often referred to as spectral noise or spectral diffusion) [27], which act on different characteristic time scales. The critical value for this time scale is the repetition of the two-pulse coherent control experiment, which is T rep = 25 ns in our case ( Fig. 3(b)). If the considered spectral jitter component is faster than the delay between the pulses, the transition energy of the emitter is likely to change between the two pulses in a single run of the experiment (Fig. 3(c)). This disturbs the quantum phase and results in an exponential loss of the coherence [43], which is usually referred to as homogeneous dephasing. If the jitter is slower than the repetition time of the experiment, it is more likely that the transition energy of the emitter changes between single runs of the experiment ( Fig. 3(d)) [44]. In this situation the successive repetition of the coherent control experiment can be treated as an ensemble average over the spectral distribution representing the jitter (see SM, Sec. S2.D.5 for details). In our case, the Gaussian dephasing dynamics correspond to a Gaussian spectral distribution of the jitter-induced energy shifts, which is usually called inhomogeneous broadening. The combined spectral jitter is schematically shown in Fig. 3(a). Figure 3(e) and Table I provide an overview of the various time scales including fast and slow spectral jitter. In addition, the table lists the determined dephasing times and their corresponding widths in optical spectra and in which figures the respective time scales are found. From homogeneous and inhomogeneous dephasing, we estimate an effective emitter coherence time, where Eq. (2) has dropped to 1/e of the initial amplitude, to This value is on the same order of magnitude as previously reported coherence times of around 80 ps in a hBN  emitter [45]. T eff ≈ 21 ps indicates that the investigated emitter experiences a strong spectral jitter on different time scales. The jitter might stem from laser induced charge fluctuations in hBN or the substrate [46], which influence the transition energy of the quantum emitter, e.g., via the Stark effect.

D. Telegraph noise
Interestingly, in the course of the experiments after cooling from room temperature to cryogenic temperatures the same emitter sometimes performed random jumps between two distinct emission energies, i.e., random telegraph noise (RTN), on a millisecond time scale as shown in Fig. 4(a) [47]. When assuming that the spectral jitter is induced by variations of the local surrounding of the emitter, heating and cooling cycles can likely result in varying noise situations. For each time step in this measurement the PL was integrated over 10 ms. As explained in the previous section (see also Table I) this spectral jitter is slow and can therefore be treated as a spectral ensemble. The RTN character is evident, because one of the emission lines is dim when the other is bright. We also find slight variations in the peak po-sitions hinting towards an additional fluctuation of the peak energies with a smaller amplitude on the order of the spectral linewidth. The amplitude of this jitter is in the range of a few tens of µeV (see SM, Sec. S1.F for details), which is in reasonable agreement with the value determined from the inhomogeneous decoherence σ j = 35 µeV in the coherent control experiment in the previous section. After time integration, the RTN results in the double-peak structure of the ZPL's PL spectrum in Fig. 4(b). The corresponding coherence dynamics in Fig. 4(c) exhibit a characteristic beat, given by the energy splitting of the RTN of ∆E ZPL = 0.7 meV. Also, this behavior is reproduced by our model when additionally considering the RTN as an ensemble average (see discussion after Eq. (10) in the Appendix). These results demonstrate that even in the presence of different sources of spectral jitter, i.e., the fast jitter leading to the homogeneous decoherence, the slow Gaussian noise resulting in the inhomogeneous decoherence, and the RTN introducing a beat, coherent control of the quantum state of an individual hBN color center is possible. The respective measurements of the coherence dynamics allow us to determine the present noise characteristics with high precision. I. List of the different time scales in the experiment. "Fast" and "slow" spectral jitter are defined with respect to the repetition time of the experiment. The decoherence times are given for a temperature of T = 8 K. The listed times are also illustrated in Fig. 3 Phonon-assisted excitation and emission have a remarkable strength in hBN color centers [21,22,29,48,49], with PSBs reaching peak intensities of almost 10% of the ZPL at T = 80 K, as discussed for the PL and PLE spectra in Fig. 1. We now detune the carrier frequency of the applied laser pulses to study the emitter's coherence when phonon-assisted excitations are addressed. As discussed throughout the SM and in the Appendix, the experimental results crucially depend on the spectral widths of the laser and the emitter. To reduce the impact of slow spectral jitters, we increase the temperature to T = 80 K, which results in a broadening of the spectral features, especially the ZPL. This generates a more robust situation as the slow jitters have a reduced relative impact on the broadening. Note, that these experiments were performed when the emitter did not show any RTN, i.e., only a single ZPL was present. The violet data in Fig. 5(a) (dark squares: experiment, bright line: simulation) shows the coherent control visibility for a resonant ZPL excitation which exhibits a dephasing time of 2/γ pd = 3.5 ps (see Eq. (8)). This more than ten-fold increase of the decoherence rate compared to T = 8 K confirms the temperature-dependent broadening of the ZPL (see Eq. (5) in Appendix). We now tune the pulsed optical excitation to the LO PSBs, LO 1 (filled circles, solid line) and LO 2 (open circles, dashed line) as depicted in the inset of Fig. 5(a). Both resonances of the PSBs in the PLE spectrum stem from flat parts of the phonon dispersion with a vanishing group velocity [32]. Consequently, the created phonons remain localized in the region of the color center's wave function. In Fig. 5(a) we plot the measured (dark blue dots) and simulated (bright blue  lines) visibility of the coherent control signals, i.e., the coherence dynamics, for optical excitations on the LO 1 and LO 2 PSB above the ZPL. As a reference, we additionally display the autocorrelation signal of the laser (green, dark diamonds: experiment, bright line: simulation, see SM, Sec. S3.B). It is measured by directly interfering the two laser pulses without interacting with the sample and calculated via Eq. (11). We find that for the LO-assisted excitations, the coherence decays significantly faster than for the resonant excitation of the ZPL (violet), but decays slower than the autocorrelation of the laser pulses.
As schematically shown in Fig. 5(b), the first laser pulse (blue arrow) creates a phonon-assisted coherence between the states |G, 0 and |X, 1 , where the second entry describes the number of excited phonons (orange ellipse). Any phonon-induced dephasing (PID) then leads to the destruction of the quantum coherence between the optical phonon, the ground state, and the excited state of the quantum emitter (dashed ellipse). This process is naturally faster than the decoherence without phonons. One obvious source of PID is the finite lifetime of optical phonons, which typically decay anharmonically into lower energy modes [50]. This process is given by the transition |X, 1 → |X, 0 . As shown in the SM, Sec. S2.E.3 a finite LO lifetime has a similar impact on the coherent control signal as a non-vanishing   dispersion of the optical modes, which is another possible source of decoherence. To reproduce the measured visibility dynamics in Fig. 5(a) we find phonon decay times of 1/γ LO1 = 0.3 ps and 1/γ LO2 = 0.4 ps within our model. LO phonon lifetimes in hBN can also be measured by Raman scattering and are on the order of a few ps [50], which is longer than the values found here. Therefore, additional decoherence mechanisms apart from their decay, e.g., a non-vanishing dispersion, are likely.

F. Acoustic phonons
Next, we reveal the impact of LA phonons on the coherence of the hBN color center. To confirm that LAphonon-assisted excitations are feasible [30,49], we perform a PLE scan at T = 80 K on a meV range around the ZPL as depicted in Fig. 6(a), where the experiment is shown as dark dots and the simulation as bright line.   We find that the PLE spectrum is asymmetric with a flatter slope on the low energy side and a plateau-like region on the high energy side. This asymmetry has two reasons: (i) at T = 80 K the LA phonon emission is more likely than LA phonon absorption, which results in a stronger LA PSB at positive detuning forming the plateau. To verify this, the dashed line in Fig. 6(a) depicts the calculated PLE spectrum without phonon coupling. (ii) The measured PLE spectrum is a convolution of the intrinsic absorption spectrum of the emitter and the laser spectrum. The asymmetry of the laser spectrum has a non-trivial impact on the final PLE spectrum, as seen in the dashed curve and discussed in detail in the SM, Sec. S4.B. The coherent control visibility with acoustic-phononassisted optical excitation is shown in Fig. 6(b) for T = 80 K. The experiment is depicted as dark dots, the simulation as bright lines, and the detuning between the laser energy and the transition energy δ = E laser − E ZPL increases from top to bottom, with negative detunings (red) in the left and positive (blue) in the right column. The resonant case (violet) is the same as in Fig. 5(a) and the considered detunings ( δ = ±4, ±5, and ±6 meV) are additionally marked as red (negative) and blue (positive) vertical lines in the PLE spectrum in Fig. 6(a). As a reference, we again plot the laser autocorrelation as green dashed lines in Fig. 6(b). For detuned (non-resonant) excitation, all coherent control signals exhibit a rapid initial visibility drop within 0.5 ps, which is followed by a slower decoherence or, interestingly, an intermediate recovery of the coherence. We find a remarkable agreement between the experimental data and the theoretical simulations. Notably, the asymmetric spectral shape of the femtosecond laser pulses has to be accurately considered in the theoretical model to achieve this almost perfect consistent match between the coherent control signal dynamics and the PL spectrum in Fig. 1(a) (see SM, Sec. S4.A for details). According to our model, the observed coherence dynamics for delays below ≈ 1 ps have two contributions (see Eq. (12) and SM, Sec. S2.D.4): (i) An initial loss of coherence resulting in the rapid drop of the coherent control visibility, followed by a quasi-stationary coherence for delays larger than ≈ 1 ps; (ii) a quantum beat resulting in an oscillation for delays below ≈ 1 ps.
Finding (i) is specific for the coupling to acoustic phonons and well known from self-assembled semiconductor quantum dots as a PID effect [51,52]. Its physical origin is the formation of a polaron, i.e., a deformed lattice configuration in the emitter region, which is created on the time scale of the optical excitation. As illustrated in Fig. 6(c), this rapid lattice deformation after the first laser pulse is accompanied by the emission of a phonon wave packet due to the linear dispersion of the LA phonons. As this phonon pulse is emitted from the color center with the speed of sound (several nm/ps in hBN), it carries part of the coherence with it, which is lost for the quantum state of the emitter. This loss of coherence is then probed by the second laser pulse in the coherent control experiment.
The amplitude of the visibility at ∆t = 2 ps reflects the coherence of the polaron state, which is generated by the first pulse. In each line of Fig. 6(b) positive (blue, right) and negative (red, left) detuning with the same absolute value are directly compared. From the visibility values at ∆t = 2 ps we find that the PID is always more pronounced for positive detunings caused by the asymmetric laser spectrum and a mismatch between phonon emission and absorption due to the thermal occupation of phonons n ph at T = 80 K (n ph ≈ 1 for modes with E ph ≈ 5 meV). For positive detuning, phonons must be generated to occupy the polaron, while for negative detunings phonons need to be absorbed. The two processes scale with (n ph + 1) ≈ 2 and n ph ≈ 1, respectively. This mismatch results in the more pronounced PID drop of the visibility for positive detunings. The quasi-stationary signal at ∆t ≈ 2 ps is a measure for the overlap of the laser spectrum with the ZPL. The visibility is generally smaller for larger absolute detuning, e.g., v PL (∆t = 2 ps) ≈ 0.5 at δ = +4 meV compared to ≈ 0.25 at δ = +6 meV, as the polaron is excited less directly and more via phononassisted processes. This leads to a stronger loss of coherence induced by the creation of the phonon wave packet when the absolute detuning |δ| is larger.
The additional oscillation of the signal, i.e., effect (ii), results in dynamics which are even faster than the laser autocorrelation (green curves in Fig. 6(b)). We can understand these oscillations when considering the theoretical expression for the coherent control signal for short delays in Eq. (12) in the Appendix. We find that the detuning of the laser and the remaining overlap with the ZPL effectively select two spectral components: the pure emitter coherence and the phonon-assisted coherence at the carrier frequency of the laser pulse. These two contributions result in a quantum beat with the frequency ω beat = |δ| that shows up as an oscillation in the coherent control signal. This directly explains why the oscillations get faster for larger detunings. The expected minimum from such a beat at T beat /2 = π/|δ| is marked by vertical dashed lines in Fig. 6(b) and matches well with the measured minima in the visibility dynamics. This beat is a clear indication of the coherent coupling between the acoustic phonon modes and the emitter.

III. CONCLUSIONS
By performing optical coherent control experiments with ultrafast laser pulses we have revealed the coherence dynamics of a single hBN color center. We found an excellent agreement with simulations, being consistent for spectrally resolved time-integrated PL and the time-resolved coherent control signal. We highlighted the impact of different sources of spectral jitter: Fast jitter on or below the picosecond time scale resulted in homogenous/exponential decoherence dynamics, while slow jitter on a nanosecond and longer time scale led to inhomogeneous/Gaussian decoherence dynamics. The presence of random telegraph noise was revealed by a characteristic beat of the coherent control signal. Furthermore, we demonstrated the impact of phonons on the emitter coherence by using phonon-assisted excitation. In the case of optical phonons we found that the increased decoherence rates are due to dephasing processes of the phonon state, which partially stem from its anharmonic decay. When the phonon-assisted excitation involved acoustic phonon modes, we found that (i) the generation of the phonon wave packet led to a characteristic drop of the coherent control signal and (ii) the simultaneous resonant and phonon-assisted excitation resulted in the selection of two energetically detuned contributions, which manifested in beat dynamics in the signal.
Hence, the coherence of hBN color centers is limited by phonon induced decoherence and spectral jitter. While the phonon impact could be controlled by deliberately shaped laser pulses [53-55], control over the crystal environment could reduce spectral jitter. An overall increased coherence time will enable advanced coherent control schemes, e.g., employing coherent optoelectronics [56].
Our demonstration of phonon-assisted coherent control of individual hBN color centers is a key step towards hybrid quantum technologies, which combine electronic and phononic excitations [35]. We envision that the impact of the phonons can further be controlled by external acoustic driving [57], which will allow to deterministically address the photon emission properties.

Experiment
Commercially available hBN multilayer nanopowder (Sigma-Aldrich, diameter < 150 nm) is annealed at 850 • C to increase the number of single-photon emitters and randomly scattered on a clean quartz substrate. The sample is then placed in a variable temperature helium bath cryostat with a low-temperature transmission microscope inside the sample chamber. The microscope consists of two objective lenses with a numerical aperture NA of 0.82, enabling efficient photon collection. The sample, as well as the microscope, is cooled by evaporating the coolant in the sample chamber. The cryostat is operated using liquid helium or nitrogen, resulting in sample temperatures of 8K and 80K respectively. Spectrally tunable laser pulses with durations of 230 fs are provided by a frequency-doubled supercontinuum from a femtosecond Er:fiber laser. To obtain two identical pulses, a single laser pulse is split at a 50:50 plate beamsplitter. The path of the split pulses can be controlled by a combination of a 204 mm long translation stage and a 100 µm closed-loop piezo scanner, resulting in a coarse and fine delay of the pulses in the range of ps and sub-fs respectively. The beam path of the two pulses is rejoined by an identical 50:50 plate beamsplitter and then focused on the single-photon emitter. The photoluminescence is detected either by a spectrometer and cryogenicallycooled CCD for spectroscopy and the dynamics measurements or a single-mode fiber-coupled Hanbury Brown and Twiss setup for photon correlation measurements. Further details can be found in the SM, Sec. S1.

Theory
In the SM in Sec. S2 we derive all required equations to model the depicted optical signals in detail. Here, we summarize the most important final results for each measurement. We model the emitter as a two-level system (TLS) and consider the coupling to phonons in the independent boson model and the coupling to light semiclassically leading to the Hamiltonian Here, ω TLS is the transition energy of the emitter in the absence of phonon coupling. The zero phonon line energy additionally includes the polaron shift ω ZPL = ω TLS − n |g n | 2 /ω n , depending on the phonon frequencies ω n and their coupling strengths g n . X ( †) and b ( †) n are the TLS and phonon annihilation (creation) operators, respectively, and E(t) describes the laser pulses, being the projection of the electric field of the laser onto the emitter's dipole moment. We also include H pd and H anh to describe the pure dephasing (pd) of the TLS and the anharmonic decay (anh) of the LO phonons, which are treated in terms of Lindblad dissipators in the model.
The PL spectrum of the emitter is calculated via where R is the spectrometer response and γ pd is the pure dephasing rate. For the PL presented in Fig. 1(a) the spectrometer response is given by a normal distribution with mean ω and standard deviation ω∆λ/λ(ω), where λ(ω) is the corresponding wavelength and ∆λ = 0.3 nm is the resolution of the spectrometer. We emphasize that the spectrometer response is crucial for accurately modeling the PL spectra as the resolution influences the relative height of ZPL and LA-PSB. For some phonon parameters, this implies that they cannot be accurately determined from the PL spectra, unless the spectrometer response is well-known, as discussed in the SM, Sec. S3.C. (5) is the phonon correlation function, which captures the information on the emitter-phonon coupling. It is given by Here, J LA (ω) = n∈{LA} |g n | 2 δ(ω n − ω) is the LA spectral density, where the sum runs over all LA modes. The LO spectral density is defined in an analogous way. n(ω) is the thermal occupation of the phonon mode with frequency ω. In the case of the LO modes we assume n(ω) = 0 since they have energies > 150 meV and we consider temperatures well below 300 K. γ(ω) is the decay rate of the LO mode with frequency ω. Equation (7) shows that the central quantity determining the phonon correlation function is the spectral density. All details on the structure of the couplings g n , i.e., its dependence on the quantum number n, is lost in this effective description and we can in principle choose a spectral density for our modeling that fits best with the experimental observation. For the LAs we choose a super-ohmic spectral density with a Lorentz-Drude cutoff. The LO spectral density is modeled as a sum of four δ-functions, i.e., it describes four dispersionless modes in the Einstein model. The PLE spectrum and the visibility dynamics of the coherent control experiment can be calculated from the same expression whereẼ(t) is the envelope of a single laser pulse and ω l its carrier frequency. In this paper we define the carrier frequency as the maximum of the detected laser spectrum, as this can be conveniently identified in experiment. The temporal shape of the laser pulse is derived from an asymmetric spectrum, modeled by a Pearson-4 distribution (see SM, Sec. S3.B for details). The PLE is then given by and the visibility for a fixed carrier frequency, depending on the delay ∆t between the pulses, is i.e., it is the absolute value of I 0 , normalized by the PLE signal, which satisfies v(∆t) ≤ 1 with equality for ∆t = 0. Slow spectral jitter, i.e., jumps of the ZPL frequency on a time scale longer than a single repetition of the experiment, can be incorporated via an ensemble average. This is done by replacing exp[−iω ZPL (∆t − τ )] under the integral in Eq. (8) by the weighted average over a probability distribution p(ω ZPL ). We show in the SM in Sec. S2.D.5 that Gaussian jitter on a spectral scale well below the width of the laser spectrum simply leads to a prefactor exp[−σ 2 j (∆t) 2 /2] in Eq. (8) with 1/σ j ≈ 19 ps throughout the paper. Telegraph noise jitter, as seen in Fig. 2, is modeled by two δ-functions in p(ω ZPL ).
As discussed in the SM in Sec. S2.D.7, the visibility of the direct laser interference, measured in the same coherent control setup, is given by the autocorrelation of the electric field of a single pulse. Here, we deal with linearly polarized light, such that the amplitude of the electric field is proportional to E(t). Therefore the visibility of the laser interference is modeled by normalizing In the SM in Sec. S2.D.4 we discuss an approximation for the visibility dynamics in Eq. (10) when exciting close to the ZPL and considering short time scales where the pure dephasing can be neglected. In this case we arrive at the expression where δ = ω l − ω ZPL is the detuning between the laser and the ZPL. The function a(∆t) is given by i.e., it is the envelope of the pulse autocorrelation, weighted by the LA-PSB contribution to the absorption spectrum at the position of the carrier frequency ω l , as where (7) is a measure of the phonon coupling strength. For strong phonon coupling or large temperatures, B 2 tends towards zero.
In the PL spectra, B 2 is also the prefactor of the pure ZPL contribution. For weak phonon coupling, B 2 tends towards unity, i.e., the spectrum consists only of a ZPL and contains no PSBs. I ideal laser (ω ZPL ) is the ideal (perfect resolution) laser spectrum at the position of the ZPL, as defined in the SM, Sec. S2.D.2. Thus, b measures the overlap of the laser pulse spectrum with the ZPL. a(∆t) measures the strength of the LA-PSB at the carrier frequency for ∆t = 0 and decays along with the pulse autocorrelation. Thus we expect dynamics that decay on the time scale of the laser autocorrelation towards a final value of b/(a(0)+b), which is a measure of how efficiently we excite the system resonantly on the ZPL compared to off-resonant excitation via PSBs. The visibility dynamics exhibits a beat with the detuning frequency |δ|, additionally to the decay, which leads to dynamics faster than the pulse autocorrelation in Fig. 6, as discussed in the main text. As a final note on the visibility dynamics on short time scales we want to emphasize that the derivation was made for the case γ pd = 0, i.e., no pure dephasing. This leads to a δ-shaped ZPL and a strong dependence of the parameter b on the ZPL frequency ω ZPL via the overlap with the laser spectrum. Therefore in this regime (γ pd → 0), spectral jitter can have a significant impact and makes the experimental identification of the discussed phenomenon challenging. Thus, rather counter-intuitive, it is more instructive to perform the coherent control experiment with detuned excitation at more elevated temperatures. In that case γ pd is larger, i.e., the system dephases more quickly, leading to a broader ZPL and therefore less sensitivity on spectral jitter.

DISCLOSURES
The authors declare no conflicts of interest.

DATA AVAILABILITY
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request. S1 Experiment S1.A Sample preparation S1.B Experimental setup S1.C Single-photon emission characteristics S1.D Data acquisition and analysis S1.E Ruling out additional zero-phonon lines S1.F Characterization of the Gaussian spectral jitter We investigate single-photon emitters in commercially available hBN nanopowder (Sigma-Aldrich) [1,2]. The particles are disk-shaped with average diameters of < 150 nm. The powder is annealed at 850 • C for 30 min in an Ar environment to increase the number of single-photon emitters [3]. Subsequently, the powder is scattered onto a 500 µm thick quartz substrate, on which the crystals stick via van der Waals forces. Figure S1 shows these randomly scattered hBN particles in a scanning electron microscope (SEM) (a), white-light (b), and photoluminescence (PL) micrograph (c) at different positions of a typical sample. The focused laser beam is used to individually address single particles or small clusters of particles.

B. Experimental setup
The sample is mounted in a variable temperature liquid helium bath optical cryostat. A low-temperature microscope inserted in the sample chamber of the bath cryostat with two high numerical aperture objective lenses (NA = 0.82) allows for high resolution optical measurements. The sample chamber is filled with helium (nitrogen) vapor, cooling the sample down to 8 K (80 K). A sketch of the experimental setup is shown in Fig. S2.
Ultrafast laser pulses are provided by a femtosecond Er:fiber laser system. The pulses are used to generate a supercontinuum in a highly nonlinear fiber and are subsequently frequency doubled using a periodically poled lithium niobate (PPLN) crystal to obtain ≈ 230 fs pulses that are tunable over the full visible spectrum with a spectral width of 5 meV [4]. Each pulse is split in a Mach-Zehnder interferometer at a 50:50 plate beamsplitter. The two pulses are delayed with respect to each other up to 1360 ps of coarse optical delay using a motorized translation stage (PI M-521.DD, 204 mm long) and up to 667 fs of fine optical delay with a piezo stage (PI P-621.1CD, 100 µm scan range) before being rejoined by an identical 50:50 plate beam splitter. The polarization of the combined laser pulses can be adjusted using a λ/2 waveplate to match the dipole orientation of the single-photon emitter, optimizing the excitation efficiency and therefore maximizing the PL intensity of the emitter. The PL of the emitter is collected and focused into a spectrometer and cryogenically cooled CCD (Princeton Instruments Acton SpectraPro SP-2500 and Pylon PYL100BRX). The resolution of the spectrometer is <0.1 nm when using 1200 lines/mm and 0.3 nm  with 300 lines/mm diffraction gratings for high resolution or large spectral range measurements, respectively. For photon correlation measurements, the photoluminescence is routed into a single-mode fiber-coupled Hanbury Brown and Twiss setup. Using two Excelitas SPCM-AQRH-14 single-photon counting modules and a PicoQuant TimeHarp 260 time-correlated single-photon counting board, we detect the arrival times of single photons in a start-stop-type experiment and measure the pulsed second-order autocorrelation function g 2 (τ ), described in the following chapter.

C. Single-photon emission characteristics
To test for single-photon emission from the emitter shown in the manuscript, we perform a Hanbury Brown and Twiss experiment with pulsed laser excitation (see Fig. S2) to obtain the second-order autocorrelation function as described in chapter S1 B. The measurement is carried out at T = 8 K and the results are shown in Fig. S3. The pulsed laser excitation is resonant to the ZPL of the emitter. All photons in the spectral range of −130 meV to −440 meV detuned from the ZPL impinge on the single-photon detectors. This includes the first and second LO phonon sideband (PSB). No background correction is carried out. A set of exponential functions is fitted to the data, where A and B are the amplitudes, t decay is the decay time, τ is the time delay, τ 0 the separation of the laser pulses. The relative amplitude at time delay τ = 0 is A/B ≈ 0.3 which is below the single-photon threshold of 0.5. The decay time of the peaks is t decay ≈ 4 ns, which is in agreement with the typical lifetime of a single-photon emitter in hBN [1].  Figure S3. Single-photon emission. Second-order autocorrelation function g 2 (τ ) of the emitter shown in the main manuscript, measured at T = 8 K. The emitter is excited resonantly. The peak height at τ = 0 time delay is reduced to 0.3 relative to the mean peak height. The data is not background corrected.

D. Data acquisition and analysis
To detect the coherent control signal of the hBN single-photon emitter, PL spectra are recorded for fine delays ∆t f of the two pulses between 0 and 10 fs, with a step size of 0.1 fs. This is done at each coarse delay step with, depending on the experiment, varying step sizes from 0.1 ps to 2 ps. In the following discussion, all intensities also depend on the coarse delay, which is kept implicit. For resonant excitation, the PL of the PSB is collected, as shown in Fig. S4(a). For excitation via the PSB, the PL of the zero-phonon line (ZPL) is collected. The PL intensity data i PL (∆t f , E) is then processed by integrating in the direction of the vertical black line shown in Fig. S4(a), i.e., over the energy E and calculating the relative deviation from the mean PL intensity to obtain the amplitude of the oscillation (Fig. S4(b)). With this definition, the coherent control visibility for a given coarse delay ∆t from Eq. (1) in the main text reads for a noise-free, perfectly sinusoidal signal.
To account for deviations from this idealized case, a fast Fourier transform (FFT) is carried out to extract the amplitude of the oscillation corresponding to the transition frequency of the emitter, which at T = 8 K is 492 THz (2.036 eV), as shown in Fig. S4(c). This procedure allows for an accurate analysis even when the oscillation amplitude becomes close to the noise level of the experiment after an almost full decoherence of the two-level system. This amplitude corresponds to the PL visibility of the coherent control experiment as described in the main text. This analysis is repeated for each coarse delay step. The resulting visibility data set is shown in Fig. S4(d). It corresponds to the envelope of the PL oscillation and is related to the coherence of the system as described in the main text and the Theory section. The PL spectrum of the PSBs of the single-photon emitter exhibits a complex structure with narrow lines (see Fig.  1 of the main text and S5, upper spectrum). To identify possible contributions from other emitters in the spectrum, we utilize the PL spectra acquired during the coherent control experiments. We measure the PL spectra after resonant two-pulse excitation at time delays ∆t of the laser pulses without temporal overlap (e.g. 10 ps). By subtracting the PL spectrum at a minimum of the oscillation in Fig. S4(b) from the spectrum at a maximum, we remove all incoherent background from the PL spectrum. Thereby, also contributions from other emitters are efficiently suppressed which are only non-resonantly excited and therefore do not show long coherence. Figure S5. Identification of other emitters in the PL spectrum. Comparison of the single-pulse PL spectrum of the PSBs (upper curves, also shown in Fig. 1 in the main text) and the difference spectrum extracted from the coherent control experiment at ∆t = 10 ps (lower curves). The data is shifted vertically for clarity. The emitter is excited resonantly. The difference spectrum is obtained by subtracting the minimum of an oscillation from the maximum. Only parts of the spectrum related to the two-level system remain after this subtraction. We can rule out, that the narrow lines LO1, LO 1 , LO 1 , and LO2 stem from another emitter. The sharp lines marked with a star are not present in the difference spectrum and therefore belong to other emitters.
The resulting spectrum is depicted in the lower part of Fig. S5. Using this analysis, we find that LO 1 , LO 1 , LO 1 , and LO 2 can be attributed to the two-level system (ZPL) under study, since they still exhibit a coherent control signal at 10 ps time delay. ZPLs from other emitters will not appear in the difference PL spectrum. In the single-pulse PL spectrum a sharp line in the second LO sideband at E − E ZPL ≈ −335 meV (marked by the star) is observed which is not reproduced by our simulation (see Fig. 1 of the main text). It is also not observed in the difference PL spectrum and therefore stems from a ZPL of another emitter. This is also the case for the smaller narrow line at -348 meV.

F. Characterization of the Gaussian spectral jitter
In the main text we identify a spectral jitter on the timescale of the repetition rate of the coherent control experiment, which has an amplitude smaller than the linewidth of the ZPL. To determine the amplitude of this jitter, we measure a series of PL spectra depicted in Fig. S6(a) where each time step is integrated over 10 ms. The time required to read out the camera after each integration is ≈ 90 ms, during which no data is acquired. As discussed in the main text and seen in the figure, the additional random telegraph noise results in a two-peak structure of the ZPL. We determine the energies of the peak maxima by fitting two Gaussians for each time step. From these 2000 energies for the peak positions we calculate the histogram depicted in Fig. S6(b) as red dots. We then determine the width of the resulting two maxima by fitting two Gaussians, which leads to the standard deviations 34 µeV for the peak at E − E ZPL ≈ +0.35 meV and 18 µeV for the peak at E − E ZPL ≈ −0.35 meV. E ZPL is in this case defined as the center between the two states. Note that each state in the time series is the ZPL at that time.
These two widths (18 µeV and 34 µeV) contain contributions from the spectral jitter and the uncertainty of the energy determination due to photon noise, which originates from the finite number of photons collected in each peak  σ noise can be estimated using the following formula, under the assumption that the measurement noise is dominated by shot noise [5] σ noise = 16 9 Here, N ≈ 700 is the number of photons detected per peak and frame, s ≈ 180 µeV is the standard deviation of the Gaussian peak fit of each time step and a = 97 µeV is the spectral width of a single pixel of the spectrometer-CCD combination. This results in a noise-induced uncertainty of the emitter energy of σ noise = 9 µeV. Using (S2), we can estimate the jitter amplitude σ jitter to 16 µeV and 33 µeV for the two peaks, respectively. Due to the finite temporal resolution of the series of spectra, the determined values are lower bounds for the spectral jitter. Nevertheless, they are in good agreement with the values obtained from the coherent control experiment shown in the main text.

S2. THEORY
A. Model of the color center As described in the main text, we model the hexagonal Boron Nitride (hBN) color center as a two-level system (TLS), coupled to phonon modes of the bulk crystal via the independent boson Hamiltonian [6][7][8][9][10]. The full Hamiltonian, that describes the color center driven by the electric field of a laser, is given by where H 0 describes the free dynamics in absence of laser pulses. H IB describes the coupling of the phonon modes with the excited state of the TLS and is modeled via Here, ω TLS is the transition energy of the TLS and X = |G X| is the annihilation operator of the excited state |X , where |G denotes the ground state. The frequencies of the phonons are given by ω n with associated bosonic annihilation operators b n and independent boson coupling parameters g n . The index n collectively denotes all quantum numbers of the phonon modes under consideration. We focus here on the longitudinal acoustic (LA) phonons and the longitudinal optical (LO) phonons coupling to the emitter. We assume that the LO modes are dispersion-less, such that we can collect all LO phonon modes with the same frequency in a single effective mode [11].
We introduce a finite lifetime of the LO modes due to anharmonic coupling to a bath of low energy phonons via H anh [12]. When speaking of 'phonons' in the following, we do not mean this bath to which the LO phonons couple. Furthermore, we account for pure dephasing of the color center due to the coupling to some unspecified bath via H pd . We assume that both of these Hamiltonians lead to phenomenological Lindblad dissipators. H I (t) describes the time-dependent interaction of the color center with the external laser field. It is in the dipole and rotating wave approximation given by with E(t) = −2 X| d · E(t) |G denoting the (effective) laser field, where d is the dipole of the TLS and E(t) is the electric field of the laser at the color center's position. Before the interaction of the color center with any laser pulse, the TLS is in its ground state and the phonons are in a thermal state ρ phon at a given temperature T , such that the initial total density matrix at t → −∞ is given by where ρ pd and ρ anh are the density matrices of the baths inducing pure dephasing and anharmonic LO decay, respectively. Since we are interested in temperatures below room temperature and all LO phonons have energies > 150 meV, they are assumed to be in their vacuum state, i.e., b n ρ phon = ρ phon b † n = 0 for n ∈ {LO}. Note, that in the present work we are interested in picosecond and sub-picosecond dynamics, such that we can neglect the decay of the excited state of the color center, which happens on a nanosecond-timescale [3].

B. Lindblad equation in the polaron frame
In the following we will discuss the effective non-unitary dynamics of the coupled TLS-phonon system in the polaron frame. One moves to the polaron frame via the unitary polaron transformation, which diagonalizes H IB and is given by [10] where are multi-mode displacement operators in phase space for the phonon modes from (S5) and D n (α) is the singlemode displacement operator of the n-th mode with coherent amplitude α [13]. With this unitary transformation, the independent boson Hamiltonian is transformed into with being the polaron shifted transition frequency of the color center. In the emission spectra discussed later, it is the frequency of the zero phonon line (ZPL) which is usually the dominant peak in the spectra. J(ω) is the so-called (phonon) spectral density, given by The coupling of the LO phonons to the bath of low energy phonons with density matrix ρ anh and the coupling of the TLS to the pure dephasing bath with density matrix ρ pd is taken to be such, that the Born-Markov approximations hold and the total density matrix approximately factorizes at all times, i.e., where ρ(t) is the reduced density matrix of the TLS-phonon system. Under these conditions and due to the fact that the polaron transformation only acts on the TLS and the phonons, the total density matrix in the polaron frame is given by Our goal is to derive an effective equation of motion for the reduced density matrix in the polaron frame ρ P . In the absence of laser pulses, the propagation of the reduced density matrix in the polaron frame from an initial time t 0 to a time t is given by Here, Tr pd,anh [. . . ], denotes the trace over the subspaces of the baths inducing pure dephasing and LO decay.
is the free time evolution operator and U 0,P (t, t 0 ) = T P U 0 (t, t 0 )T † P is its version in the polaron frame. Finally, we introduce the phenomenological Lindblad dissipators, that describe pure dephasing and LO decay, by stating that the relation between ρ P at different times in the absence of laser pulses is obtained by solving [12] where is the Lindblad superoperator acting on the polaron system and is the Lindblad superoperator acting on the phonon system. The phonon and polaron subsystems are thus not coupled by the time evolution induced by H 0,P = T P H 0 T † P . The Lindblad dissipators for pure dephasing of the TLS and decay of the LO phonons are given by respectively. Here, γ pd is the pure dephasing rate and γ (n) LO is the decay rate of the respective LO mode n. The subscript n ∈ {LO} under the sum indicates that the dissipator only acts on the LO modes.

C. Spectra
The central formula used in this section to calculate the spectra of the laser and the emitter is presented in Ref. [14] as the time-dependent physical spectrum of light. It is given by where r(τ, ω) is the response function of the spectrometer, i.e., the outgoing field, when the incoming field is a δ-pulse, and ω is the frequency, at which the spectrum is measured. E (±) is the positive/negative frequency part of the electric field operator at the position of the spectrometer, projected onto the detector's dipole.

Photoluminescence of the emitter
To compare between theory and experiment, we need to know the spectral densities of the LO and LA phonons. Information on these can be gained from photoluminescence (PL) spectra of the emitter. The time-dependent emission spectrum of the emitter is given by [12,14,15] We consider the case in which the emitter is excited by a laser pulse with a duration of a few hundred fs. The subsequent radiative decay happens on a nanosecond-timescale, whereas the phonons thermalize on a much faster timescale of a few picoseconds [16]. After the pulse interaction, only the part of the density matrix proportional to |X X| contributes to the correlation function X † (t 1 )X(t 2 ) . In this subspace the phonons thermalize into a displaced thermal phonon state B − ρ phon B + [12], such that the reduced density matrix |X X| B − ρ phon B + is stationary under the effective equations of motion in (S16). For the time-integrated PL spectrum, the initial non-equilibrium dynamics are irrelevant compared with the stationary dynamics due to the long radiative lifetime. Thus in calculating the PL spectrum we assume a stationary emitter with the total density matrix |X X| B − ρ phon B + ⊗ ρ anh ⊗ ρ pd = B − X † ρ tot (−∞)XB + and time evolution due to H 0 . The correlation function X † (t 1 )X(t 2 ) is calculated in Sec. S2 E and is equal tõ C(t 1 , t 2 ) =C(t 1 − t 2 ) from (S73) as shown later. For a stationary emitter, i.e., long after the interaction with any laser pulses, the time-integrated spectrum is given by

S10
where we omitted all irrelevant prefactors. This can be rewritten as In the following section ((S26)) we will see that R(ω, Ω), the spectral response of the spectrometer, is given by δ(ω −Ω) in the case of perfect resolution, such that is the ideal emission spectrum of the emitter.

Ideal and detected laser spectrum
To calculate the laser spectrum, we start again with (S21) and use the fact, that for coherent laser light it is Employing Parseval's theorem, we can determine the detected, time-integrated laser spectrum We see, that the detected spectrum is given by the absolute square of the Fourier transform of the electric field, i.e., the ideal laser spectrum, if R(ω, Ω) = δ(ω − Ω), which suggests to call this the case of perfect resolution. The spectral response R(ω, Ω) can be determined by detecting the spectrum of a spectrally narrow light source approaching the limit of having a single frequency component ∼ exp(−iω l t), such that the detected spectrum is proportional to R(ω, ω l ). If we assume a homogeneous response over the frequency interval of interest, i.e., R(ω, Ω) = R(ω − Ω), we can thus determine the spectral response entirely. The ideal laser spectrum has a direct connection to the autocorrelation of the electric field of the laser. This can be seen by calculating the Fourier transform of the ideal spectrum where the -operator denotes the autocorrelation of the electric field E (+) (τ 1 ) .

D. Modeling the coherent control experiment
In the coherent control experiment, two laser pulses emitted from the same source and with a relative delay ∆t interact with the color center. The first pulse creates a coherence, which is then probed by the second pulse. The PL after the second pulse is the observable in the experiment. An optical excitation of the coupled phonon-TLS systems induces non-equilibrium dynamics that last for several picoseconds. In the following, the entire system is in a quasi-equilibrium, where excited and ground state subspaces are in their respective thermal equilibrium and the excited state decays on a nanosecond time scale. During the quasi-equilibrium stage, the time-dependent PL spectrum, integrated over an arbitrary frequency interval, yields the occupation of the excited state [12]. Since the decay dynamics takes much longer than the initial non-equilibrium dynamics, the time-integrated PL is basically proportional to the occupation of the emitter after the laser pulse interaction. For the simulation of the coherent control signal (Ramsey fringes), this implies that we need to calculate the occupation of the excited state |X after the interaction with the two laser pulses. In our case, since we neglect the decay of the emitter during the short non-equilibrium stage of the dynamics, this means that we can just calculate the occupation of the excited state at time t = ∞. Consequently, the observable is the TLS's excited state occupation where U (t, t 0 ) is the time evolution operator associated with the full Hamiltonian H(t) from (S4), i.e., andT is the time-ordering operator. The effective laser field of the double pulse reads and introduces a dependence of the final state occupation f on the delay ∆t between the pulses via the interaction H I (t) from (S6). The carrier frequency of the laser is given by ω l andẼ(t) is the potentially still complex envelope of a single pulse.

Perturbation expansion
In the following we perform a perturbation expansion of our observable in (S29) in terms of the exciting laser field. We will keep the calculation to the lowest non-vanishing order, which is the main contribution at small pulse S12 areas, i.e., low intensities of the exciting laser pulses. As it turns out, this contribution is sufficient to reproduce the experimental data. We start by writing down a perturbation expansion of the full time evolution operator in terms of the interaction Since U 0 does not change the occupation of the emitter, i.e., G| U 0 |X = 0, and ρ tot (−∞) ∼ |G G| describes the emitter in its ground state, the lowest order in (S32) is irrelevant and the first contribution to the observable is given by Again, since U 0 does not alter the occupation of the emitter, this reduces to where the correlation function is calculated in Sec. S2 E (see (S72)). f (2) denotes the second order perturbation theory contribution to the full observable f .

Representation of the observable in terms of the spectra of the emitter and the laser pulse
In this section we show that the second-order contribution f (2) is essentially determined by the spectrum of the applied laser pulse and by the absorption spectrum of the emitter. From (S34) we see that the second-order contribution to the full coherent control signal can be written as

S13
Introducing the variables τ = τ 2 − τ 1 and T = τ 2 /2 + τ 1 /2, such that dτ 1 dτ 2 = dT dτ , we obtain Now we introduce the representation of a δ-function By looking at (S27), we see that this expression contains the ideal spectrum of the double pulse used in the coherent control experiment. Note that for linearly polarized light, as is used in the experiment, we can assume E(t) ∼ E (+) (t) .
To make the dependence on the delay ∆t obvious, we will now express the double pulse spectrum via the single pulse spectrum. Using the representation of the laser field in terms of the two laser pulses from (S31), we find .

(S39)
The expression inside the absolute value is proportional to the ideal (R = δ) spectrum of a single laser pulse, as can be seen by examining (S27) together with (S31), and is therefore named I ideal laser (Ω). With this intermediate result, we can write Now we turn to the integral with the correlation function C(τ ). From Eqs. (S109), (S110), and (S25) we know, that is the ideal emission spectrum of the emitter. Therefore, is the ideal emission spectrum, mirrored at the ZPL, or in other terms: the absorption spectrum [1,6]. Consequently, the second order contribution to the observable for the coherent control experiment can be written as and depends on the ideal spectrum of a single laser pulse and on the absorption spectrum of the color center. It is not sensitive to chirp, which would be a frequency-dependent phase for the Fourier transform of E(t) and thus cancels in (S27). At higher orders of the expansion we expect chirp to have an impact. Since the laser spectrum is centered around the carrier frequency ω l , it is obvious from (S43) that the observable will be oscillatory with respect to the delay ∆t with a frequency determined by the weight of the overlap of the ideal laser spectrum and the absorption spectrum. Due to the similarity with Ramsey interferometry, these oscillations are often called Ramsey fringes [17].

Relation between visibility of the Ramsey fringes and the phonon correlation function
One interesting quantity to obtain from the coherent control experiment is the phonon correlation function G −+ (τ ) (see (S106) in Sec. S2 E 4), which contains information on the emitter-phonon coupling. We want to gain information on this quantity from our observable f (∆t). To this aim, let us look at the contribution of lowest order in the electric field of the laser pulse from (S43). The spectrum of the laser I ideal laser is peaked around the carrier frequency ω l . The absorption spectrum I abs is peaked at the ZPL frequency ω ZPL , which roughly corresponds to 2eV but also contains dominant peaks at the LO phonon sidebands (PSBs). The time-dependence of f (2) (∆t) is dominated by a high-frequency oscillation, as can be seen by shifting the frequency integration by the mean frequency of the laser and absorption spectrum overlap such that f (2) (∆t) = 1 2π 2 R dΩ cos 2 (Ω + Ω)∆t 2 I ideal laser (Ω + Ω)I abs (Ω + Ω)

(S45)
The quantity I 0 (Ω) is a distribution with mean frequency 0 and the width ∆Ω of this distribution is given by the width of the overlap between laser spectrum and absorption spectrum, which is in our case on the order of a few meV. The mean frequencyΩ however, is on the order of 2 eV/ , as it is mainly determined by the frequencies ω ZPL and ω l . If the laser excites resonantly on the ZPL, we have Ω ≈ ω ZPL . If the laser excites on a LO-PSB for a LO phonon with frequency ω LO , we have Ω ≈ ω ZPL + ω LO . Therefore f (2) (∆t) shows an oscillatory behavior with the mean frequency Ω. The envelope of this oscillation varies on a timescale which is long compared to the oscillation period 2π/Ω, as it is determined by the comparably small width of the distribution I 0 in the meV range. As we will show in the following the visibility of the oscillations has a direct relation with the correlation function G −+ (τ ). To demonstrate this we rewrite f (2) , using cos 2 (x) = 1 2 [1 + cos(2x)], as . (S46) The first term does not depend on the delay ∆t, being the average of the oscillations, whereas the second term consists of a fast oscillatory contribution exp iΩ∆t and a slowly varying oneĨ 0 (∆t). According to the harmonic addition theorem [18], we can rewrite the second term as Re e iΩ∆tĨ 0 (∆t) = Re Ĩ 0 (∆t) cos Ω∆t − Im Ĩ 0 (∆t) sin Ω∆t = sgn Re Ĩ 0 (∆t) Ĩ 0 (∆t) cos Thus, sinceĨ 0 (∆t) is slowly varying, the visibility v of the oscillation of f (2) is given by where f (2) max (∆t) = max τ ∈[0,2π/Ω] f (2) (∆t + τ ), and analog for f (2) min (∆t). The visibility has the property v(0) = 1 and v(∆t) ≤ v(0), as it is indeed the absolute value of the characteristic function exp(iΩ∆t) of the distribution I 0 (Ω) which can be seen in (S48). Furthermore it is a direct measure of which is the absolute value of the Fourier transform of the overlap of the ideal laser spectrum and the absorption spectrum.
In the case of δ-shaped laser pulses, i.e., a flat laser spectrum, this quantity is proportional to the absolute value of the Fourier transform of the absorption spectrum This in turn is given by the ZPL-mirrored PL spectrum (see (S42)) in the case of a perfect spectrometer, i.e., where we used (S25) and (S109) to express the PL spectrum in terms of the phonon correlation function G −+ . Thus, we see that the visibility is a measure of the pure dephasing and of the phonon correlation function, which contains the information about the color-center-phonon coupling and the LO decay. The case of finite temporal pulse widths is slightly more complicated, leading to a finite time resolution. In this case, the visibility is no direct measure of the phonon correlation function since the finite pulse width leads to smeared-out visibility dynamics. Using the convolution theorem for Fourier transformations and the fact that the inverse Fourier transform of the laser spectrum is proportional to the autocorrelation of a single laser pulse (see (S28)) we arrive at which is the convolution between the laser pulse autocorrelation and the correlation function of the phonons G −+ times the correlation function e −iωZPLτ e − γ pd 2 |τ | of the polaron. Thus, a finite laser pulse width leads to a finite time resolution due to the convolution of the quantities of interest, i.e., the correlation functions, with the autocorrelation of the laser pulse. We note here, that there is at this level of approximation, i.e., lowest order in the electric field contribution, no special emphasis on the delays for which the pulses overlap. The convolution with the autocorrelation of the pulse leads to the same smearing out for every value of the delay ∆t.

Approximation for the visibility on short timescales for small detunings
In the context of slightly detuned excitation around the ZPL aiming for the influence of LA phonons on the visibility dynamics on short timescales, it is useful to have an approximative formula at hand to understand the data. We S16 derive this here by first separating the absorption spectrum into a absorption via the ZPL and an phonon-assisted absorption via the PSBs .

(S56)
We are interested in excitation of the emitter close to the ZPL, i.e., ω l ≈ ω ZPL . The spectral width of the laser pulse shall be of a few meV, such that for our purposes the influence of LO phonons in I abs,PSB can be neglected, since they have energies of > 150 meV. The LA PSB in hBN is rather flat on a scale of a few meV, such that we make the crude approximation I ideal laser (Ω)I abs,PSB (Ω) ≈ I ideal laser (Ω)I abs,PSB (ω l ). The visibility is then of the form v(∆t) ≈ R dΩ e iΩ∆t I ideal laser (Ω)I abs,PSB (ω l ) + 2πB 2 e iωZPL∆t I ideal laser (ω ZPL ) R dΩ I ideal laser (Ω)I abs,PSB (ω l ) + 2πB 2 I ideal laser (ω ZPL ) , where we inserted I abs,ZPL (Ω) = 2πB 2 δ(Ω − ω ZPL ) since γ pd = 0 is assumed. Using now the fact that the Fourier transform of the ideal laser spectrum is the autocorrelation of the pulse (see (S28)), we arrive at where we took the complex conjugate of the expression in the absolute value. To understand this formula a little better, we take a look at two aspects. (A) The dynamics on timescales, where the absolute value of the autocorrelation |[E s E s ] (∆t)| is more or less constant and (B) the dynamics for large delays when the autocorrelation has decayed. In the experiment, the autocorrelation decays within 0.5 ps to 1 ps. Up to this time, it is approximately given by This is strictly true for symmetric laser spectra and only approximately true for asymmetric spectra. The carrier frequency of the single pulse is also approximately the carrier frequency of the autocorrelation. Focusing now on aspect (A), we approximate the autocorrelation further as . (S59) Thus we see that it is of the form |ae −iδ∆t +b| a+b , where δ = ω l − ω ZPL is the detuning between the carrier frequency of the laser and the frequency of the ZPL. The visibility therefore exhibits a beat with the beat frequency determined by the detuning. The origin of this beat is the separate excitation of the polaron at frequency ω ZPL and of polaron-phonon states at frequency ω l = ω ZPL + δ, where the detuning δ together with the spectral width of the pulse determines the involved phonon frequencies. Now we turn to aspect (B), which concerns the long time behavior. By this we mean the behavior after the autocorrelation of the pulse has vanished but before the signal is diminished by pure dephasing, which we have neglected S17 here. We will denote this limit by ∆t(large) and the visibility then approaches the stationary value v(∆t) .

(S60)
The stationary value of the visibility is thus of the form b a+b . We can take a look at two limiting cases. The first is a b, such that the visibility approaches zero as stationary value. This limiting case can be reached by two different routes. The first is I ideal laser (ω ZPL ) → 0, i.e., an excitation that is strongly detuned such that the overlap of the laser with the ZPL is negligible. In this case the polaron is not excited directly and no coherence remains after the autocorrelation of the laser has vanished. The other route is B 2 → 0, i.e., strong phonon coupling leading to dominant LA sidebands. In this case also the absorption is dominantly via the PSB and not the ZPL, such that the polaron is only excited via phonon-assisted absorption resulting in a strong phonon induced dephasing. So the case v(∆t) 0 is a signature of negligible direct excitation of the polaron. Keep in mind that we neglected the pure dephasing here, which would lead to a vanishing of the visibility in any case. The other limiting case of the stationary value of the visibility is b a, such that the visibility approaches unity. For this limiting case we need a strong overlap of the laser spectrum with the ZPL and weak LA PSBs. In this way we excite the polaron alone, whose coherence later decays by pure dephasing. The above discussion shows that the short time visibility dynamics for excitation close to the ZPL are very sensitive to the detuning, which determines the overlap of the laser pulse spectrum with the ZPL. It is also very sensitive to the form of the laser spectrum, as this form is scanned by the ZPL when varying the detuning. Furthermore the stationary value of the visibility is determined by the relative strength between PSB excitation and ZPL excitation and thus the form of the LA PSB is important for the correct detuning dependence of the visibility.

Slow spectral jitter
The color centers in hBN are sometimes subject to telegraph-noise type jitter and are usually subject to slow gaussian jitter [19]. Here we discuss shortly how to incorporate spectral jitter into the calculation of the coherent control visibility. We consider only the case of slow jitter, i.e., the system jumps between different values for the ZPL frequency ω ZPL on a long timescale, such that every repetition of the coherent control experiment is performed on a fixed ZPL frequency and the emitter jumps between frequencies only in between the experiments. Thus, in the experiment the contributions for different ZPL frequencies are just integrated, such that the observable becomes (see (S43)) where the dependence of the absorption spectrum on the ZPL frequency is explicitly given in the argument. p(ω ZPL ) is the probability distribution of finding the color center with the frequency ω ZPL in a repetition of the experiment.
Defining the average absorption spectrum as we see that (S43) and (S61) are equivalent, but I abs (Ω, ω ZPL ) is replaced with its ω ZPL -averaged value. Thus we can calculate the visibility in the same way as in (S48) but with the averaged absorption spectrum. In this case, we need to calculate and normalize it. In the case of random telegraph noise, p(ω ZPL ) is a sum of two delta functions. In the case of Gaussian jitter, p(ω ZPL ) is a Gaussian.
In the following we explicitly consider a Gaussian jitter around ω ZPL , whose frequency spread is much smaller than the typical width of the laser spectrum. The probability distribution is given by where σ j is the spectral standard deviation. For the case of spectrally narrow Gaussian jitter we see from (S63) that the visibility is given by normalizing In the last step we used that the timescale 1/σ j of the Gaussian is much longer than the typical timescale of the laser pulse autocorrelation, i.e., Therefore, a spectrally narrow Gaussian jitter leads to a Gaussian envelope of the visibility on a long timescale and can be incorporated just by an additional Gaussian factor. The result still holds in the case of additional telegraph noise, which can here be incorporated by averaging I abs Ω, ω ZPL over another probability distribution for ω (0) ZPL . As discussed in the main text the telegraph noise results in a beat of the visibility with a frequency corresponding to the frequency jump.

Photoluminescence excitation
Closely related, at least on a conceptional level, to the coherent control experiment is the photoluminescence excitation (PLE) experiment, in which the emitter is excited by a laser pulse and the integrated photoluminescence is detected. By tuning the carrier frequency of the laser pulse one can gain information on the absorption properties of the color center. The observable is again the emitter occupation after the pulse and we can interpret the single pulse as a double pulse with zero delay, such that for the PLE we need to calculate f (∆t = 0) for different settings of the carrier frequency of the laser. In lowest perturbation order, the PLE signal is therefore (see (S43)) where we explicitly write the dependence on the carrier frequency ω l . Since the ideal laser spectrum is a function of the difference ω l − Ω, we see that the PLE is determined by the convolution of the absorption spectrum with the laser spectrum. From Eqs. (S46) and (S51) we see, that We thus need to perform essentially the same numerical calculations for the visibility dynamics and the PLE signal, namely: determine I 0 (∆t). The value at ∆t = 0 yields the PLE, whereas taking the absolute value and normalizing yields the visibility, as we see from (S48).

Laser interference
To obtain direct experimental information on the autocorrelation of the laser pulse, one can send a double pulse directly onto the spectrometer. The time and frequency integrated spectrum of the double pulse is given by (see S19 Eqs. (S26) and (S39)) . (S69) where we assumed a homogeneous spectral response R(ω, Ω) = R(ω − Ω) in the relevant frequency interval and E (+) (t) ∼ E(t). We see that this signal has the same form as f (2) (∆t) from (S43), but with I abs (Ω) ∼ 1. Therefore the time and frequency integrated double pulse spectrum also oscillates with the delay ∆t, like f (2) (∆t). The oscillation frequency is however here directly given by ω l and not determined by the overlap with the absorption spectrum of the emitter. We can retrace the steps leading to the visibility of the oscillations and find that it is here given by normalizing the expression Thus, in this way we can gain information on the absolute value of the autocorrelation of the single laser pulse used in the coherent control experiment. Apart from its spectrum, this gives additional valuable information on the laser pulse and aids in finding a sufficiently good theoretical fit.

E. Correlation functions
In the calculation of the optical spectra of the color center and of the observable of the coherent control experiment we encountered two related correlation functions (see Eqs. (S22) and (S34)), that are calculated here. We will see that the phonon spectral density is the central quantity to determine these correlation functions and thus also of central importance for the calculation of the considered optical signals. First we note that the color center is assumed to be initially in its ground state |G and the phonons to be in a thermal state ρ phon at a given temperature T . Thus, the initial total density matrix is of the form which is not affected by going into the polaron frame. The LO phonons are assumed to be in their vacuum state, i.e., b n ρ phon = ρ phon b † n = 0 for n ∈ {LO}, since they have energies > 150 meV and we are interested in temperatures below room temperature. The two correlation functions that we encountered are Since it holds that C(τ 2 , τ 1 ) * = C(τ 1 , τ 2 ) (same forC) we only need to calculate them for τ 2 ≥ τ 1 , while the other case can be inferred by complex conjugation. Both correlation functions can be evaluated by going first into the polaron frame and then employing the quantum regression theorem together with (S16) [12,20], which we will do here explicitly.
First we bring C in a form more similar toC by noting that H 0 and thus U 0 does not alter the occupation of the emitter and commutes with X † X. This can be seen by observing that H IB commutes explicitly with X † X, H pd leads to pure dephasing, which does not influence the occupation and H anh acts only on LO phonon degrees of freedom. Thus with X † XU 0 X † = U 0 X † and employing the cyclic property of the trace we can write

S20
Now we go to the polaron frame, which yields where ρ P (−∞) = |G G| ⊗ ρ phon is the initial reduced density matrix of the TLS-phonon system in the polaron frame. Now we will focus on the case τ 2 ≥ τ 1 and perform similar steps for both correlation functions. First we propagate ρ P and X † ρ P X from −∞ to τ 1 with (S16). This propagation leaves the phonon density matrix unaffected since the LA phonons are in a thermal state and the LO phonons in their vacuum state. Thus, L phon [ρ phon ] = 0. Likewise, L pol [|G G|] = L pol [|X X|] = 0, such that the reduced density matrices ρ P (−∞) and X † ρ P (−∞)X remain unchanged.
Then we need to propagate X † B + ρ P (−∞) = X † B + ρ phon from τ 1 to τ 2 in the case of C and XB − ρ phon in the case ofC. Observing that the polaron part of the propagation just yields a factor exp ±iω ZPL − γ pd 2 (τ 2 − τ 1 ) . The phonon part is more involved, but can be rewritten by remembering that the LO phonons are in their vacuum state, i.e., ρ phon b † n = 0 for n ∈ {LO}. Thus, This can be solved by introducing an effective non-unitary time evolution operator Thus the solution of (S16) at time τ 2 , when the initial reduced density matrix at τ 1 is given by XB − ρ phon , is while for the other case of X † B + ρ phon at τ 1 it is given by With this, the correlation functions take the form Here, G ∓,± (τ ) are phonon correlation functions, given by The properties C(τ 2 , τ 1 ) * = C(τ 1 , τ 2 ) andC(τ 2 , τ 1 ) * =C(τ 1 , τ 2 ) lead to G ∓,± (τ ) * = G ∓,± (−τ ), which can be seen in (S82) for the special case of γ pd = 0 and will be taken as the definition of G ∓,± (τ ) in the case of τ < 0. This can be seen by substituting ω = ω − ω n and shifting the lower limit of the integration to −∞. The exponential decay term then emerges as the Fourier transform of the Lorentzian factor in the spectral density. Thus, on this level, there is no difference between a decaying dispersion-less mode and a mode with a broadened spectral density, e.g., due to a non-vanishing dispersion.

Full phonon correlation function
The full phonon correlation function is a product of the LO correlation function and the LA correlation function. Note that for the thermal expectation value of the multi-mode displacement operators we have B = Tr [B ± ρ phon ] = B LA B LO , such that with . (S107) Since we will assume dispersionless LO modes, the LO spectral density J LO (ω) will be a sum of δ-functions. The LA spectral density J LA (ω) will be discussed in Sec. S2 F.

Phonon correlation function G+−(τ )
The phonon correlation function G +− (τ ), which appears in the correlation functionC in (S83), is given by This correlation function has the same structure as G −+ (τ ) in (S85), however the operators B ± are interchanged. If we change the sign of all coupling constants g n in the independent boson Hamiltonian H IB in (S5), we see from (S9), that this leads to an interchange of B + with B − . Since changing the sign of the couplings g n does not change the form of G −+ (τ ), which depends only on the spectral density and thus on |g n |, we conclude that G −+ (τ ) = G +− (τ ).

Full correlation functions -final result
The correlation functions C(τ 2 , τ 1 ) andC(τ 2 , τ 1 ) are now given by where G −+ is given in (S106) together with (S107) . We thus see that both correlation functions are completely determined by the phonon spectral density, the decay rates of the LO phonons, the polaron shifted transition frequency, and the pure dephasing rate of the exciton.
F. Modeling the LA phonon spectral density So far we have not considered a specific model for the LA spectral density. The derivation of the correlation function in Sec. S2 E does not require any structure of the coupling constants g n and we have seen that the phonon quantity of S24 interest for the coherent control experiment and the PL spectra is the spectral density J LA (ω), which we can model to fit the experiment best. Therefore, to keep the model as simple as possible, we employ a phenomenological spectral density of the form which is a standard super-ohmic spectral density with a Lorentz-Drude cutoff [20] and a sharp high-energy cutoff at a frequency Ω c describing the maximum acoustic phonon energy. By this cutoff also the polaron shift in (S11) is well defined. The frequency ω c essentially determines the width of the LA PSB in the PL spectra and α determines its height.

S3. IMPACT OF THE SYSTEM PARAMETERS
In this section we discuss the impact of the fitted system parameters in the theoretical results. Some of the intricacies of the fitting process, especially the influence of the spectrometer resolution and the spectral shape of the laser pulses, are highlighted.

A. Spectrometer response
To determine the spectral response R(ω, Ω) of the spectrometer, we use a spectrally narrow light source (light from a green laser diode), as discussed after (S26). The spectrum of this laser diode is proportional to R(ω, Ω) with Ω being the optical carrier frequency, which is given by the maximum of the laser spectrum. In the experiment, two different diffraction gratings are used, as discussed in Sec. S1 B. We refer to the spectrometer when using the 300 lines/mm grating as "300-grating" and when using the 1200 lines/mm grating as "1200-grating".  The spectra, depending on the wavelength of the detected light, are shown in Fig. S7 as dots. Figure S7(a) shows the case of the 300-grating and Fig. S7(b) the case of the 1200-grating. The spectra are approximated by Gaussians (lines), with a standard deviation of ∆λ 300 ≈ 0.3 nm in the case of the 300-grating and ∆λ 1200 ≈ 0.075 nm in the case of the 1200-grating. The spectrometer response in the spectral domain is thus given by with λ(ω) = 2πc/ω and analog for Λ(Ω). Since usually λ ∆λ, we can expand the exponent around λ ≈ Λ, or equivalently ω ≈ Ω into a Taylor series yielding Thus the spectrometer response can be approximated by the Gaussian . (S114) If we assume a constant resolution ∆λ in wavelength space, the resolution ∆ω in the frequency domain depends on the detected frequency ω.

B. Modeling the laser pulses
To find an accurate theoretical description for the laser pulse, we (i) employ (S26) for the detected laser spectrum using the spectrometer response from (S114) (keeping in mind that E (+) (t) ∼ E(t) for linearly polarized light) and (ii) use (S70) for the laser interference, detected in the coherent control setup. To retrieve the most suitable model for the laser pulses, we simultaneously use both, the spectral (i) and the temporal domain (ii) to get a reasonable fit for the laser pulse. We also found that the PL visibility dynamics in the coherent control experiment critically depend on the spectral shape of the laser pulses, when exciting close to the ZPL on the LA sideband (see Sec. S2 D 4 and the main text). This is discussed in more detail in Sec. S4 A. We achieve the best overall results by using an asymmetric Pearson-4 distribution [21] for the Fourier transform E s (ω) of the laser field of a single pulse E s (t) =Ẽ(t)e −iω l t with This distribution converges to a Gaussian with standard deviation 1/σ l for α l = √ 2m l /σ l and m l → ∞. Another interesting limiting case is m l = 1 and ν l = 0, where E s (ω) is a Lorentzian. For m l close to unity, the parameter α l is a measure for the width of the distribution. The parameter ν l is a measure for the asymmetry of the spectrum. For ν l < 0 the spectrum has a longer tail to the right, i.e., a positive skewness. The location parameter ω 0 is a measure for the position of the maximum. The true maximum however is given by ω 0 − α l ν l /(2m l ).
We define the carrier frequency of the laser ω l as the maximum of the distribution, as this frequency is easy to  Figure S8. Characterization of the laser pulse. (a) Measurement (dark) and simulation (bright) of the laser spectrum. The spectrum is modeled via (S115) with α l = √ 2m l /σ l , ω0 = ω l + α l ν l /(2m l ), σ l = 350 fs, m l = 1.2 and ν l = −0.5 √ m l . The spectrum is measured with the 1200-grating, such that the resolution for the spectrometer response in (S114) is ∆λ1200 = 0.075nm. (b) Measurement (dark dots) and simulation (bright line) of the laser interference visibility, i.e., the pulse autocorrelation. identify in the experiment. For the fit of the spectrum and the laser autocorrelation, we use (S115) together with α l = √ 2m l /σ l and ω 0 = ω l + α l ν l /(2m l ). A fit, leading to an overall good agreement with the experiment is achieved for the parameters σ l = 350 fs, m l = 1.2 and ν l = −0.5 √ m l . This choice of parameters leads to a value of α l ≈ (230 fs) −1 , which we refer to as a pulse of 230 fs duration in the main text. In Fig. S8(a) and (b) the simulated (bright) and detected (dark) spectrum and autocorrelation of the laser pulse are shown, respectively. For the detection of the spectrum, the 1200-grating is used and in the simulation we use the respective response in (S114) with ∆λ = ∆λ 1200 = 0.075 nm.
C. Fitting PL spectra of the emitter By fitting the detected PL spectra of the emitter for temperatures of 8 K and 80 K we can determine the spectral density of the acoustic phonons and the coupling strengths and lifetimes of the dispersionless optical modes. These parameters have to be chosen such that also the PL visibility dynamics in the coherent control experiment are consistently reproduced. The PL at 80 K is shown in the main text in Fig. 1. In addition we show the measured (dark red) and simulated (bright red) PL at T = 8 K in Fig. S9. In both cases, a constant background was subtracted from the experimental data. The PL spectra for both temperatures are well reproduced by the theory using the same set of parameters for the phonons, demonstrating consistency of the approach. The spectra are detected with the 300-grating, such that we model it via (S24) together with (S114) and ∆λ = ∆λ 300 = 0.3 nm. All phonon parameters are summarized in Tab. II and Tab. III. Note, that the spectral resolution of the spectrometer for a detection at Ω ≈ 2 eV is ∆ω = ∆λ Ω 2 /(2πc) ≈ 1 meV. Therefore, in the PL spectra we can disregard the Gaussian spectral jitter on an energy scale of ∼ 10 µeV, which is discussed in the main text and in Sec. S1 F. To illustrate the effect of the spectrometer resolution on the PL, in Fig. S9 we additionally show simulations of the PL spectrum for resolutions between ∆λ = 0.2 nm and ∆λ = 0.4 nm as a blue band. For a spectrometer with worse resolution (0.4 nm, upper edge of the blue band), especially the LA PSB gains intensity relative to the ZPL compared to the case of better resolution (0.2 nm, lower edge of blue band). This can be understood intuitively by considering the limit ∆λ → ∞, i.e., a spectrometer that cannot destinguish between different wavelengths. In that case the spectrum should be constant, i.e., the PSBs and the ZPL would have the same height and cannot be distinguished. In the other extreme case of ∆λ → 0, i.e., an ideal spectral resolution, we would be able to definitely distinguish between PSB and ZPL. The height of the LA sideband is mainly determined by the parameter α in (S111), which is directly connected to the deformation potentials when considering a more microscopic approach to the modeling of the LA spectral density [1]. This implies that one cannot accurately determine deformation potential constants from PL spectra, if the spectrometer resolution is not well known. exp. hom. inhom. hom.+inhom. Figure S10. Coherent control visibility for resonant ZPL excitation at T = 8 K. The data (dots) are the same as in Fig. 2 in the main text. The lines are fits to the data using (S116). The blue line is the best fit for the case σj = 0, i.e., no Gaussian jitter. The green line is the best fit for the case γ pd = 0, i.e., no pure dephasing. The red line includes both pure dephasing and Gaussian jitter.
In the context of Fig. 2 in the main text, the impact of Gaussian spectral jitter is discussed. In Fig. S10 we show the same data (dots), together with three fits (lines). All fits are simplified versions of the formula presented in (S65). Since the excitation is in resonance with the ZPL, we can as a first approximation neglect the influence of phonons and since the dephasing occurs on a much longer timescale than the laser autocorrelation, we can approximate (S65) by The simulation in the main text of course fully includes the effect of phonons and of the finite width of the laser pulse. The blue line in Fig. S10 now represents the best fit, when σ j = 0 is chosen, the green line represents the best fit for γ pd = 0, and the red line is the best fit for non-vanishing σ j and γ pd . We clearly see that the fit with Gaussian jitter and pure dephasing (red line) describes the data most accurately with σ j ≈ 36 µeV and 2/γ pd ≈ 55 ps. The presence of Gaussian jitter is also supported by the discussion in Sec. S1 F and the fitted value of σ j agrees with the values determined there.
E. Summary of the parameters and relevant formulas TABLE I. Parameters for modeling the laser pulses together with (S115), yielding the laser spectrum and laser interference dynamics presented in Fig. S8.
The spectrometer response R(ω, Ω) is modeled by (S114) with a resolution of ∆λ 1200 = 0.075 nm in the case of the 1200-grating and a resolution of ∆λ 300 = 0.3 nm in the case of the 300-grating. It is used in (S24) to calculate the PL of the emitter and in (S26) to calculate the laser spectrum. Whenever a measured PL is displayed together with a simulated one, the detection was performed with the 300-grating. The laser spectrum was always detected with the 1200-grating. The Fourier transform of the laser field of a single pulse is modeled by (S115). The parameters are given in Tab. I. With this, we can determine the laser spectrum using (S26), the laser interference using (S70), and most importantly the expressionĨ 0 (∆t) in (S54), which yields the normalized coherent control visibility and the PLE depending on the carrier frequency ω l for ∆t = 0. The LA spectral density is modeled by (S111) with parameters given in Tab. II. The LO spectral density is modeled by four dispersionless modes (sum of four delta-functions, see (S101)) with frequencies ω n and coupling constants g n given in Tab. III. The lifetimes 1/γ n = 1/γ(ω n ) are also presented there. These parameters are used in (S106) together with (S107) to calculate the phonon correlation function G −+ , which is relevant to all optical signals retrieved from the emitter. The ZPL frequency ω ZPL and the pure dephasing rates γ pd for 8 K and 80 K are given in Tab. IV. The spectral width σ j of the slow Gaussian jitter, which modifies all visibility calculations according to (S65) is also presented there, together with the parameters for the probability distribution describing the impact of random telegraph noise on the visibility dynamics in Fig. 4 in the main text. This probability distribution is modeled by a sum of two-delta functions and incorporated into the calculation of the visibility dynamics according to the discussion in Sec. S2 D 5. It is not necessary to use a normalized probability distribution, as the visibility is automatically normalized. Thus, A and B give the relative strength of the two frequency components of the random telegraph noise, which are separated by an energy of ∆E ZPL . All parameters given here are considered throughout the main text, leading to a consistent description of the experiment, including PL, PLE, and PL visibility signals. As mentioned in Sec. S2 D 4, the coherent control visibility for the case of slightly detuned excitation around the ZPL addressing the acoustic PSB, crucially depends on the form of the laser spectrum. Here, we present three auxiliary simulations to Fig. 6 in the main text.
1. Asymmetric Pearson-4 laser spectrum Figure S11 shows the same situation as Fig. 6 in the main text, however plotted on a slightly longer timescale and for more values of the detuning δ = ω l − ω ZPL . In the main text, the data was shown up to 2 ps to emphasize the non-trivial phonon induced dephasing dynamics occurring in the first picosecond of the delay. The experimental data is shown as dots, the simulation as lines. The left column shows the case of resonant excitation (violet) and negative detunings (red). The right column shows the case of positive detunings (blue). Simulation and experimental data match remarkably well, also for the cases not shown in the main text, i.e., δ = +2, ±3, +7 meV. By plotting data and simulation on a longer timescale of 3 ps, we highlight that the chosen value of the pure dephasing rate γ pd accurately describes the long-time dephasing behavior.  Figure S11. Coherent control visibility for excitation close to the ZPL, as in Fig. 6 of the main text. Experimental data as dots, simulation as lines. The case of resonant excitation (violet) and negative detuning (red) is shown in the left column, the case of positive detuning (blue) in the right.

Symmetric Pearson-4 laser spectrum
Now we consider the case of a symmetric spectrum of the Pearson-4 type from (S115) with the same parameters as in the main text (see Tab. I), but with ν l = 0. In Fig. S12(a) and (b) the simulated (bright) and detected (dark) spectrum and autocorrelation of the laser pulse are shown, respectively. For the measurement of the spectrum, the 1200-grating is used, such that in the simulation we use the response in (S114) with ∆λ = ∆λ 1200 = 0.075 nm. Compared to the case of an asymmetric spectrum considered in the main text (see Fig. S8), this laser pulse model yields a worse fit to the measured spectrum. The laser interference however still agrees well with the experiment. sim. Figure S12. Characterization of the laser pulses with a symmetric spectrum. (a) Measurement (dark) and simulation (bright) of the laser spectrum. The spectrum is modeled via (S115) with α l = √ 2m l /σ l , ω0 = ω l + α l ν l /(2m l ), σ l = 350 fs, m l = 1.2, and ν l = 0. The spectrum is measured with the 1200-grating, such that the resolution for the spectrometer response from (S114) is ∆λ1200 = 0.075 nm. (b) Measurement (dark dotted line) and simulation (bright line) of the laser interference visibility, i.e., the pulse autocorrelation.  Figure S13. Coherent control visibility for excitation close to the ZPL, as in Fig. 6 of the main text, but with a symmetric (ν l = 0) Pearson-4 pulse. Experimental data as dots, simulation as lines. The case of resonant excitation (violet) and negative detuning (red) is shown in the left column, the case of positive detuning (blue) in the right. Figure S13 shows the same data as Fig. S11, but the simulations include the symmetric Pearson-4 pulse. We see that the simulation in the case of negative detuning (left) tends to lie below the experimental data at ∆t = 3 ps, while it tends to lie above the data for positive detuning (right). When comparing the PID drop of the signal, i.e., the final PL visibility value at ∆t = 3 ps, for detunings of the same absolute value, we find almost the same values. For the asymmetric Pearson-4 pulse in Fig. S11 the difference between ±δ is larger. The reason for this is the asymmetry of the laser spectrum used in Fig. S11. As discussed in Sec. S2 D 4, the value of the visibility for large times is determined by the overlap of the laser with the ZPL and the height of the LA PSB at the position of the maximum of the laser spectrum. In the case of an asymmetric laser spectrum with more spectral weight towards higher frequencies (see Fig. S8), the overlap of the laser with the ZPL is enhanced for negative detunings resulting in larger visibility values at ∆t = 3 ps. When the asymmetry is lifted by reducing the high frequency flank of the spectrum the overlap with the ZPL reduces and the visibility values at ∆t = 3 ps shrink. The effect is less pronounced for positive detunings, because the lower frequency flank of the laser spectrum remains nearly the same when making the pulse symmetric (see Fig. S8(a)). Still in the case of a completely symmetric laser spectrum, the visibility dynamics for +δ and −δ behave differently. For positive detunings, phonon emission plays a more dominant role than phonon absorption and vice versa for negative detuning. Phonon absorption however can only take place if a thermal occupation n ph is present and the process scales with ∼ n ph . Phonon emission on the other hand includes spontaneous emission as well as induced emission and the full process scales with n ph + 1. Therefore an asymmetry, which can best be seen when comparing the complete dynamics in Fig. S13, e.g., for the case of δ = ±5 meV, remains even in the case of a symmetric laser spectrum.

Gaussian laser spectrum
A commonly used laser pulse spectrum to simulate the optical excitation of quantum system is given by a Gaussian with (S118) To test its performance in our situation we keep the temporal and spectral width comparable to the case considered in the main text and choose α l = (250 fs) −1 . In Fig. S14(a) and (b) we show the simulated (bright) and detected (dark) spectrum and autocorrelation of the laser pulse, respectively. For the detection of the spectrum, the 1200-grating is used, such that in the simulation we consider the response in (S114) with ∆λ = ∆λ 1200 = 0.075 nm. Neither the spectrum, nor the autocorrelation of the laser are well reproduced by this model, compared to the case of an asymmetric Pearson-4 spectrum in Fig. S8.    Figure S15. Coherent control visibility for excitation close to the ZPL, as in Fig. 6 of the main text, but with a Gaussian pulse (see (S118) and Fig. S14). Experimental data as dots, simulation as lines. The case of resonant excitation (violet) and negative detuning (red) is shown in the left column, the case of positive detuning (blue) in the right.
Figure S15 again shows the same data as Fig. S11, but using the Gaussian laser spectrum. While the simulations with the Pearson-4 spectra in Fig. S11 and Fig. S13 differed only slightly on a quantitative level, the visibility dynamics with the Gaussian pulse show a completely different qualitative behavior. Some of the cases, e.g., δ = −3, −4 meV fit well with the experimental data, while others, e.g., δ = −6 meV do not fit at all. The value of the visibility at ∆t = 3ps drops much faster than in the experiment, when going from δ = ±3 meV to = ±6 meV. This can be explained by the form of the Gaussian spectrum and the form of the Pearson-4 spectrum. The Gaussian spectrum drops ∼ exp −x 2 when moving from the maximum outwards. The Pearson-4 spectrum for our choice of parameters, specifically m l = 1.2, is close to a Lorentzian spectrum, which has more spectral weight far from the maximum. Consequently, the overlap with the ZPL is still significant for large absolute detunings in the case of the Pearson-4 spectrum and vanishes more quickly with the detuning for the Gaussian pulse.

B. Impact of LAs and laser asymmetry on the PLE
In this section we discuss the impact of the laser pulse spectrum on the PLE signal shown in Fig. 6(a) in the main text. Figure S16 shows the PLE spectrum of the emitter at T = 80 K with the experimental data as dots and the simulations as lines. (a) shows simulations including the coupling to acoustic phonons and no phonons are considered in (b). The general impact of the phonon coupling is discussed in the main text. Here, we focus on the additional impact of the spectral shape of the laser. The solid lines present the asymmetric Pearson-4 pulse with the parameters in Tab. I and the dashed lines the symmetric pulse with ν l = 0. In Fig. S16(a) the simulated PLE with LA phonons is asymmetric between positive and negative detuning E − E ZPL . For negative detunings, the PLE vanishes, whereas it reaches a plateau for positive detunings. The reason is again the asymmetry between absorption and emission of phonons. When going further from the ZPL in the direction of negative detunings, phonons with the energy ≈ |E − E ZPL | need to be absorbed in order to yield a PLE signal. The thermal occupation decreases for increasing phonon energy, such that the PLE signal decreases as well for increasing negative detuning. For positive detuning, spontaneous emission, as well as induced emission of phonons can generate a PLE signal. While the probability for induced emission decreases when increasing the detuning, spontaneous emission  processes can always take place, leading to the plateau for positive detuning.
If we now remove any influence of the LA phonons in Fig. S16(b), e.g., by setting α = 0 in (S111), the absorption spectrum of the emitter close to the ZPL is just a Lorentzian. In that case the PLE, being the convolution of the mirrored laser spectrum with the absorption spectrum (see (S67)), is just a broadened version of the laser spectrum, mirrored with respect to E − E ZPL = 0. Therefore an asymmetry of the PLE remains when removing the influence of LA phonons (solid line), since the laser spectrum is asymmetric (see Fig. S8). When additionally making the pulse spectrum symmetric, as depicted by the dashed line in Fig. S16(b), also the PLE spectrum becomes symmetric.
Thus we see that an asymmetry of the PLE in general stems from two influences, namely (i) the asymmetry between phonon emission and absorption and (ii) a possible asymmetry of the laser spectrum. Both have to be considered accurately to reproduce the measurement properly.