k-Wave User Forum » Topic: Transducer simulation with Imager

k-Wave User Forum » Topic: Transducer simulation with Imager http://www.k-wave.org/forum/topic/transducer-simulation-with-imager Support for the k-Wave MATLAB toolbox en-US Wed, 13 May 2026 03:46:09 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "Transducer simulation with Imager" http://www.k-wave.org/forum/topic/transducer-simulation-with-imager#post-4134 Sun, 24 Nov 2013 15:10:28 +0000 Bradley Treeby 4134@http://www.k-wave.org/forum/ Hi nealmp, The line of code giving the error is used in the "Simulating B-mode Ultrasound Images Example" to increase the resolution of a 2D ultrasound image consisting of a series of scan lines. In your case, you only have a single scan line, so the 2D interpolation will generate an error. Brad. nealmp on "Transducer simulation with Imager" http://www.k-wave.org/forum/topic/transducer-simulation-with-imager#post-4131 Fri, 22 Nov 2013 19:00:04 +0000 nealmp 4131@http://www.k-wave.org/forum/ So my objective with the simulation that I pulled together is to map out an image of a simulated environment using acoustic waves. In the code, I have a transducer, a couple of sensors and an object. I'm actually having trouble with a couple of things in the simulation, the first of which being that I'm not getting a reflection off my object. But my main concern involves the Imaging portion of the simulation. The error I'm hit with is displayed after the code, telling me basically that my vectors do not define a grid matching kgrid size. I'm not sure how to go about solving this error. Any help is appreciated. clear all; close all; clc % simulation settings DATA_CAST = 'single'; % ========================================================================= % DEFINE THE K-WAVE GRID % ========================================================================= % set the size of the perfectly matched layer (PML) PML_X_SIZE = 20; % [grid points] PML_Y_SIZE = 10; % [grid points] PML_Z_SIZE = 10; % [grid points] % set total number of grid points not including the PML Nx = 64 - 2*PML_X_SIZE; % [grid points] Ny = 64 - 2*PML_Y_SIZE; % [grid points] Nz = 32 - 2*PML_Z_SIZE; % [grid points] % set desired grid size in the x-direction not including the PML x = 40e-3; % [m] % calculate the spacing between the grid points dx = x/Nx; % [m] dy = dx; % [m] dz = dx; % [m] % create the k-space grid kgrid = makeGrid(Nx, dx, Ny, dy, Nz, dz); % ========================================================================= % DEFINE THE MEDIUM PARAMETERS % ========================================================================= % define the properties of the propagation medium medium.sound_speed = 1540; % [m/s] medium.density = 1000*ones(Nx,Ny,Nz); % [kg/m^3] medium.alpha_coeff = 0.75; % [dB/(MHz^y cm)] medium.alpha_power = 1.5; medium.BonA = 6; object=makeBall(kgrid.Nx,kgrid.Ny,kgrid.Nz,kgrid.Nx/2,kgrid.Ny/2,kgrid.Nz/2,1); medium.sound_speed(object==1) = 150; medium.density(object==1) =720; for i=1:1:length(medium.alpha_coeff) if(i==length(medium.alpha_coeff)) medium.alpha_coeff = medium.alpha_coeff(i); end end for i=1:1:length(medium.sound_speed) if(i==length(medium.sound_speed)) medium.sound_speed = medium.sound_speed(i); end end for i=1:1:length(medium.density) if(i==length(medium.density)) medium.density = medium.density(i); end end % create the time array t_end = 300e-6; % [s] % kgrid.t_array = makeTime(kgrid, medium.sound_speed, [], t_end); [kgrid.t_array, dt] = makeTime(kgrid, medium.sound_speed); kgrid.t_array = 0:dt:t_end; % ========================================================================= % DEFINE THE INPUT SIGNAL % ========================================================================= % define properties of the input signal source_strength = 1e6; % [Pa] tone_burst_freq = 0.5e6; % [Hz] tone_burst_cycles = 5; % create the input signal using toneBurst input_signal = toneBurst(1/kgrid.dt, tone_burst_freq, tone_burst_cycles); % scale the source magnitude by the source_strength divided by the % impedance (the source is assigned to the particle velocity) input_signal = (source_strength./(medium.sound_speed*medium.density)).*input_signal; % ========================================================================= % DEFINE THE ULTRASOUND TRANSDUCER % ========================================================================= % physical properties of the transducer transducer.number_elements = 25; % total number of transducer elements transducer.element_width = 1; % width of each element [grid points] transducer.element_length = 6; % length of each element [grid points] transducer.element_spacing = 0; % spacing (kerf width) between the elements [grid points] transducer.radius = inf; % radius of curvature of the transducer [m] % calculate the width of the transducer in grid points transducer_width = transducer.number_elements*transducer.element_width ... + (transducer.number_elements - 1)*transducer.element_spacing; % use this to position the transducer in the middle of the computational grid transducer.position = round([1, Ny/2 - transducer_width/2, Nz/2 - transducer.element_length/2]); % properties used to derive the beamforming delays transducer.sound_speed = 1540; % sound speed [m/s] transducer.focus_distance = Nx/4; % focus distance [m] transducer.elevation_focus_distance = 19e-3;% focus distance in the elevation plane [m] transducer.steering_angle = 0; % steering angle [degrees] % apodization transducer.transmit_apodization = 'Rectangular'; transducer.receive_apodization = 'Rectangular'; % define the transducer elements that are currently active transducer.active_elements = zeros(transducer.number_elements, 1); transducer.active_elements(6:20) = 1; % append input signal used to drive the transducer transducer.input_signal = input_signal; % create the transducer using the defined settings transducer = makeTransducer(kgrid, transducer); % print out transducer properties transducer.properties; % ========================================================================= % DEFINE SENSOR MASK % ========================================================================= %%%%ORIGINAL CODE%%%%%%% % % create a binary sensor mask with four detection positions % sensor.mask = zeros(Nx, Ny, Nz); % sensor.mask([Nx/4, Nx/2, 3*Nx/4], Ny/2 - transducer_width/2, Nz/2 - transducer.element_length/2) = 1; % create a binary sensor mask with four detection positions sensor.mask = zeros(Nx, Ny, Nz); sensor.mask(1, [5, 40], Nz/2) = 1; %sensor.mask(Nx/2, Ny/2, Nz/2) = 1; % ========================================================================= % RUN % ========================================================================= % set the input settings input_args = {'DisplayMask', transducer.all_elements_mask | sensor.mask | object... 'PMLInside', false, 'PlotPML', false, 'PMLSize', [PML_X_SIZE, PML_Y_SIZE, PML_Z_SIZE], ... 'PMLInside', false, 'DataCast', DATA_CAST, 'PlotScale', [-source_strength/2, source_strength/2]}; % run the simulation [sensor_data] = kspaceFirstOrder3D(kgrid, medium, transducer, sensor, input_args{:}); % calculate the amplitude spectrum of the input signal and the signal % recorded each of the sensor positions [f_input, as_input] = spect([input_signal, zeros(1, 2*length(input_signal))], 1/kgrid.dt); [f, as_1] = spect(sensor_data(1, :), 1/kgrid.dt); [f, as_2] = spect(sensor_data(2, :), 1/kgrid.dt); %[f, as_3] = spect(sensor_data(3, :), 1/kgrid.dt); % ========================================================================= % VISUALISATION % ========================================================================= % plot the input signal and its frequency spectrum figure; subplot(2, 1, 1), plot((0:kgrid.dt:(length(input_signal)-1)*kgrid.dt)*1e6, input_signal, 'k-'); xlabel('Time [\mus]'); ylabel('Input Particle Velocity [m/s]'); subplot(2, 1, 2), plot(f_input/1e6, as_input./max(as_input(:)), 'k-'); hold on; line([tone_burst_freq, tone_burst_freq]/1e6, [0 1], 'Color', 'k', 'LineStyle', '--'); xlabel('Frequency [MHz]'); ylabel('Relative Amplitude Spectrum [au]'); f_max = medium.sound_speed / (2*dx); set(gca, 'XLim', [0 f_max/1e6]); % plot the recorded time series figure; stackedPlot(kgrid.t_array*1e6, {'Sensor Position 1', 'Sensor Position 2', 'Sensor Position 3'}, sensor_data); xlabel('Time [\mus]'); % plot the corresponding amplitude spectrums figure; plot(f/1e6, as_1./max(as_1(:)), 'k-', f/1e6, as_2./max(as_1(:)), 'b-'); legend('Sensor Position 1', 'Sensor Position 2', 'Sensor Position 3'); xlabel('Frequency [MHz]'); ylabel('Normalised Amplitude Spectrum [au]'); f_max = medium.sound_speed / (2*dx); set(gca, 'XLim', [0 f_max/1e6]); % ========================================================================= % IMAGE PROCESSING % ========================================================================= % create radius variable assuming that t0 corresponds to the middle of the % input signal t0 = length(input_signal)*kgrid.dt/2; r = ( (1:length(kgrid.t_array))*kgrid.dt/2 - t0); % [m] % create time gain compensation function based on attenuation value, % transmit frequency, and round trip distance tgc_alpha = 0.25; % [dB/(MHz cm)] tgc = exp(2*tgc_alpha*tone_burst_freq/1e6*r*100); % apply the time gain compensation to each of the scan lines scan_lines = bsxfun(@times, tgc, transducer.scan_line(sensor_data)); % filter the scan lines using both the transmit frequency and the second harmonic scan_lines_fund = gaussianFilter(scan_lines, 1/kgrid.dt, tone_burst_freq, 100, true); scan_lines_harm = gaussianFilter(scan_lines, 1/kgrid.dt, 2*tone_burst_freq, 30, true); % envelope detection scan_lines_fund = envelopeDetection(scan_lines_fund); scan_lines_harm = envelopeDetection(scan_lines_harm); % normalised log compression compression_ratio = 3; scan_lines_fund = logCompression(scan_lines_fund, compression_ratio, true); scan_lines_harm = logCompression(scan_lines_harm, compression_ratio, true); length(scan_lines_fund) length(kgrid) % upsample the image using linear interpolation scale_factor = 2; scan_lines_fund = interp2(1:kgrid.Nt, (1:2).', scan_lines_fund, 1:kgrid.Nt, (1:1/scale_factor:2).'); scan_lines_harm = interp2(1:kgrid.Nt, (1:2).', scan_lines_harm, 1:kgrid.Nt, (1:1/scale_factor:2).'); mat2gray(scan_lines_fund) %end program The Error: Error using griddedInterpolant The grid vectors do not define a grid of points that match the given values. Error in interp2/makegriddedinterp (line 214) F = griddedInterpolant(varargin{:}); Error in interp2 (line 127) F = makegriddedinterp({X, Y}, V, method,extrap); Error in example_us_defining_transducer_test (line 319) scan_lines_fund = interp2(1:kgrid.Nt, (1:2).', scan_lines_fund, 1:kgrid.Nt, (1:1/scale_factor:2).');