%23%20gallery_category%3A%20FDN%20Design%20%26%20Analysis%0A%0Aimport%20marimo%0A%0A__generated_with%20%3D%20%220.23.13%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%0A%20%20%20%20return%20(mo%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%20Frequency-dependent%20pole%20boundaries%0A%0A%20%20%20%20FDN%20with%20frequency-dependent%20absorption%20filters%2C%20but%20*not*%20with%20homogeneous%0A%20%20%20%20(delay-proportional)%20decay.%20Still%2C%20boundaries%20for%20the%20pole%20magnitudes%20can%20be%0A%20%20%20%20computed%20from%20the%20singular%20values%20of%20the%20loop%20transfer%20function%20and%20tested%0A%20%20%20%20against%20the%20actual%20poles.%0A%0A%20%20%20%20The%20loop%20here%20is%20%24P(z)%20%3D%20%5Cmathrm%7Bdiag%7D(z%5E%7Bm%7D)%20-%20A%5C%2C%5Cmathrm%7Bdiag%7D(h(z))%24%20with%0A%20%20%20%20a%20two-tap%20FIR%20absorption%20filter%20%24h(z)%20%3D%200.65%20%2B%200.3%20z%5E%7B-1%7D%24%20on%20every%20delay%0A%20%20%20%20line%20and%20a%20non-orthogonal%20feedback%20matrix%20%24A%20%3D%20Q%2F1.5%24.%0A%0A%20%20%20%20Reference%3A%20*Schlecht%2C%20S.%2C%20Habets%2C%20E.%20(2019).%20Modal%20Decomposition%20of%20Feedback%0A%20%20%20%20Delay%20Networks.%20IEEE%20Transactions%20on%20Signal%20Processing%2067(20)%2C%205340-5351.*%0A%20%20%20%20%5Bdoi%3A10.1109%2Ftsp.2019.2937286%5D(https%3A%2F%2Fdx.doi.org%2F10.1109%2Ftsp.2019.2937286)%0A%0A%20%20%20%20Original%20MATLAB%3A%20%60example_poleBoundaries.m%60%2C%20Sebastian%20J.%20Schlecht%2C%0A%20%20%20%2023%20April%202018.%20Delays%20are%20scaled%20down%20relative%20to%20MATLAB%20so%20the%0A%20%20%20%20root-finding%20stays%20fast.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20from%20types%20import%20SimpleNamespace%0A%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20import%20plotly.graph_objects%20as%20go%0A%20%20%20%20import%20plotly.io%20as%20pio%0A%20%20%20%20import%20torch%0A%0A%20%20%20%20import%20pyFDN%0A%0A%20%20%20%20pio.renderers.default%20%3D%20%22sphinx_gallery%22%0A%20%20%20%20return%20SimpleNamespace%2C%20go%2C%20np%2C%20pyFDN%2C%20torch%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Define%20FDN%20with%20FIR%20absorption%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(SimpleNamespace%2C%20np%2C%20pyFDN)%3A%0A%20%20%20%20np.random.seed(6)%0A%20%20%20%20fs%20%3D%2048000%0A%20%20%20%20num_delays%20%3D%208%0A%20%20%20%20delays%20%3D%20np.random.randint(50%2C%20301%2C%20num_delays)%0A%20%20%20%20input_gain%20%3D%20np.eye(num_delays%2C%201)%0A%20%20%20%20output_gain%20%3D%20np.eye(1%2C%20num_delays)%0A%20%20%20%20direct%20%3D%20np.random.randn(1%2C%201)%0A%0A%20%20%20%20feedback_matrix%20%3D%20pyFDN.random_orthogonal(num_delays)%20%2F%201.5%0A%0A%20%20%20%20%23%20two-tap%20FIR%20absorption%20per%20delay%20line%3A%20h(z)%20%3D%200.65%20%2B%200.3%20z%5E%7B-1%7D%0A%20%20%20%20absorption%20%3D%20SimpleNamespace(%0A%20%20%20%20%20%20%20%20b%3Dnp.zeros((num_delays%2C%201%2C%202))%2C%20a%3Dnp.zeros((num_delays%2C%201%2C%202))%0A%20%20%20%20)%0A%20%20%20%20absorption.a%5B%3A%2C%200%2C%200%5D%20%3D%201.0%0A%20%20%20%20absorption.b%5B%3A%2C%200%2C%200%5D%20%3D%200.65%0A%20%20%20%20absorption.b%5B%3A%2C%200%2C%201%5D%20%3D%200.3%0A%0A%20%20%20%20print(f%22Delays%3A%20%7Bdelays%7D%20(sum%20%3D%20%7Bdelays.sum()%7D)%22)%0A%20%20%20%20return%20(%0A%20%20%20%20%20%20%20%20absorption%2C%0A%20%20%20%20%20%20%20%20delays%2C%0A%20%20%20%20%20%20%20%20direct%2C%0A%20%20%20%20%20%20%20%20feedback_matrix%2C%0A%20%20%20%20%20%20%20%20fs%2C%0A%20%20%20%20%20%20%20%20input_gain%2C%0A%20%20%20%20%20%20%20%20output_gain%2C%0A%20%20%20%20)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Pole%20boundaries%20and%20modal%20decomposition%0A%0A%20%20%20%20%60pole_boundaries%60%20combines%20the%20singular%20values%20of%20the%20feedback%20matrix%20with%0A%20%20%20%20the%20absorption%20magnitude%20responses%20and%20group%20delays.%20For%20the%20poles%2C%20the%0A%20%20%20%20FIR%20absorption%20is%20placed%20as%20an%20SOS%20filter%20behind%20the%20delays%20in%20a%20FLAMO%0A%20%20%20%20model%20(loop%3A%20delay%20%E2%86%92%20absorption%20%E2%86%92%20%24A%24)%20via%20%60dss_to_flamo%60%2C%20and%0A%20%20%20%20%60flamo_to_pr%60%20extracts%20the%20poles%20with%20Ehrlich%E2%80%93Aberth%20refinement%20in%20the%0A%20%20%20%20%24w%20%3D%20z%5E%7B-1%7D%24%20domain.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(%0A%20%20%20%20absorption%2C%0A%20%20%20%20delays%2C%0A%20%20%20%20direct%2C%0A%20%20%20%20feedback_matrix%2C%0A%20%20%20%20fs%2C%0A%20%20%20%20input_gain%2C%0A%20%20%20%20np%2C%0A%20%20%20%20output_gain%2C%0A%20%20%20%20pyFDN%2C%0A%20%20%20%20torch%2C%0A)%3A%0A%20%20%20%20min_curve%2C%20max_curve%2C%20f_bounds%20%3D%20pyFDN.pole_boundaries(%0A%20%20%20%20%20%20%20%20delays%2C%20absorption%2C%20feedback_matrix%5B%3A%2C%20%3A%2C%20None%5D%2C%20fs%0A%20%20%20%20)%0A%0A%20%20%20%20%23%20absorption%20FIR%20as%20one%20SOS%20section%20per%20delay%20line%3A%20%5Bb0%2C%20b1%2C%20b2%2C%20a0%2C%20a1%2C%20a2%5D%0A%20%20%20%20sos_loop%20%3D%20np.zeros((1%2C%206%2C%20delays.size))%0A%20%20%20%20sos_loop%5B0%2C%200%2C%20%3A%5D%20%3D%20absorption.b%5B%3A%2C%200%2C%200%5D%0A%20%20%20%20sos_loop%5B0%2C%201%2C%20%3A%5D%20%3D%20absorption.b%5B%3A%2C%200%2C%201%5D%0A%20%20%20%20sos_loop%5B0%2C%203%2C%20%3A%5D%20%3D%201.0%0A%0A%20%20%20%20model%20%3D%20pyFDN.dss_to_flamo(%0A%20%20%20%20%20%20%20%20A%3Dfeedback_matrix%2C%0A%20%20%20%20%20%20%20%20B%3Dinput_gain%2C%0A%20%20%20%20%20%20%20%20C%3Doutput_gain%2C%0A%20%20%20%20%20%20%20%20D%3Ddirect%2C%0A%20%20%20%20%20%20%20%20m%3Ddelays%2C%0A%20%20%20%20%20%20%20%20Fs%3Dfs%2C%0A%20%20%20%20%20%20%20%20shell%3DFalse%2C%0A%20%20%20%20%20%20%20%20sos_filter%3Dsos_loop%2C%0A%20%20%20%20%20%20%20%20dtype%3Dtorch.float64%2C%0A%20%20%20%20)%0A%20%20%20%20_residues%2C%20poles%2C%20_direct_term%2C%20_is_pair%2C%20_meta%20%3D%20pyFDN.flamo_to_pr(%0A%20%20%20%20%20%20%20%20model%2C%0A%20%20%20%20%20%20%20%20quality_threshold%3D1e-10%2C%0A%20%20%20%20%20%20%20%20refinement_tol%3D1e-10%2C%0A%20%20%20%20%20%20%20%20maximum_iterations%3D80%2C%0A%20%20%20%20%20%20%20%20reject_unstable_poles%3DTrue%2C%0A%20%20%20%20%20%20%20%20deflation_type%3D%22fullDeflation%22%2C%0A%20%20%20%20%20%20%20%20verbose%3DFalse%2C%0A%20%20%20%20)%0A%20%20%20%20print(f%22Number%20of%20FDN%20poles%3A%20%7Bpoles.size%7D%20(conjugate%20pairs%20reduced)%22)%0A%20%20%20%20return%20f_bounds%2C%20max_curve%2C%20min_curve%2C%20poles%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Poles%20between%20the%20boundaries%0A%0A%20%20%20%20Pole%20magnitudes%20converted%20to%20T60%20over%20frequency.%20All%20poles%20lie%20between%20the%0A%20%20%20%20minimum%20and%20maximum%20boundary%20curves.%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(f_bounds%2C%20fs%2C%20go%2C%20max_curve%2C%20min_curve%2C%20np%2C%20poles%2C%20pyFDN)%3A%0A%20%20%20%20pole_freq%20%3D%20pyFDN.rad_to_hertz(np.angle(poles)%2C%20fs)%0A%0A%20%20%20%20fig%20%3D%20go.Figure()%0A%20%20%20%20fig.add_trace(%0A%20%20%20%20%20%20%20%20go.Scatter(%0A%20%20%20%20%20%20%20%20%20%20%20%20x%3Dpole_freq%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y%3DpyFDN.slope_to_rt(pyFDN.lin_to_db(np.abs(poles))%2C%20fs)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20mode%3D%22markers%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20marker%3D%7B%22size%22%3A%204%7D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20name%3D%22Poles%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20)%0A%20%20%20%20for%20_curve%2C%20_name%20in%20%5B%0A%20%20%20%20%20%20%20%20(min_curve%2C%20%22Minimum%20boundary%22)%2C%0A%20%20%20%20%20%20%20%20(max_curve%2C%20%22Maximum%20boundary%22)%2C%0A%20%20%20%20%5D%3A%0A%20%20%20%20%20%20%20%20fig.add_trace(%0A%20%20%20%20%20%20%20%20%20%20%20%20go.Scatter(%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20x%3Df_bounds%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20y%3DpyFDN.slope_to_rt(pyFDN.lin_to_db(_curve)%2C%20fs)%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20mode%3D%22lines%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20line%3D%7B%22width%22%3A%203%7D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20name%3D_name%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20fig.update_layout(%0A%20%20%20%20%20%20%20%20title%3D%22Pole%20T60%20and%20frequency-dependent%20boundaries%22%2C%0A%20%20%20%20%20%20%20%20xaxis%3D%7B%22title%22%3A%20%22Frequency%20(Hz)%22%7D%2C%0A%20%20%20%20%20%20%20%20yaxis%3D%7B%22title%22%3A%20%22Pole%20RT%20(s)%22%7D%2C%0A%20%20%20%20%20%20%20%20template%3D%22plotly_white%22%2C%0A%20%20%20%20%20%20%20%20height%3D460%2C%0A%20%20%20%20)%0A%20%20%20%20fig.show()%0A%20%20%20%20return%20(pole_freq%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22%0A%20%20%20%20%23%23%20Verify%20bounds%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(f_bounds%2C%20max_curve%2C%20min_curve%2C%20np%2C%20pole_freq%2C%20poles%2C%20pyFDN)%3A%0A%20%20%20%20upper_ok%2C%20_%20%3D%20pyFDN.is_bounding_curve(%0A%20%20%20%20%20%20%20%20pole_freq%2C%20np.abs(poles)%2C%20f_bounds%2C%20max_curve%2C%20%22upper%22%0A%20%20%20%20)%0A%20%20%20%20lower_ok%2C%20_%20%3D%20pyFDN.is_bounding_curve(%0A%20%20%20%20%20%20%20%20pole_freq%2C%20np.abs(poles)%2C%20f_bounds%2C%20min_curve%2C%20%22lower%22%0A%20%20%20%20)%0A%20%20%20%20print(f%22All%20poles%20below%20the%20maximum%20boundary%3A%20%7Bbool(upper_ok)%7D%22)%0A%20%20%20%20print(f%22All%20poles%20above%20the%20minimum%20boundary%3A%20%7Bbool(lower_ok)%7D%22)%0A%20%20%20%20assert%20upper_ok%20and%20lower_ok%0A%20%20%20%20return%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A
d500e73c5ba32f35b065793c9905302c