%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%20Time-domain%20FDN%20vs%20FLAMO%20with%20GEQ%20absorption%0A%0A%20%20%20%20The%20same%20FDN%20with%20frequency-dependent%20absorption%20is%20rendered%20by%20two%0A%20%20%20%20independent%20implementations%20and%20the%20impulse%20responses%20are%20compared%3A%0A%0A%20%20%20%201.%20**%60process_fdn%60**%20%E2%80%94%20block%20time-domain%20recursion%3B%20the%20per-delay-line%0A%20%20%20%20%20%20%20SOS%20cascades%20run%20in%20an%20%60SOSFilterBank%60%20and%20the%20FIR%20feedback%20matrix%20in%0A%20%20%20%20%20%20%20an%20%60FIRMatrixFilter%60%2C%20both%20with%20persistent%20state.%0A%20%20%20%202.%20**%60dss_to_flamo%60**%20%E2%80%94%20FLAMO%20frequency-domain%20model%20with%20the%20same%20SOS%0A%20%20%20%20%20%20%20cascades%20as%20%60parallelSOSFilter%60%20and%20the%20FIR%20feedback%20matrix%20as%20a%0A%20%20%20%20%20%20%20%60Filter%60%20module%20in%20the%20loop.%0A%0A%20%20%20%20The%20feedback%20matrix%20is%20a%20paraunitary%20scattering%20matrix%20from%0A%20%20%20%20%60filter_matrix_gallery%60%3B%20the%20absorption%20is%20a%2010-band%20graphic%20EQ%0A%20%20%20%20(%60absorption_geq%60%2C%2011%20biquad%20sections%20per%20delay%20line)%20targeting%20a%0A%20%20%20%20frequency-dependent%20reverberation%20time.%20The%20two%20impulse%20responses%20must%0A%20%20%20%20match%20to%20numerical%20precision.%0A%0A%20%20%20%20Reference%3A%20*Schlecht%2C%20S.%2C%20Habets%2C%20E.%20(2020).%20Accurate%20reverberation%20time%0A%20%20%20%20control%20in%20feedback%20delay%20networks.%20Proc.%20Int.%20Conf.%20Digital%20Audio%20Effects%0A%20%20%20%20(DAFx).*%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%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%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%20FDN%20parameters%20and%20target%20reverberation%20time%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(np%2C%20pyFDN)%3A%0A%20%20%20%20np.random.seed(5)%0A%20%20%20%20fs%20%3D%2048000%0A%20%20%20%20num_delays%20%3D%204%0A%20%20%20%20ir_len%20%3D%20fs%20%20%23%201%20second%0A%0A%20%20%20%20delays%20%3D%20np.sort(np.random.randint(500%2C%202001%2C%20num_delays))%0A%20%20%20%20feedback_matrix%20%3D%20pyFDN.filter_matrix_gallery(%0A%20%20%20%20%20%20%20%20num_delays%2C%20%22Velvet%22%2C%20num_stages%3D3%2C%20sparsity%3D3%0A%20%20%20%20)%0A%20%20%20%20input_gain%20%3D%20np.ones((num_delays%2C%201))%20%2F%20num_delays%0A%20%20%20%20output_gain%20%3D%20np.ones((1%2C%20num_delays))%0A%20%20%20%20direct%20%3D%20np.zeros((1%2C%201))%0A%0A%20%20%20%20%23%20Target%20RT%20at%20the%2010%20GEQ%20bands%20(seconds)%2C%20decaying%20towards%20high%20frequencies%0A%20%20%20%20target_rt%20%3D%20np.array(%5B1.0%2C%200.9%2C%200.8%2C%200.7%2C%200.6%2C%200.5%2C%200.4%2C%200.3%2C%200.25%2C%200.2%5D)%0A%0A%20%20%20%20print(f%22Delays%3A%20%7Bdelays%7D%22)%0A%20%20%20%20print(f%22Feedback%20matrix%3A%20%7Bfeedback_matrix.shape%5B2%5D%7D%20taps%22)%0A%20%20%20%20print(f%22Target%20RT%3A%20%7Btarget_rt%7D%22)%0A%20%20%20%20return%20(%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%20ir_len%2C%0A%20%20%20%20%20%20%20%20num_delays%2C%0A%20%20%20%20%20%20%20%20output_gain%2C%0A%20%20%20%20%20%20%20%20target_rt%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%20Design%20GEQ%20absorption%20filters%0A%0A%20%20%20%20%60absorption_geq%60%20converts%20the%20target%20T60%20to%20a%20per-sample%20dB%20slope%2C%20fits%20a%0A%20%20%20%20graphic%20EQ%2C%20and%20returns%20one%20SOS%20cascade%20per%20delay%20line%2C%20shape%0A%20%20%20%20(N%2C%2011%2C%206).%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(delays%2C%20fs%2C%20go%2C%20np%2C%20num_delays%2C%20pyFDN%2C%20target_rt)%3A%0A%20%20%20%20sos_absorption%20%3D%20pyFDN.absorption_geq(target_rt%2C%20delays%2C%20fs)%0A%20%20%20%20print(f%22Absorption%20SOS%20shape%3A%20%7Bsos_absorption.shape%7D%22)%0A%0A%20%20%20%20%23%20constract%20the%20SOSFilter%0A%20%20%20%20absorption%20%3D%20pyFDN.SOSFilterBank(sos_absorption%2C%20len(delays))%0A%0A%20%20%20%20fig_mag%20%3D%20go.Figure()%0A%20%20%20%20for%20_i%20in%20range(num_delays)%3A%0A%20%20%20%20%20%20%20%20_%2C%20_H_bands%2C%20_W_bands%20%3D%20pyFDN.probe_sos(%0A%20%20%20%20%20%20%20%20%20%20%20%20sos_absorption%5B...%2C%20_i%5D%2C%20np.array(%5B%5D)%2C%20fft_len%3D2**14%2C%20fs%3Dfs%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20_mag_db%20%3D%20pyFDN.lin_to_db(np.abs(np.prod(_H_bands%2C%20axis%3D1)))%0A%20%20%20%20%20%20%20%20fig_mag.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%3D_W_bands%5B%3A%2C%200%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20y%3D_mag_db%20%2F%20delays%5B_i%5D%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%20name%3Df%22delay%3D%7Bdelays%5B_i%5D%7D%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20line%3D%7B%22width%22%3A%201.2%7D%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_mag.update_layout(%0A%20%20%20%20%20%20%20%20title%3D%22Per-delay%20absorption%20magnitude%20(one%20application)%22%2C%0A%20%20%20%20%20%20%20%20xaxis%3D%7B%0A%20%20%20%20%20%20%20%20%20%20%20%20%22title%22%3A%20%22Frequency%20(Hz)%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22type%22%3A%20%22log%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%22range%22%3A%20%5Bnp.log10(50)%2C%20np.log10(fs%20%2F%202)%5D%2C%0A%20%20%20%20%20%20%20%20%7D%2C%0A%20%20%20%20%20%20%20%20yaxis%3D%7B%22title%22%3A%20%22Magnitude%20(dB%2Fsample)%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%3D400%2C%0A%20%20%20%20)%0A%20%20%20%20fig_mag.show()%0A%20%20%20%20return%20absorption%2C%20sos_absorption%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%20Render%20with%20both%20implementations%0A%0A%20%20%20%20%60process_fdn%60%20filters%20the%20delay%20outputs%20block%20by%20block%20(%60SOSFilterBank%60)%0A%20%20%20%20and%20runs%20the%20FIR%20feedback%20matrix%20in%20the%20time-domain%20recursion%0A%20%20%20%20(%60FIRMatrixFilter%60)%3B%20the%20FLAMO%20model%20places%20the%20same%20cascades%20as%20a%0A%20%20%20%20%60parallelSOSFilter%60%20behind%20the%20delays%20and%20the%20FIR%20matrix%20as%20a%20%60Filter%60%0A%20%20%20%20feedback%20module.%20FLAMO%20renders%20circularly%20with%20period%20%60nfft%60%2C%20so%20%60nfft%60%0A%20%20%20%20is%20chosen%20long%20enough%20for%20the%20tail%20to%20decay%20below%20numerical%20precision.%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%20ir_len%2C%0A%20%20%20%20np%2C%0A%20%20%20%20output_gain%2C%0A%20%20%20%20pyFDN%2C%0A%20%20%20%20sos_absorption%2C%0A%20%20%20%20torch%2C%0A)%3A%0A%20%20%20%20impulse%20%3D%20np.zeros(ir_len)%0A%20%20%20%20impulse%5B0%5D%20%3D%201.0%0A%20%20%20%20ir_td%20%3D%20pyFDN.process_fdn(%0A%20%20%20%20%20%20%20%20impulse%2C%0A%20%20%20%20%20%20%20%20delays%2C%0A%20%20%20%20%20%20%20%20feedback_matrix%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%20%20%20%20direct%2C%0A%20%20%20%20%20%20%20%20absorption%3Dabsorption%2C%0A%20%20%20%20)%0A%0A%20%20%20%20model%20%3D%20pyFDN.dss_to_flamo(%0A%20%20%20%20%20%20%20%20feedback_matrix%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%20%20%20%20direct%2C%0A%20%20%20%20%20%20%20%20delays%2C%0A%20%20%20%20%20%20%20%20fs%2C%0A%20%20%20%20%20%20%20%20nfft%3D2**17%2C%0A%20%20%20%20%20%20%20%20sos_filter%3Dsos_absorption%2C%20%20%23%20canonical%20(n_sections%2C%206%2C%20N)%20bank%0A%20%20%20%20%20%20%20%20shell%3DTrue%2C%0A%20%20%20%20%20%20%20%20dtype%3Dtorch.float64%2C%0A%20%20%20%20)%0A%20%20%20%20ir_flamo%20%3D%20pyFDN.flamo_time_response(model).squeeze().astype(np.float64)%5B%3Air_len%5D%0A%0A%20%20%20%20difference%20%3D%20ir_td%20-%20ir_flamo%0A%20%20%20%20max_deviation%20%3D%20np.max(np.abs(difference))%0A%20%20%20%20print(f%22Max%20%7CIR_process%20-%20IR_flamo%7C%20%3D%20%7Bmax_deviation%3A.3e%7D%22)%0A%20%20%20%20assert%20pyFDN.is_almost_zero(difference%2C%20tol%3D1e-9)%0A%20%20%20%20return%20difference%2C%20ir_flamo%2C%20ir_td%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%20Impulse%20responses%0A%0A%20%20%20%20The%20two%20impulse%20responses%20overlap%20to%20numerical%20precision%20(mu-law%20encoded%0A%20%20%20%20for%20visibility%20of%20the%20tail).%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(ir_flamo%2C%20ir_td%2C%20np%2C%20pyFDN)%3A%0A%20%20%20%20t_axis%20%3D%20np.arange(len(ir_td))%0A%20%20%20%20pyFDN.plot_impulse_response(%0A%20%20%20%20%20%20%20%20ir_td%2C%0A%20%20%20%20%20%20%20%20ir_flamo%2C%0A%20%20%20%20%20%20%20%20labels%3D%5B%22process_fdn%20(time%20domain)%22%2C%20%22FLAMO%20(frequency%20domain)%22%5D%2C%0A%20%20%20%20%20%20%20%20title%3D%22Impulse%20response%3A%20process_fdn%20vs%20FLAMO%22%2C%0A%20%20%20%20)%0A%20%20%20%20return%20(t_axis%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%20Error%20over%20time%0A%20%20%20%20%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(difference%2C%20fs%2C%20go%2C%20pyFDN%2C%20t_axis)%3A%0A%20%20%20%20fig_err%20%3D%20go.Figure()%0A%20%20%20%20fig_err.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%3Dt_axis%20%2F%20fs%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20y%3DpyFDN.lin_to_db(difference)%2C%0A%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%20name%3D%22%7CIR_process%20-%20IR_flamo%7C%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20line%3D%7B%22width%22%3A%200.8%7D%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20)%0A%20%20%20%20fig_err.update_layout(%0A%20%20%20%20%20%20%20%20title%3D%22Difference%20between%20the%20two%20implementations%22%2C%0A%20%20%20%20%20%20%20%20xaxis%3D%7B%22title%22%3A%20%22Time%20(s)%22%7D%2C%0A%20%20%20%20%20%20%20%20yaxis%3D%7B%22title%22%3A%20%22Error%20(dB)%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%3D360%2C%0A%20%20%20%20)%0A%20%20%20%20fig_err.show()%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
1f07e2fa66a59934fb34d4663649b994