MUSE crowded field 3D spectroscopy in NGC 300 : IV. Planetary nebula luminosity function

We perform a deep survey of planetary nebulae (PNe) in the spiral galaxy NGC 300 to construct its planetary nebula luminosity function (PNLF). We aim to derive the distance using the PNLF and to probe the characteristics of the most luminous PNe. We analyse 44 fields observed with MUSE at the VLT, covering a total area of $\sim11$ kpc$^2$. We find [OIII]5007 sources using the differential emission line filter (DELF) technique. We identify PNe through spectral classification using the aid of the BPT-diagram. The PNLF distance is derived using the maximum likelihood estimation technique. For the more luminous PNe, we also measure their extinction using the Balmer decrement. We estimate the luminosity and effective temperature of the central stars of the luminous PNe, based on estimates of the excitation class and the assumption of optically thick nebulae. We identify 107 PNe and derive a most-likely distance modulus $(m-M)_0 = 26.48^{+0.11}_{-0.26}$ ($d = 1.98^{+0.10}_{-0.23}$ Mpc). We find that the PNe at the PNLF cut-off exhibit relatively low extinction, with some high extinction cases caused by local dust lanes. We present the lower limit luminosities and effective temperatures of the central stars for some of the brighter PNe. We also identify a few Type I PNe that come from a young population with progenitor masses $>2.5 \, M_\odot$, however do not populate the PNLF cut-off. The spatial resolution and spectral information of MUSE allow precise PN classification and photometry. These capabilities also enable us to resolve possible contamination by diffuse gas and dust, improving the accuracy of the PNLF distance to NGC 300.


