simus fails for subaperture transmit (please review the fix mentioned)

[hyper-trophy]

I am trying simulating focused pressure fields using a sub-aperture of a linear array, simus fails when txdelay has nan entries.
Here's the code to recreate the issue

import pymust, tqdm
import numpy as np, matplotlib.pyplot as plt

# grid of reflectors
nx, nz = (7, 10)
x = np.linspace(-32e-3/2, 32e-3/2, nx)
z = np.linspace(1e-3, 4e-2, nz)
xv, zv = np.meshgrid(x, z)
scatterer_coords = np.stack([xv, zv], axis=-1).reshape(-1,2)
xs, zs = scatterer_coords[:,0], scatterer_coords[:,1]
RC = np.ones_like(xs)*np.random.uniform(0,1,xs.shape)
plt.scatter(xs, zs, c=RC)

param = pymust.getparam('L12-3V')

# delay profile, using nan as delay for inactive elements
param_subaperture = param.copy()
param_subaperture['Nelements'] = 32
xf, zf = 0, 2e-2 # focus wrt subaperture
txdelay_subaperture = pymust.txdelay(xf, zf, param_subaperture)
txdelay = np.empty(param['Nelements'])
txdelay.fill(np.nan)
txdelay[0:32] = np.squeeze(txdelay_subaperture)
plt.bar(np.linspace(1, param['Nelements'], param['Nelements']), txdelay)

# pfield works with txdelays having nan entries
nx, nz = (100, 100)
x_extent = np.linspace(-2e-2,2e-2,nx)
y_extent = np.linspace(0,4e-2,nz)
xi,zi = np.meshgrid(x_extent, y_extent);
yi = np.zeros(xi.shape)
pressure_field, _, _ = pymust.pfield(xi,yi,zi,txdelay[None],param);

plt.imshow(20 * np.log10(pressure_field/ np.max(pressure_field)), cmap='hot', extent=[xi[0,0], xi[0, -1], zi[-1,-1], zi[0, 0]])
plt.colorbar()
plt.clim(-30, 0)

# simus fails when txdelay has nan entries. 
receives, _ = pymust.simus(xs,zs,RC,txdelay[None], param)

Error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[48], line 1
----> 1 receives, _ = pymust.simus(xs,zs,RC,txdelay[None], param)

File [~/miniconda3/envs/pymust/lib/python3.13/site-packages/pymust/simus.py:358](~/miniconda3/envs/pymust/lib/python3.13/site-packages/pymust/simus.py#line=357), in simus(*varargin)
    356 df = 1/2/(2*maxD/param.c + np.max(delaysTX.flatten() + param.RXdelay.flatten())) # % to avoid aliasing in the time domain
    357 df = df*options.FrequencyStep
--> 358 Nf = 2*int(np.ceil(param.fc/df))+1 # % number of frequency samples
    359 #%-- Run PFIELD to calculate the RF spectra
    360 RFspectrum = np.zeros((Nf,NumberOfElements), dtype = np.complex64)# % will contain the RF spectra

ValueError: cannot convert float NaN to integer

This change seems to fix the issue

    valid_tx_delays = np.array([e for e in delaysTX.flatten() if not np.isnan(e)])
    df = 1/2/(2*maxD/param.c + np.max(np.concat((valid_tx_delays,param.RXdelay.flatten())))) # % to avoid aliasing in the time domain
    # df = 1/2/(2*maxD/param.c + np.max(delaysTX.flatten() + param.RXdelay.flatten())) # % to avoid aliasing in the time domain
    df = df*options.FrequencyStep
    Nf = 2*int(np.ceil(param.fc/df))+1 # % number of frequency samples

[gbernardino]

Closed with commit #48