LDOS and flux monitor resonant frequency disrepancy

I am simulating a circular Bragg grating cavity in cylindrical coordinates with a Au layer underneath the system to act as a backreflector. LDOS function gives me a resonant emission at about 50nm lower wavelength that I can see with flux monitors encompassing the structure. This shift seems to reduce significantly when I am not using the gold layer. I tried different simulation settings such as resolution, simulation time, Courant factor, and none can get rid of this difference entirely. My code is as followed:

import meep as mp
import math
import numpy as np
import matplotlib.pyplot as plt
import argparse
import cmath 
from meep.materials import SiO2, Au

Ac = mp.Medium(index=1.8)
TiO2 = mp.Medium(index=2.67)
Air = mp.Medium(index=1)

l = 0.4     
diam = 4*l   
ff = 0.67   
N = 10      

dpml = 0.5              # thickness of PML
pad_r = 2.0               # padding between structure and edge of PML
pad_z = 2.0
dau = 0.15              # thickness of gold backreflector
ddev = 0.22             # thickness of device layer
dsio2 = 0.44            # thickness of SiO2 spacer
dac = 0.2               # thickness of Ac layer on top

sr = (N+1-ff)*l+diam/2+pad_r+dpml          # cell size in x and y direction
sz = 2*(dpml+pad_z)+ddev +dac+dsio2 +dau   # cell size in z direction

cell_size = mp.Vector3(sr, 0, sz)      # simulation cell size
res = 100                               # resolution

geometry = []

geometry.append(mp.Block(size=mp.Vector3(1e20,0,ddev),            # Add TiO2 device layer
                     center=mp.Vector3(0,0,0),
                     material=TiO2))

geometry.append(mp.Block(size=mp.Vector3(1e20,0,dac),             # Add Ac on top of structure
                     center=mp.Vector3(0,0,ddev/2+dac/2),
                     material=Ac))

geometry.append(mp.Block(size=mp.Vector3(1e20,0,dsio2),           # Add SiO2 spacer layer 
                     center=mp.Vector3(0,0,-(ddev+dsio2)/2),   
                     material=SiO2))

geometry.append(mp.Block(size=mp.Vector3(1e20,0,dau),             # Add Au back reflector
                     center=mp.Vector3(0,0,-(ddev+dau)/2-dsio2),
                     material=Au))

for n in range (0,N+1,+1):                                        # Add N+1 etches (N gratings)
    geometry.append(mp.Block(size=mp.Vector3((1-ff)*l, 0, ddev),center=mp.Vector3(diam/2+(1-ff)*l/2+n*l,0,0), material=Air))  

wvl_min = 0.7           # min wavelength
wvl_max = 0.9            # max wavelength
wvl_cen = 0.5*(wvl_min+wvl_max)
fmin = 1/wvl_max        # min frequency
fmax = 1/wvl_min        # max frequency
fcen = 0.5*(fmin+fmax)  # center frequency
df = fmax-fmin          # frequency width
src = mp.Source(mp.GaussianSource(fcen, fwidth=df), component=mp.Er, center=mp.Vector3(0,0,ddev/2+0.15))  #Add dipole-like source

sim = mp.Simulation(cell_size= cell_size,
                    sources=[src],
                    resolution=res, 
                    geometry=geometry,
                    boundary_layers=[mp.PML(dpml)],
                    dimensions = mp.CYLINDRICAL,
                    m=+1,
                    Courant = 0.5,
                    progress_interval=60,
                    eps_averaging = False,
                    force_complex_fields=True)

nfreq = 100
box_dis = 0.2 

box_z1 = sim.add_flux(fcen, df, nfreq, mp.FluxRegion(center=mp.Vector3(0,0,ddev/2+0.15-box_dis),size=mp.Vector3(2*box_dis,0,0)))
box_z2 = sim.add_flux(fcen, df, nfreq, mp.FluxRegion(center=mp.Vector3(0,0,ddev/2+0.15+box_dis),size=mp.Vector3(2*box_dis,0,0)))
box_r = sim.add_flux(fcen, df, nfreq, mp.FluxRegion(center=mp.Vector3(box_dis,0,ddev/2+0.15),size=mp.Vector3(z=2*box_dis)))

time_after_sources = 150

harminv_instance = mp.Harminv(mp.Er, mp.Vector3(0,0,ddev/2+0.15), fcen, df)
ldos_instance = mp.Ldos(fcen, df, nfreq)

sim.run(mp.after_sources(harminv_instance), mp.dft_ldos(ldos=ldos_instance),until_after_sources=time_after_sources)

ldos_results = np.transpose(
                    np.array([mp.get_ldos_freqs(ldos_instance), sim.ldos_data]))

plt.figure(dpi=100)
plt.plot(1 / ldos_results[:, 0], ldos_results[:, 1], 'b-')
plt.show()

flux_freqs = np.array(mp.get_flux_freqs(box_z2))
flux_up = np.array(mp.get_fluxes(box_z2))
flux_bot = np.array(mp.get_fluxes(box_z1))
flux_side = np.array(mp.get_fluxes(box_r))
flux_total= -flux_bot+flux_up-flux_side
flux_upwards= np.array(mp.get_fluxes(box_z2))

flux_wvl=1/flux_freqs
plt.figure(dpi=100)
plt.plot(flux_wvl, flux_total, 'r-', label='Total emission')

plt.plot(flux_wvl, flux_upwards, 'b-',label='Upward emission')
plt.legend(loc='upper right')
plt.xlabel('Wavelength (µm)')
plt.ylabel('Arbitrary intensity')

q_results = []
for mode in harminv_instance.modes:
    q_results.append([1000/mode.freq, mode.decay, mode.Q, abs(mode.amp)])

q_results = np.array(q_results)

Can the metal layer causing such problems?

NanoComp / meep

LDOS and flux monitor resonant frequency disrepancy #1327