Introduction
The planetary nebula luminosity function (PNLF) is a distance determination method with a precision and accuracy that is comparable to those of the tip of the red giant branch (TRGB) and Cepheid methods (Jacoby 1989;Ciardullo et al. 1989;Ciardullo 2010Ciardullo , 2012Roth et al. 2021). Using evidence from narrowband photometric surveys in [O iii]λ5007, Ciardullo et al. (1989) have shown that the magnitude distribution of planetary nebulae (PNe) for a given galaxy follows an empirical power law defined as where the brightest PN at the cut-off has an absolute magnitude of M * = −4.53 ± 0.06 with a possible minor dependency on metallicity (Jacoby 1989;Dopita et al. 1992;Ciardullo et al. 2002;Ciardullo 2012). While a number of different formulations have been developed to model the various shapes of the PNLF at fainter magnitudes (Rodríguez-González et al. 2015;Longobardi et al. 2013;Bhattacharya et al. 2019Bhattacharya et al. , 2021, such faint-end variation do not affect the definition of the PNLF's bright end cut-off, which is the critical feature for distance determinations (Spriggs et al. 2021;Ciardullo 2022). Until the early 2010s, most PNLF distance measurements were obtained using 4-meter class telescopes and narrow-band interference filters, and as a result, the method has been traditionally limited to distances of ∼ 20 Mpc (Jacoby et al. 1990;Ciardullo 2010Ciardullo , 2012Ciardullo , 2022. Although 8-meter class telescopes were available and even observed PNe at the Coma cluster (∼ 100 Mpc, Gerhard et al. 2005), most of the instruments had wider bandpass filters, which increased the inclusion of sky background signal. This limited the PN detection sensitivity, that was necessary to significantly improve the distance range of the PNLF (Ciardullo 2022). This situation has now changed due to the use of the Multi Unit Spectroscopic Explorer (MUSE, Bacon et al. 2010) integral-field spectrograph on the 8.2-meter Very Large Telescope to survey PNe in distant systems (Spriggs et al. 2020(Spriggs et al. , 2021Roth et al. 2021;Scheuermann et al. 2022). In fact, Roth et al. (2021) have shown that by using a differential emission line filter technique on MUSE data, PNLF measurements are now possible out to distances of ∼ 40 Mpc under excellent seeing condition and with the aid of the adaptive optics system. This is mainly due to the narrow effective bandpass of MUSE, Article number, page 1 of 20 arXiv:2301.03437v1 [astro-ph.GA] 9 Jan 2023 A&A proofs: manuscript no. pnlf_ngc300 that is five times narrower than the typical narrow-band filters, which can substantially suppressed the background sky noise (Roth et al. 2021).
Previous PNLF studies of late type galaxies using [O iii]λ5007 narrow-band filters were also hampered by the possible confusion with supernova remnants (SNRs) or H ii regions (Herrmann et al. 2008;Herrmann & Ciardullo 2009;Frew & Parker 2010). While Hα narrow-band image allows the exclusion of H ii regions, it cannot be used to exclude the SNRs, whose classification typically rely on the [S ii]λ6716, 6731 lines. In M31 and M33,  found that the SNR contamination does not change the shape nor the position of the PNLF cutoff. In contrary, Kreckel et al. (2017) have shown with MUSE observations of NGC 628 that the presence of SNR contaminants can affect the bright cut-off. However, Scheuermann et al. (2022) demonstrated that the latter study had an issue with the background subtraction in Hα, which affected the classification. Their reanalysis concluded that the PNLF cutoff was indeed unaffected by the contaminants. Nevertheless, the possible contamination by SNRs is relevant for the investigation of the faint end of the PNLF. The impressive spatial resolution and spectroscopic capability of the MUSE instrument allows the instant identification of interlopers, even in star forming disk galaxies Roth et al. 2018Roth et al. , 2021Scheuermann et al. 2022).
NGC 300 is a spiral galaxy in the foreground of the Sculptor group. Being fairly isolated from its neighbouring galaxies (Karachentsev et al. 2003) and close in distance (Gieren et al. 2005;Rizzi et al. 2006Rizzi et al. , 2007 makes it interesting for studying star formation histories and galactic evolution (Muñoz-Mateos et al. 2007;Kudritzki et al. 2008;Bernard-Salas et al. 2009;Gogarten et al. 2010;Jang et al. 2020). A previous PNLF study, Soffner et al. (1996) identified 34 PNe, a small number that was not ideal to make a proper PNLF and therefore opted the cumulative PNLF to derive the distance by using the LMC as a yardstick. More recently, Peña et al. (2012) observed 104 PN candidates using narrow-band imaging from the central and the eastern outskirt region to construct the PNLF, with a follow-up spectroscopy for the brighter candidates (Stasińska et al. 2013). In Paper I (Roth et al. 2018), seven 1 × 1 MUSE fields in the central region of NGC 300 were observed with the goal of resolving stellar populations in crowded fields of nearby galaxies, from which they discovered 45 PN candidates. Again, this number was too small to create a useful PNLF, since the sample spans a very wide magnitude range of 22 m 5007 29. In the present work, using publicly available archival data from McLeod et al. (2020McLeod et al. ( , 2021) -or ML20 -and 2 additional MUSE-GTO fields, we expand the observed area from 7 to 44 MUSE fields in order to detect more PNe and obtain a PNLF distance to NGC 300 using integral field spectroscopy.
Our lack of a complete understanding of the underlying physics behind the invariance of the PNLF cut-off has prevented the PNLF technique to become a primary standard candle (Ciardullo 2010(Ciardullo , 2012. Although simulations have provided an impression of the physical properties of the most luminous PNe (Jacoby 1989;Dopita & Meatheringham 1990, 1991Mendez & Soffner 1997;Méndez et al. 2008b;Schönberner et al. 2007Schönberner et al. , 2010Valenzuela et al. 2019), an observational characterisation is still limited to the LMC (Dopita & Meatheringham 1991;Dopita et al. 1992;Reid & Parker 2010a,b), and M31 (Kwitter et al. 2012;Galera-Rosillo et al. 2022). If the most luminous PNe at the PNLF cut-off have indeed originated from a single-star stellar evolution, then placing the central stars in the HR-diagram will provide insights into the underlying stellar population, and also the nature of the cut-off itself. Using the data quality that MUSE offers, we aim to constrain the luminosity and effective temperature of the central stars for some of the bright PNe to understand their origin and expand our understanding of PNe beyond the Local Group.
The structure of this paper is as follows: details on observations and data reduction are described in Section 2. The data analysis regarding the PN detection and classification, the literature comparison of the PN number, the [O iii]λ5007 photometry, and the measurement of the Balmer decrement is explained in Section 3. The resulting luminosity function and the distance measurement are presented in Section 4. The discussion and the implications of this work follow in Section 5. Lastly, the conclusions are given in Section 6.

Observations and data reduction
The data for this project were acquired using the MUSE spectrograph on the 8.2-meter Very Large Telescope (Bacon et al. 2010). 9 fields were obtained as part of the MUSE guaranteed time observation (GTO) program 1 , while 35 fields were taken from the ESO Archive 2 . The area covered by these observations is shown in Figure 1.
The initial MUSE-GTO data (field A, B, C, D, E, I, J) were obtained between the years 2014 -2016 using the extended wide field mode with a spatial coverage of 1 ×1 and spectral coverage of 4650 − 9350 Å with 1.25 Å sampling. First results from the 7 fields in the centre area were reported in Paper I, covering the nucleus, part of the spiral arm that extends from the nucleus to the north-west, and inter-arm regions of the galaxy. In late 2018, additional fields P and Q were observed to cover the part of the outer spiral arm. Moreover, fields A, B, and C were re-observed with adaptive optics support to obtain better image quality. In the adaptive optics mode, a notch filter at 5750 − 6100 Å blocks the laser light, which is, however, not affecting the emission lines of interest. Most of the fields were obtained with an exposure time of 6 × 900 s, with the exception of field J (4 × 900 s), field C (8 × 900 s), and field B, D (11 × 900 s).
The 9 MUSE-GTO fields were reduced using the MUSE pipeline (Weilbacher et al. 2020) within the MUSE-WISE environment (Vriend 2015), as explained in more detail in Paper III (Micheva et al. 2022). In addition to the field distortion correction produced by the pipeline, the astrometry is also calibrated using the Gaia DR2 catalogue (Gaia Collaboration 2018), providing absolute astrometry within 0 . 1. Sky subtraction was performed using an offset field outside the galaxy. Since we will perform photometry in [O iii]λ5007, we measure the seeing quality at this wavelength based on the FWHM of 3 to 4 stars for each field. These stars, which are typically giants or supergiants in the disk of NGC 300, have apparent magnitudes of F606W 21 in the HST ACS magnitude system (Roth et al. 2018). Unlike the situation in more distant systems, such as Fornax cluster ellipticals (Sextl et al. 2021), confusion with globular clusters is not a concern. For the MUSE GTO data, the image quality ranges from 0 . 6 − 0 . 8 FWHM, as presented in Appendix C.
Another 35 fields, the ML20 data, publicly available at the ESO Archive, were originally observed to study stellar feedback in NGC 300 (McLeod et al. 2020(McLeod et al. , 2021. The data were obtained using the nominal wide field mode, which has the same spatial Fig. 1. MUSE fields of NGC 300 are marked with red, green, and cyan. The red fields are MUSE-GTO data from the Paper I pilot study. The green fields labelled P and Q are additional MUSE-GTO observation for the outer spiral arm. The cyan fields indicate ML20 fields. The magenta fields are the previous PNe survey area of Peña et al. (2012) using the FORS2 instrument. Image: NGC 300 in Hα taken with the Wide Field Imager (ESO) -Program ID 065.N-0076 . coverage of 1 × 1 , but with a slightly shorter wavelength coverage of 4750 − 9350 Å. The observations were conducted in the period of 2016 -2018 without the support of adaptive optics, and each field was observed with an exposure time of 3×900s. These data were reduced using the fully automated MUSE pipeline (Weilbacher et al. 2020) with default parameters, as provided in the ESO Archive. The astrometry of the ML20 data only relied on the distortion correction within each field, which limited the absolute positional accuracy of the object catalogue to ∼ 3 . Moreover, the sky subtraction was performed using a reference region within each field instead of an offset field; this resulted in sky oversubtraction, especially in areas where diffuse gas is prominent. However, since we perform local sky subtraction for flux measurements of individual objects (see Section 3), the effect cancels out. Based on our measurements, the seeing quality of the ML20 data in [O iii]λ5007 ranges between 0 . 8 − 1 . 5 FWHM. These [O iii]λ5007 seeing measurements are presented in Appendix C.

PN detection and classification
To find PN candidates, we employed the differential emission line filter (DELF) method described by Roth et al. (2021). This is performed by extracting 15 datacube layers around the wavelength of redshifted [O iii]λ5007 (the systematic velocity of NGC 300 is v sys = 144 km/s; Lauberts & Valentijn 1989) and treating each layer as an on-band image; this 18.75 Å range accounts for the different line-of-sight velocities (LOSV) within the galaxy. Then, an intermediate broadband continuum image is constructed from the wavelength range between λ5063−5188 Å, which is free from strong absorption line features; this is used as the off-band image. By subtracting the scaled off-band image (see scaling factor in Equation 8, Roth et al. 2021) from the on-band images, we obtain a series of continuum-free differential images. Using the DS9 software (Joye & Mandel 2003), the differential images are visually inspected to find the [O iii]λ5007 sources. After experimenting unsuccessfully with DAOPHOT FIND (Stetson 1987) to identify PN candidates, which turned out to be unable to cope with the spatially variable emission line background in [O iii]λ5007, we resorted to the datacube layer blinking technique, that is described in Roth et al. (2021).
The typical physical size of planetary nebulae is of the of order ∼ 0.3 pc (Osterbrock & Ferland 2006). If we assume a distance of 1.88 Mpc (Gieren et al. 2005) and a scale of ∼ 9 pc/ , we expect the PNe in NGC 300 to appear as point sources. After marking the coordinates of the point sources in [O iii]λ5007, we apply aperture photometry (Stetson 1987) at each wavelength layer along the datacube for these objects to obtain their spectra. We employ an aperture diameter of 3 spaxels (0 . 6), an inner sky annulus of 12 spaxels, and an outer sky annulus of 15 spaxels. Although the seeing conditions of the datacubes vary between ∼ 0 . 6 − 1 . 5, we extract the spectra using the same parameters. The small aperture of 3 spaxels is chosen to minimise contamination of background gas or nearby H ii regions. Then, the line fluxes are extracted using Gaussian fitting with the LM-FIT routine in Python (Newville et al. 2016), keeping in mind that the MUSE data has a wavelength sampling of 1.25 Å. Since the MUSE-GTO and the ML20 data have overlapping areas, the classifications were done independently for each data set. Crossmatching was performed after the PNe were identified.
To classify the sources into PNe, H ii regions, and supernova remnants (SNR), we employed the BPT-diagram (Baldwin et al. 1981) that is based on the line ratio of [O iii]λ5007/Hβ and [S ii]λλ6716, 6731/Hα. Besides the application of classifying active galaxies (Kewley et al. 2001(Kewley et al. , 2006, the BPT-diagram has been demonstrated to effectively discriminate the PNe from their mimics (Kniazev et al. 2008;Frew & Parker 2010;Sabin et al. 2013;Roth et al. 2021). As the classification rely on emission line ratios with very similar wavelengths, we can assume that the relative line fluxes have negligible extinction and seeing dependency on wavelength.
Our BPT-diagrams for the MUSE-GTO and ML20 data are presented in Figure 2. To separate the SNRs, we adopt the value log [S ii]λλ6716, 6731/Hα ≥ −0.5 from Roth et al. (2021). To discriminate PNe from H ii regions, we employ the theoretical line by Kewley et al. (2001), which was originally intended to differentiate starburst galaxies. We also consider the line ratio of [S ii]λ6731/6716 as a proxy for density, since bright PNe are expected to be denser than both H ii regions and SNRs (Osterbrock & Ferland 2006). While the brightest [O iii]λ5007 sources have sufficient line fluxes for the BPT-diagram classification, fainter sources may lack the weaker emission lines, i.e. the Hβ line or the [S ii]λ6716, 6731 lines. In such cases, we assume lower limits for the line fluxes and classify the sources as PNe if the [O iii]λ5007 line is stronger than the Hα line, which may introduce deviation from the separation lines in the diagram. We also found some faint objects, which only have the detection of the [O iii]λ5007 and the [N ii]λ6584 line without Hα detection, which are possibly Type I PNe (Frew & Parker 2010). Moreover, we also put remarks for PNe, which are only classified solely based on [O iii]λ5007 detection. This is the case for few of our faintest PN candidates. Nevertheless, such cases will not affect the distance determination because the PNLF cut-off is only defined by the brightest PNe.
In the MUSE-GTO data, we classified 37 PNe, 62 H ii regions, and 59 SNRs. In ML20 data, we classified 85 PNe, 176 H ii regions, and 105 SNRs. To cross-match the PNe candidates in the overlap area between the two dataset, we attempted an automated algorithm by comparing the sky coordinates. However, since the astrometric accuracy of both data differs, our attempt was not successful. Therefore, the cross-matching was performed visually using the DS9 software (Joye & Mandel 2003). The final PN number from the MUSE-GTO and ML20 fields: 105 PNe in the central region, and 2 in the P and Q fields at higher galactocentric distance. The PN catalogue is presented in Appendix D.

PN number comparison
Previous PN surveys of NGC 300 were conducted by Soffner et al. (1996) MUSE-GTO data (upper) and ML20 data (lower). The orange dot-dashed line is taken from Kewley et al. (2001) and the purple dashed line is defined by Roth et al. (2021). Open circles indicate PN candidates, which only have [O iii]λ5007 as the diagnostic line for the diagram. The deviation from the separation lines is explained in the text. more area and contains more PNe than the other studies. PE12 observed NGC 300 with the FORS2 imager (Appenzeller et al. 1998) in two 6.8 × 6.8 fields, one in the centre, and another in the eastern outskirts of the galaxy. The study employed the on/off-band technique to detect PN candidates and classified objects based on the criterion of whether or not a central star was present in their 5105 Å image. The expectation was that central stars of PNe would be too faint to be detected in the visual, while the ionising O stars in H ii regions could be seen. For brighter candidates with m 5007 < 25, they also performed additional spectroscopy using the MXU-mode with the same instrument (Stasińska et al. 2013). Since our observations cover a smaller area on the sky, we only made the comparison for the intersecting 5 × 7 region in the centre of the galaxy, which is also indicated in Figure 1.
In the overlapping region at the centre of the galaxy, we identify 105 PNe, compared to 58 in the PE12 sample. Moreover, although we recover all 58 sources found by PE12, our classification indicates several discrepancies. While 43 of the PE12 sources are confirmed as PNe, we classify 9 objects as compact H ii regions and 3 as SNRs. These misclassifications could have happened due to the fact that PE12 only have the spectral classification for candidates with m 5007 < 25, while the fainter objects completely relied on the detection of a central star. This approach also lacked the ability to identify SNRs amongst the fainter candidates, as such objects can be discriminated through the detection of the [S ii]λ6716, 6731 lines (Frew & Parker 2010;Sabin et al. 2013). Since all of our candidates are classified on the basis of their spectral properties, we believe that our classification is more reliable. Moreover, in terms of the number of detection, we also demonstrate that the MUSE observations are more sensitive and able to reach fainter magnitudes.

[O III]λ5007 photometry
The [O iii]λ5007 fluxes were obtained using DAOPHOT aperture photometry (Stetson 1987), applied to the PNe candidates in the 15 differential layers for each datacube. Then, the magnitudes were computed using the V-band equivalent conversion (Jacoby 1989) defined as where the flux is in erg cm −2 s −1 . Here, the aperture radius was adjusted to a value of approximately the FWHM of the PSF in a given exposure to accommodate the respective seeing condition. The inner and outer sky annulus were fixed to 12 and 15 spaxels, respectively. Most of the flux of the PSF was obtained by adding the 5 bins closest to the Gaussian peak and the remaining flux is recovered through the use of an aperture correction based on the information of a PSF reference. This correction is crucial to obtain accurate fluxes, which however can be a challenge when there is no reference available especially with the small field of view of MUSE. The aperture correction method is explained in Appendix A. The photometric uncertainty was calculated from the Gaussian fit errors, convoluted with an assumed flux calibration error of 5% (Weilbacher et al. 2020). In high surface brightness regions of distant galaxies, double-peaked profiles are occasionally found and indicate the presence of two superposed PNe with different radial velocities (Roth et al. 2021). Unsurprisingly, we do not find such cases in our sample. Since NGC 300 is a quite nearby galaxy, spatial coincidences are less likely. Additionally, the five datacube layers containing the total flux for [O iii]λ5007 were also inspected to insure the PN candidate was not extended or contaminated by surrounding gas emission.
As an internal test of our photometry, we used the regions of field overlap to compare our PN measurements made in the MUSE-GTO fields to those from the ML20 data. This comparison is presented in the upper panels of Figure 3. We find that the ML20 observations obtained under poor seeing condition tend to be systematically fainter than the MUSE-GTO data, while the photometry of the same object from different datacubes with similar image quality gives identical results (exception for the faintest PN in the comparison). Thus, the difference in seeing conditions can introduce a magnitude error; this is most likely due to the choice of a too small aperture for the asymptotic assumption for the aperture correction. Because the disk of NGC 300 contains a large amount of diffuse emission-line gas, we chose to not to increase this radius. However, in order to obtain the same photometric quality between the two data sets, we applied and additional corrections of 0.2 mag for objects with seeing FWHM ∼ 1 . 2 (∼ 6 spaxels), and 0.3 mag for seeing FWHM ∼ 1 . 4 (∼ 7 spaxels). We found these values based on empirical trial and error to achieve the minimum average discrepancy between the two sets of magnitudes.
The comparison after applying the correction is presented in the lower panels of Figure 3. The average discrepancy is now 0.05 mag, which is still within the typical measurement error of 0.06 mag. For PNe with m 5007 ∼ 27 or fainter, the correction has no meaningful implication, because the candidates are close to the detection limit. In the pilot study, Roth et al. (2018) performed a completeness simulation for the MUSE-GTO data of given exposure time, with seeing quality ranging from 0 . 6−1 . 2. This was conducted by embedding artificial PNe with different magnitudes into the real datacubes. It was found that for a seeing of 1 . 2, the expected completeness is 90% at m 5007 = 27. Although we are able to detect a PN as faint as m 5007 = 28.91, the seeing quality on average for all our data is 1 . 0, with 23% of the fields exhibiting larger than 1 . 2 FWHM. This shows that completeness of 90% at m 5007 = 27 is only achieved for 77% of our fields. However, since the emphasis of this work is on the bright candidates that define the PNLF cut-off for the distance determination, our results are not suffering from sample incompleteness at the faint end. For the final [O iii]λ5007 magnitudes, we preferred the MUSE-GTO data, if available, and otherwise we employed the ML20 data. We also applied the correction for fields outside the overlapping area with seeing FWHM > 1 . 2. In total, 11 PNe from 4 fields were corrected in this manner.
To test the accuracy of our photometry, we compared our magnitudes to the results from the literature. Figure 4 shows a comparison with SO96 and PE12. While our data is in reasonable agreement with SO96 within 0.01 mag on average, there is a systematic offset with regard to PE12. We find that our magnitudes are systematically brighter by an average of 0.71 mag. PE12 obtained instrumental [O iii]λ5007 magnitudes for the FORS2 on-band image using aperture photometry with the aperture diameter of 5 pixels (1 . 25), based on the average PSF FWHM of 2.9 pixels. To obtain the apparent m 5007 magnitudes, they calibrated the instrumental measurements using an empirical relation derived from the objects' spectroscopic fluxes, which are only available for the brightest PNe in their sample. We can try to understand what may be the reason for the discrepancy. First of all, we note that flux calibration is an established MUSE procedure in operation at the VLT and part of the data reduction pipeline. According to Weilbacher et al. (2020), flux calibration has been measured to be accurate to within 3-5%. If a significant number of our MUSE exposures would have been affected by non-photometric observing conditions -for which we have no evidence, we would expect a scattered, but not the tightly constrained linear correlation, that we see in Figure 4, in particular for the brightest 3 magnitudes. Secondly, Roth et al. (2018) have tested synthetic MUSE datacube broadband photometry of stars against published HST photometry for the same GTO datacube subset that has been used in our work, showing no hint of an offset to within a magnitude of F606W=22.7. Thirdly, the agreement with SO96, who obtained their data with narrow-band imaging at the ESO NTT, i.e. a different instrument at a different telescope, gives us reasonable confidence that our flux calibration cannot be off by as much as a factor of almost 2. Finally, we can follow the argument put forward by Roth et al. (2004), that by definition, integral field spectroscopy is an ideal tool for spectrophotometry, as is does not suffer from any kind of slit effects.
We can speculate, though, that the spectrophotometry from PE12 might have been affected in various ways to give rise to the observed systematic offset. In case of PE12, the calibration relies on spectroscopic fluxes, which were obtained using a slit spectrograph, and thus may suffered from slit-losses that were estimated to be 10 − 15%. However, from our comparison, we infer that the loss might be underestimated since 0.71 mag difference is equivalent to a loss of ∼ 48%. To test this, we performed a slit-loss simulation based on a model of the PSF with the quoted seeing conditions of 0 . 7 − 0 . 9, and a slit size of 1 . The simulation is done in the R-band, as the seeing measurement is typically done in this band. Based on these parameters, our simulation predicts that the slit-losses should be between 8 − 20%. However, since the PSF FWHM is expected to be larger in the blue wavelength region, the slit-loss in [O iii]λ5007 will be larger than that in the R-band. Moreover, additional losses can be introduced by positioning and guiding errors, as investigated by Jacoby & Kaler (1993). Spectrophotometry with a slit spectrograph also requires the slit to be oriented at the parallactic angle to minimise the effect of atmospheric dispersion (Filippenko 1982;Jacoby & Kaler 1993). Since our observations are performed with an IFU, we are not affected by any of these problems. While it is possible that our use of small apertures has caused some flux to be lost, we are able to compensate for this loss using aperture corrections as described above. Such a procedure is not easily performed for data observed with a slit spectrograph.
Besides the issue of slit-losses, the follow-up spectroscopy by PE12 is limited to PNe with m 5007 < 25, which is less than 40% of their whole sample. This implies that most of their PNe are dependant on the measurement accuracy of the brighter PNe, which are likely to be affected by systematic errors. Moreover, the use a 5 pixel aperture to measure the flux for a 2.9 pixel FWHM PSF incurs a risk of including light from background contamination. In such crowded fields with a variable background and ubiquitous diffuse emission-line gas, the aperture might unexpectedly collect [O iii]λ5007 flux of the ambient interstellar medium. Without proper background inspection and subtraction, this might lead to an overestimation of brightness. The inclusion of background emission likely explains the scatter for m 5007 > 25 in Figure 4. The spatial resolution of MUSE allows us to carefully check and analyse the condition of the background on a case-by-case basis and provide more accurate photometry. The variable background is also the main consideration to opt for a smaller aperture size for our flux measurements, and to rely on the aperture correction to deliver the final values.

Balmer decrement
Measurement of extinction using the Balmer decrement with Hα/Hβ ratio have been demonstrated on MUSE data for different objects, i.e. Pillars of Creation in M16 (McLeod et al. 2015), core of R136 in the LMC (Castro et al. 2018), faint H ii regions in NGC 300 (Micheva et al. 2022). To obtain this, we employed the spectra extracted for the classification, as explained in Section 3.1, using the aperture of 3 spaxels, with the inner and outer sky annulus of 12 and 15 spaxels, respectively. We also apply aperture correction for the Balmer lines, which can be referred to in Appendix A.
However, not all of our PNe candidates are detected at these two wavelengths. In order to filter out the candidates, we put a threshold of F(Hα) = 2 × 10 −17 erg cm −2 s −1 . For typical PNe with electron temperature T e = 10.000 K, the expected Balmer ratio is Hα/Hβ = 2.86 (Osterbrock & Ferland 2006), which corresponds to F(Hβ) ∼ 8.75 × 10 −18 erg cm −2 s −1 for the Hα threshold. Any Hβ flux lower than the threshold is too close to the detection limit. For such cases, we assume the upper limit of Hβ flux derived from the Hα line, which consequently also assumes no extinction. To avoid possible biased exclusion of high extinction PNe, we flag the objects with upper limit Hβ flux. If the Hα flux is below the threshold, then the extinction measurement is not performed. Based on the Hα threshold criterion, we have a complete sample for PNe down to m 5007 = 24.5. If we extend the sample to fainter magnitudes, we reach 87% completeness until m 5007 = 26.0 and 64% completeness until m 5007 = 27.0. To calculate the extinction, we then used the Balmer decrement defined as where k(λ) is the wavelength dependant extinction constant. For the foreground extinction, we employed the extinction curve of Cardelli et al. (1989)  , we employed the average LMC extinction curve to obtain the extinction of our PNe. The uncertainty of our extinction measurement is highly dependent to the aperture correction method. Therefore, we quote an estimated error based on the comparison of extinction calculated from different aperture correction methods, as explained in Appendix B. In spiral galaxies, Herrmann & Ciardullo (2009) found that the typical extinction for the PNe in [O iii]λ5007 is A 5007 ∼ 0.7. However, we discovered three high extinction cases with A 5007 > 1.5, including one with an extreme value of A 5007 ∼ 3.3. While it is possible that some PNe exhibit high intrinsic extinction, as high metallicity populations and massive progenitors tend to produce dustier PNe (Stanghellini et al. 2012), we suspect that the Balmer decrement might not always be accurate due to the local contamination. To investigate this further, and to highlight possible pitfalls that may play a role in studies based on slit spectroscopy, we examined spatial maps of the high-extinction objects in the wavelengths of Hβ, [O iii]λ5007, and Hα, using the p3d software . We found that these PNe candidates are co-spatial with nearby H ii regions. In Figure 5, a PN is shown to be an isolated point source in [O iii]λ5007. However, in the spatial map of Hα, the point source is entirely embedded inside the extended emission surface brightness distribution of an unrelated nebula. This clearly shows that the Balmer line flux of the PN candidate is contaminated, and an accurate extinction measurement of the PN itself cannot be obtained. All three of the objects in question show similar patterns of contamination. We therefore excluded them from the sample.
We compared our PN extinction measurements with the results from Stasińska et al. (2013) -also referred as ST13 -who observed PNe in NGC 300 with the FORS2-MXU instrument at the VLT (Appenzeller et al. 1998). They used 3 grisms 600B, 600RI, and 300T to cover spectral ranges of 3600 − 5100Å, 5000 − 7500 Å, and 6500 − 9500 Å, respectively. To avoid uncertainties from the flux calibration of different bands, they employed the Hγ and Hβ lines from the 600B grism spectra to measure the Balmer decrement. The extinction values for 18 PNe in common are presented in Table 1.
We have contemplated several reasons to explain the discrepancy. Firstly, slit losses that occur for the measurement of [O iii]λ5007 most likely also occur for the Hγ and Hβ lines. The short baseline between the lines is also very sensitive to systematic errors, making it difficult to derive accurate extinction values. Moreover, higher order Balmer lines are typically weak. The accuracy of measuring their flux depends on the precision to which the background of stellar absorption line spectra can be subtracted (Jacoby & Kaler 1993;Roth et al. 2004). Using the PMAS instrument, Roth et al. (2004) compared the accuracy of the Hβ line flux obtained with the IFU and simulated slits. They demonstrated that the orientation of the slit can introduce a different sampling of the background, leading to systematic differences of the derived flux measurements. Since ST13 performed their measurements with slit spectroscopy, they were susceptible to this type of error.
For cases where the internal extinction of ST13 is reported higher than our values, we find that 3 of their PNe are discarded from our sample (L6-5, H1-1, and H12-1 in Table 1), either because of severe contamination as shown in Figure 5, or by their low excitation, which we consider typical for compact H ii regions (Frew & Parker 2010). We also found that in some of these cases, the PN is embedded in diffuse gas. In our sample, if the diffuse gas is assumed to be uniformly distributed, the flux excess can be corrected using the background sky annulus. It remains an open question whether the background correction of diffuse gas was accurately accounted for in the slit spectroscopy of ST13, but we conclude that a careful consideration of back-A&A proofs: manuscript no. pnlf_ngc300

The PNLF
The PN luminosity function of this work is presented in Figure  6. It exhibits the dip between 1 and 3 magnitudes below the cutoff. Such dip is typically observed in star-forming galaxies ( To determine the distance, we employed the maximum likelihood technique, where the empirical PNLF is treated as a probability function (Ciardullo et al. 1989), assuming M * = −4.53 ± 0.06 and a fixed slope parameter of 0.307. When the number of PNe at the bright end cut-off is less than ∼ 50, distance determinations based on χ 2 minimisation depend significantly on the details of how the PNLF is binned. Such methods are not recommended (Ciardullo et al. 1989;Roth et al. 2021). Although our observation extends to m 5007 ∼ 29, the PNLF fit is only performed for the sample brighter than the dip until m 5007 = 23.6, since equation (1) does not consider the dip feature, which nevertheless is insignificant for the distance determination (Spriggs et al. 2021;Ciardullo 2022). By taking the foreground extinction of E(B − V) = 0.011 (Schlafly & Finkbeiner 2011) into account, the most-likely distance modulus is (m − M) 0 = 26.48 +0.11 −0.26 with the uncertainties representing the statistical error of the fit and the M * uncertainty.
We also calculated the distance using PE12 photometry to make a comparison. Assuming a completeness limit of m 5007 = −0.20 , a value that is significantly larger than our MUSE distance. We argue that this is due to the systematically fainter magnitudes of PE12 photometry, as discussed in Section 3.3. It is important to mention that PE12 measured a modulus distance of (m − M) 0 = 26.29 +0.12 −0.22 , a value much smaller than our maximum likelihood values, including our distance measurement using the PE12 data. This difference is possibly due to their use of the Levenberg-Marquardt χ 2 fitting technique, which is dependant on the binning method to construct the luminosity function. Since the sample size is limited, they employed rather wide magnitude bins of 1.16 mag, which did result in a luminosity function shape that closely resembles the empirical law. However, in the luminosity function of PE12, their first magnitude bin is located at m 5007 ∼ 22, despite the fact that their brightest PN has a magnitude of m 5007 = 22.99. Thus, when the fit is performed, this systematic shift to brighter magnitudes results in a smaller distance modulus. This demonstrates that the choice of bin size can produce unintended systematical shifts of the luminosity function when the number of PNe in the top ∼ 0.5 mag of the luminosity function is small. Similarly, PE12's choice of bin size also smeared out detail in the PNLF's shape, as they did not report the observation of the PNLF dip. In Figure 6, we present the PNLF that we plot with the original data from PE12 using higher binning resolution than the original work. In fact, the dip is present in the PE12 data, confirming that it is not an artefact in our measurements.
Finally, PE12 employed a larger extinction correction for the photometry with A 5007 = 0.2 compared to our value of A 5007 = 0.05. PE12 assumed this as the intermediate value be-tween found by Gieren et al. (2005) with E(B − V) = 0.096 (A 5007 = 0.3) and Schlegel et al. (1998) with E(B − V) = 0.013 (A 5007 = 0.05). In the case of Gieren et al. (2005), the extinction value is the sum of both the foreground extinction of Schlegel et al. (1998) and internal extinction derived from the Cepheids. For Cepheid distances, the internal extinction correction is necessary since they are originated from Population I stars, that are typically surrounded by galactic dust. However, this is less true for the PNe, so the foreground extinction correction for the PNLF distance is sufficient (Ciardullo 2010). This implies that the extinction correction A 5007 = 0.2 by PE12 is overestimated and also contributes to the smaller distance modulus. Therefore, the discrepancy between our calculation of the PE12 data and the original calculation is traced back to the issue of binning a limited sample and also the extinction correction.

PNLF Distance
To demonstrate the accuracy of our distance measurement, we compare our result with previous distances in the literature derived using Cepheids and tip of the red giant branch (TRGB), taken from NED, in Figure 7. We can see that most of the distances are within the uncertainties of our PNLF result. One aspect that may introduce the discrepancy is the correction of extinction. For instance, the Cepheid distance of Gieren et al. (2005) and the TRGB distance of Rizzi et al. (2006), both part of the Auracaria project, are corrected with foreground and internal extinction of E(B − V) = 0.096. However, Rizzi et al. (2007) argue that the extinction derived from dusty young Cepheids by Gieren et al. (2005) are not representative for the whole galaxy; their result only applies the foreground component of extinction. This shows the importance of having the same zero-point when comparing different distances derived from different methods.
Moreover, we also show that the MUSE observation, combined with the differential emission line filter (DELF, Roth et al. 2021) and maximum likelihood technique (Ciardullo et al. 1989), has improved the accuracy of PNLF method, as shown in Figure 7. The early result by Soffner et al. (1996) is based on the limited sample of only 34 PNe, from which they only construct a cumulative PNLF and employed the distance modulus of the LMC as a yardstick. Peña et al. (2012) identified a significantly larger sample of 104 PNe, but as shown in Section 3.3, their data may suffer from slit-losses and contamination. The systematically fainter PN magnitudes then led to larger distance modulus, as described in Section 4. Since the cut-off of the PNLF of NGC 300 is defined by a very small number of PNe, minimisation fitting methods become too dependant on the binning (Ciardullo et al. 1989). In a study by Jacoby (1997), a correction for PNLF distance based on the number of PN sample is suggested. For a PNLF cut-off sample < 20 PNe, they estimated a distance correction of ∼ 0.1 mag (see Figure 5 in Jacoby 1997). However, since the Cepheid and TRGB distance also varies with standard deviation of 0.1 mag, there are no solid distance reference to test if the correction is appropriate. Nevertheless, we have shown that the PNLF distance derived with the maximum likelihood technique is more robust. We take this as a motivation to improve PNLF distance measurements for nearby galaxies with our method.

Local dust effect on PN extinction
Dust formation plays important role in the early stages of PN evolution since it occurs at the surface of the progenitor AGB star and presumably plays an important role in the envelope ejection (Herwig 2005;Stanghellini et al. 2012). Infrared studies in the Milky Way and the LMC have revealed that the dust production is dependant on metallicity, with dustier systems found in higher metallicity environments (Stanghellini et al. , 2012Bernard-Salas et al. 2009). Although it cannot tell the properties of the dust, Balmer decrement extinction measurements can also probe the presence of dust in PNe. In the study of , a comparison was made between the PN extinction distribution in the bulge of M31 and several other galaxies: the LMC (Reid & Parker 2010a), NGC 4697 (Méndez et al. 2008a), and NGC 5128 (Walsh et al. 2012). Despite the limited samples involved, the authors found that the average extinction of PNe in each galaxy roughly follows the metallicity of the system.
To investigate such trends in NGC 300, we plot the extinction distribution in [O iii]λ5007 for our PNe until m 5007 = 23.6 (15 PNe) in Figure 8. These PNe are the ones employed for the maximum likelihood distance measurement. We find that in general these bright PNe have low extinction in [O iii]λ5007. The average extinction value for this sample is A 5007 = 0.31 (c(Hβ) ∼ 0.09), which is lower than the average of the bright PNe sample in the LMC with A 5007 = 0.57 (Reid & Parker 2010a;). However, we refrain from further interpreting the extinction distribution with the PN dust production, since the distribution is likely to be affected by local dust clouds, which can vary from one object to another. Such a problem has been reported in NGC 5128, where the high extinction of some PNe was attributed to local dust clouds rather than the PNe themselves (Walsh et al. 2012).
At the distance of NGC 300, our MUSE observations offer a spatial resolution of between 6 and 14 pc. This resolution should be sufficient to visually resolve the spatial variation of dust extinction (Kreckel et al. 2013;Tomičić et al. 2017). To test this, we inspected several objects with high extinction. As an example, we present the spatial map in [O iii]λ5007 and RGB colours, which is constructed from Johnson-VRI filters, for the PN with the highest extinction value (PN A-23, A 5007 = 1.18) in Figure  9. Although it is not obvious in the [O iii]λ5007 image, the RGB image shows a dust lane patch, extending from the lower left corner to the centre. Since PN A-23 is in proximity to the dust lane, we suggest that the measured high extinction of this object is composed of both local dust within the galaxy and circumnebular extinction associated with the PN itself.
Based on the comparison study between Balmer decrement extinction and infrared dust distribution in M31, Tomičić et al. (2017) concluded that vertical distribution of diffuse interstellar gas (DIG) and dust can vary in different locations of the galaxy and thus cause differing amounts of extinction. For NGC 300, variation of extinction also has been reported by Roussel et al. (2005). Therefore, there is currently no guarantee that the measured extinction of individual PNe is free from local effects, which is confirmed with our images. We must therefore refrain from making conclusions based on the extinction values alone, until the different components of the extinction can be quantitatively resolved.
Although the extinction of individual PNe might be affected by local dust lanes, such effects are less significant for the luminosity function as a whole. The effect of dust scale height in the PNLF distances of late-type disk galaxies has been discussed by Feldmeier et al. (1997). They modelled PNLF with varying ex-A&A proofs: manuscript no. pnlf_ngc300 Fig. 7. Distance modulus difference between our PNLF result with Cepheids (blue triangles) and TRGB (red squares), obtained from NED and sorted based on publication date. The Cepheid distances are from Willick & Batra (2001); Paturel et al. (2002); Gieren et al. (2004Gieren et al. ( , 2005; Saha et al. (2006); Bono et al. (2010); Bhardwaj et al. (2016). The TRGB distances are from Butler et al. (2004); Sakai et al. (2004);Tikhonov et al. (2005); Tully et al. (2006); Rizzi et al. (2006Rizzi et al. ( , 2007; Jacobs et al. (2009);Dalcanton et al. (2009). Previous PNLF distances of Soffner et al. (1996) and Peña et al. (2012) is also presented (green circles). The green shadow indicated the uncertainty of our PNLF distance.  tinction in [O iii]λ5007 and concluded that the inferred distance modulus should always be within 0.1 mag of the derived distance without extinction. A similar result also obtained by Rekola et al. (2005), who modelled the PNLF with different scale heights of dust in the starburst galaxy NGC 253. They found that even when the disk was optically thick with 1 mag of extinction, the PNLF distance is robust to within 0.1 mag. Both studies suggest that the brighter PNe tend to be located above the dust layer from the point of view of the observer, or for other reasons suffer little extinction from within the galaxy. With these arguments, we do not expect the occurrence of dust lane extinction to significantly affect our distance result and a correction for internal extinction is at this point not necessary.

PN parent populations
To gain a better understanding of the parent population of the PNe, we estimate the luminosity and the effective temperature of the central stars of the planetary nebula (CSPNs). These parameters are calculated for PNe until m 5007 = 26 with measurable extinction, which corresponds to 87% of the objects within this magnitude limit.
Simulation studies suggest that the maximum conversion efficiency of a central star luminosity into nebular [O iii]λ5007 emission is ∼ 11% (Jacoby 1989;Dopita et al. 1992;Schönberner et al. 2007Schönberner et al. , 2010Gesicki et al. 2018). This occurs under the ideal assumption of optically thick nebula and assumes that [O iii]λ5007 acts as the sole coolant. If the PNe is optically thin, then the efficiency of [O iii]λ5007 production is less, and the luminosity inferred for a PN's central star will be underestimated (Mendez et al. 1992). A high abundance of nitrogen, such as that typically found in Type I PNe (Peimbert & Torres-Peimbert 1983;Phillips 2005), can also increase cooling, and lead to an underestimation of central star luminosity. Moreover, the assumption of lower limit extinction for some cases can also underestimate the luminosity. Therefore, we only consider our luminosity estimates as the lower limits.
To estimate the central stars' effective temperatures, we employed the excitation class method based on the PNe in the LMC (Dopita & Meatheringham 1990;Reid & Parker 2010a). For optically thick PNe, the excitation class temperatures are found to have an empirical correlation with temperature as derived from photo-ionisation modelling (Dopita & Meatheringham 1991;Dopita et al. 1992;Reid & Parker 2010b). To employ this method, we also assume that the metallicity difference between the LMC and NGC 300 is negligible (Bresolin et al. 2009). The revised excitation classes by Reid & Parker (2010b) are defined as with E low employed for low excitation PNe (0 < E < 5) and E high for medium-to high excitation PNe (5 ≤ E < 12). The empirical relation between the excitation class and effective temperature for optically thick PNe is then defined by Reid & Parker (2010b) as Since only the extended mode used in the MUSE-GTO dataset has the wavelength coverage to include the He ii λ4686 line, for uniformity, we determine the excitation class using just the [O iii]λ5007 line and Hβ line. This implies that the effective temperatures for medium-and high excitation PNe with E ≥ 5 (or log T eff ≥ 4.98) are only lower limits. This includes the measurements, which have uncertainties beyond the condition for lowexcitation PNe. Based on 6 PNe in our sample that have He ii λ4686 detection, we estimate that the effective temperatures can be underestimated by 1.5 − 3 times if we only rely on the E low .
On the other hand, if the nebula is optically thin, based on the study in the LMC, the excitation class temperatures can be overestimated by at least 50% compared to the Zanstra temperatures (Villaver et al. 2007;Reid & Parker 2010b). Both luminosity and effective temperature estimates rely on the optical thickness of the PNe. To obtain this, we adopt the criterion of [N ii]λ6584/Hα ≤ 0.3 as the condition for optically thin PNe Reid & Parker 2010b). Since this criterion is not based on nebular modelling in our sample, we use the indication as a more likely condition rather than a definite indicator to explain the possible limitation in our estimations. The estimated stellar parameters are presented in Figure 10. We include the post-AGB tracks from Miller Bertolami (2016) with a stellar metallicity of Z = 0.01, which is the closest to the observed value at the central area of NGC 300 with Z = 0.007 (Kudritzki et al. 2008;Gogarten et al. 2010).
A stellar population study by Jang et al. (2020) using the Hubble Space Telescope found young stars of ∼ 300 Myr, AGB stars with an age between 1 − 3 Gyr and significant number of RGB stars older than 3 Gyr. From a single stellar evolution perspective, the stellar population of NGC 300 can produce a PN central star mass of ∼ 0.7 M from a progenitor mass of 3.0 M , which would have a main sequence lifetime of τ MS > 320 Myr (Miller Bertolami 2016). This implies that, theoretically, central stars within any of the stellar tracks in Figure 10 can be expected. Unfortunately, since most of our luminosities and effective temperatures are lower limits, we are unable to put more constraints on the central star masses at this point.
Within our sample, we also identified several objects as Type I PNe (Peimbert & Torres-Peimbert 1983;Phillips 2005), highly enriched in nitrogen, and classified using [N ii]λ6584/Hα > 1. These objects likely to arise from younger and more massive stars. For a progenitor mass above ∼ 2.5 M , the convective envelope in the thermal pulsing AGB phase is likely to extend to the hydrogen-shell burning layer and produce "hot bottom burning" (HBB). This can dredge up the products of the CNO cycle to the surface, to be later expelled by the stellar wind, therefore increasing the nitrogen-to-oxygen ratio in the nebula (Henry et al. 2018). Observations of PNe in M31 by Fang et al. (2018) put a lower limit of ∼ 2.0 M for HBB. A more thorough analysis, performed for a Type I PN in the M31 young open cluster B477-D075, yields a HBB lower mass limit of ∼ 3.4 M (Davis et al. 2019). This suggests that our approximation for the central star luminosities of Type I PNe is greatly underestimated. For this particular case, we argue that the assumption of [O iii]λ5007 as the only coolant is not true. Since the nitrogen-to-oxygen ratio is high, the nitrogen contribution as additional coolant cannot be neglected (Jacoby 1989), causing the underestimation. This might also explain why we did not see the Type I PNe at our PNLF cut-off, although they are expected to have more massive cores than the typical PNe.
More implications regarding the underlying stellar population can also be inferred from the faint end of the PNLF (Ciardullo 2010). However, the current observational study is still limited to a relatively small sample. Recently, based on a very deep survey in M31, Bhattacharya et al. (2021) found that the steep rise in the number of PN fainter M * + 5 mag is caused by the increased mass fraction of a population older than 5 Gyr. For NGC 300, this implies that the photometry should be complete for m 5007 > 27. Since our PNLF completeness breaks after m 5007 = 27.5, we are unable to provide any insights on this matter at the moment.

Insights on the most luminous PNe
Numerous simulation studies have been conducted to investigate the nature of the PNe at the cut-off of the luminosity function (Jacoby 1989;Schönberner et al. 2007Schönberner et al. , 2010Méndez et al. 2008b;Gesicki et al. 2018;Valenzuela et al. 2019). We review some of them and compare it to our estimated properties to investigate the nature of the most luminous PNe in NGC 300. Using the most recent post-AGB models by Miller Bertolami (2016), simulations of the [O iii]λ5007 fluxes for different progenitor mass have been performed by Gesicki et al. (2018). They found that progenitors with the mass range between 1.5 − 3.0 M are able reach the cutoff absolute magnitude M * = −4.5, assuming that the fluxes at the stellar evolution stages are maximised -also known as maximum nebula hypothesis. It is important to note that the timescale of the 3.0 M track is too short and less likely to be observed. Additionally, they also performed a simulation with an intermediate nebula hypothesis, where the PNe are predominately opaque; this model suggests that the brightest PNe in the luminosity function will have the luminosity log L/L = 3.75 ± 0.13. For comparison, the intrinsically most luminous PN in our sample, PN H9-1, has a lower limit luminosity of log L/L > 3.53. It is also indicated as more likely optically thick, which is in agreement with the simulation. Again, we note that this assumes the ideal 11% maximum efficiency. For example, based on the chemical abundance analysis, the bright PNe in M31 exhibit less conversion efficiency (Jacoby & Ciardullo 1999;Kwitter et al. 2012). Therefore, the actual central star luminosity is likely to be brighter.
Simulation of [O iii]λ5007 flux evolution has also been conducted by Schönberner et al. (2007), who employed 1dimensional radiative-hydrodynamical simulations for the nebulae. They calculated that the most luminous PNe that populate the PNLF cut-off will achieve their maximum luminosity at log T eff = 5.00 K and spend ∼ 500 years in this phase. For PN H9-1, the lower limit temperature is log T eff > 4.97. They also suggest that UV-to [O iii]λ5007 flux conversion process happens most efficiently for central star mass of ∼ 0.62 M , if the nebular shell remains optically thick during the evolution. Referring to the post-AGB models by Miller Bertolami (2016), the initial mass of the progenitor star would be ∼ 2.5 M , and the object would spend less than 1000 years before entering the white dwarf cooling sequence.
Similarly, hydrodynamical models have been used to investigate PNe in nearby galaxies by Schönberner et al. (2010). In these simulations, it was found that central star masses greater than 0.65 M do not exist at the PNLF cut-off. This also supports the result from Gesicki et al. (2018), in that a progenitor mass of 3.0 M for PNe is not expected. This also agrees with the progenitor masses between 2.0 − 2.5 M for the bright PNe in NGC 300 predicted by Stasińska et al. (2013), despite the concerns we mentioned regarding their spectroscopic fluxes in Section 3.3. They derived the progenitor masses using the stellar tracks of Bloecker (1995), which evolve slower than the recent models of Miller Bertolami (2016). This implies the possibility of less massive progenitors if the new stellar tracks are adopted, which is however beyond the scope of our current study.
Recently, the properties of luminous PNe near the PNLF cutoff of M31 have been studied by  for the bulge, and by Galera-Rosillo et al. (2022) for the disk. For the disk, it was found that the four brightest PNe have an average progenitor mass of 1.5 M , which is lower than the values predicted by Schönberner et al. (2007), but still in agreement with Gesicki et al. (2018). Galera-Rosillo et al. (2022) also measure a relatively low average extinction of the PNe with c(Hβ) ∼ 0.1. This means that the PNe originated from an older stellar population, although the disk of M31 also exhibits star forming regions.
In contrast, in the older population of the bulge of M31,  found that the brightest PN have a central star mass > 0.66 M , which means progenitor masses of > 2.5 M . This is found for cases with high extinction, one even reaching c(Hβ) ∼ 0.6, with the average of c(Hβ) ∼ 0.3 for 23 PNe. Current simulations do not predict such massive central stars to be observable, if they exist at all (Schönberner et al. 2010;Gesicki et al. 2018). In old systems, the most luminous PNe are suggested to be products from of blue stragglers -stars that result from a merger during the main sequence (Ciardullo et al. 2005), or symbiotic nebula (Soker 2006). However, both scenarios still do not predict such massive central stars to exist. While the bright Hα background might overestimate the measured extinction, which can lead to the overestimation of luminosity and subsequently the progenitor mass,  in fact did their measurement with an IFU. Their sky subtraction was based on a PSF model, which was claimed to be accurate within 10%.
Lately, Ueta & Otsuka (2021) suggested that the extinction measurement should be solved iteratively, considering the dependency of Hα/Hβ ratio on the electron temperature (T e ) and electron density (n e ). Assuming those two parameter as con-stants would increase the uncertainty of the extinction, and subsequently the stellar parameters. They demonstrate this approach by reanalysing the M31 disk PNe, worked by Galera-Rosillo et al. (2022), and found that the iterative approach yields an average progenitor mass of 2.2 M , instead of 1.5 M for the four brightest PNe (Ueta & Otsuka 2022). While the extinction does not necessarily affect the T e and n e , it may compromise the ionic and elemental abundance analysis (Ueta & Otsuka 2021, 2022. Since we also assume constant T e and n e for our parameters, we are not excluded from this problem. However, as we are missing the diagnostic lines in the blue spectral region and the ones within MUSE wavelength coverage are below the detection limit, we are unable to put constraints on the T e and n e . It would be interesting to repeat the exercise of modelling PN spectra on the basis of improved IFU observations that we believe are superior to slit-based spectroscopy in controlling systematic errors, with the more recent stellar evolution tracks and more careful plasma diagnostics. The future BlueMUSE instrument for the VLT (Richard et al. 2019) will offer the capability with a wavelength coverage down to the atmospheric limit in the UV, which includes the necessary nebular lines for such study.

Conclusions
We analyse 44 fields, obtained with the MUSE instrument to find PNe and construct the PNLF. Using the differential emission line filter (DELF, Roth et al. 2021), we identified more than 500 point sources in [O iii]λ5007, 107 of which were designated as PNe based on spectral classification with the aid of the BPTdiagram (Baldwin et al. 1981). The [O iii]λ5007 magnitudes for the PNe were obtained using DAOPHOT aperture photometry (Stetson 1987) with aperture corrections. With the sample completeness at m 5007 = 27 for most fields, we constructed the PNLF, which exhibits the dip that has been observed in other star forming galaxies (Jacoby & De Marco 2002;Ciardullo 2010;Reid & Parker 2010a). To derive the distance, we employed the maximum likelihood estimation method (Ciardullo et al. 1989) to yield a most likely distance modulus (m − M) 0 = 26.48 +0.11 −0.26 (d = 1.98 +0.10 −0.23 Mpc). For PNe, that are isolated from surrounding emission line sources, and that exhibit bright enough Balmer lines, we measured their extinction. We estimated parameters of the central stars using the extinction corrected fluxes in an attempt to track their origin from the underlying stellar population. We discuss the accuracy of our distance measurement, the effect of local dust for our PNe extinction measurements, and the properties of the most luminous PNe in our sample. The conclusions are as follows: 1. The PNLF distance measurement to NGC 300 is improved with our method and is in excellent agreement with Cepheids and TRGB distances. This is due to the spectral information and spatial resolution of MUSE, that provides a higher PN detection per area, better classification, and accurate photometry. 2. With a limited sample, distance determination based on the minimisation technique is very dependent on the binning. Although coarse binning might provide a better apparent shape of the luminosity function for fitting, it can introduce an unintended systematic shift. Moreover, the details of the PNLF shape, which can provide insights on the stellar population, are also smeared out.
were able to resolve several PNe that are likely obstructed by dust lanes. Any attempt to link the internal extinction and the underlying stellar population requires a quantitative technique to separate the local and internal PNe extinction. 4. We found a few Type I PNe, that evolved from main sequence mass > 2.5 M . Their luminosities are likely underestimated due to the high abundance of nitrogen that serves as a competing coolant with oxygen. They do not populate our PNLF cut-off.
With these results, and other works reported in the literature, we feel encouraged to further develop the IFU observing technique with MUSE to study extragalactic PNe. One of the inherent parameters that we have as yet not utilised is the radial velocity of individual PNe that comes for free as a by-product of the analysis. It will be interesting to find out whether the kinematics can provide hints as to the membership in different populations in NGC 300. Such study was recently done for other disc galaxies: NGC 628 (Aniyan et al. 2018), NGC 6946 (Aniyan et al. 2021), and M31 (Bhattacharya et al. 2019). In the interest of understanding the physical parameters of the PNe, we are currently dependant on the ideal assumption of [O iii]λ5007 maximum conversion and excitation classes to derive the central star parameters, which is not ideal, especially when most cases have no He ii λ4686 coverage. Better constraints on the luminosities and effective temperatures are obtainable through nebular abundance modelling. However, our current wavelength coverage of the MUSE instrument limit us to explore this possibility. Future IFUs, that are optimised in the blue wavelength, such as Blue-MUSE (Richard et al. 2019), will play an important role and allow us to gain more understanding about PNe in the nearby galaxies beyond the Local Group, getting us closer to comprehend the underlying physics behind the constancy of PNLF cutoff across galaxies.