k-Wave User Forum » Topic: Source Positioning

k-Wave User Forum » Topic: Source Positioning http://www.k-wave.org/forum/topic/source-positioning Support for the k-Wave MATLAB toolbox en-US Wed, 13 May 2026 00:06:49 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-6420 Wed, 18 Apr 2018 09:19:53 +0000 Bradley Treeby 6420@http://www.k-wave.org/forum/ Hi cahyakirana, If you're running simulations in 2D, take a look at the <code>makeMultiArc</code> function, it should do what you need. Brad. cahyakirana on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-6406 Wed, 11 Apr 2018 11:46:37 +0000 cahyakirana 6406@http://www.k-wave.org/forum/ hello, i'm super new to this toolbox, and i need to simulate 20 transducers to see the field pattern based on the transducers' configuration (e.g 5 x 4 arrays or bee-hive configuration). how can i assign the position of each transducer? thank you. % Simulating Transducer Field Patterns Example % % This example demonstrates the use of k-Wave to compute the field pattern % generated by a curved single element transducer in two dimensions. It % builds on the Monopole Point Source In A Homogeneous Propagation Medium % Example. % % author: Bradley Treeby % date: 10th December 2009 % last update: 4th May 2017 % % This function is part of the k-Wave Toolbox (<a href="http://www.k-wave.org" rel="nofollow">http://www.k-wave.org</a>) % Copyright (C) 2009-2017 Bradley Treeby % This file is part of k-Wave. k-Wave is free software: you can % redistribute it and/or modify it under the terms of the GNU Lesser % General Public License as published by the Free Software Foundation, % either version 3 of the License, or (at your option) any later version. % % k-Wave is distributed in the hope that it will be useful, but WITHOUT ANY % WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS % FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for % more details. % % You should have received a copy of the GNU Lesser General Public License % along with k-Wave. If not, see <http://www.gnu.org/licenses/>. clearvars; % ========================================================================= % SIMULATION % ========================================================================= % create the computational grid Nx = 216; % number of grid points in the x (row) direction Ny = 216; % number of grid points in the y (column) direction dx = 50e-3/Nx; % grid point spacing in the x direction [m] dy = dx; % grid point spacing in the y direction [m] kgrid = kWaveGrid(Nx, dx, Ny, dy); % define the properties of the propagation medium medium.sound_speed = 1500; % [m/s] medium.alpha_coeff = 0.75; % [dB/(MHz^y cm)] medium.alpha_power = 1.5; % create the time array kgrid.makeTime(medium.sound_speed); % define a curved transducer element cx = ; cy = ; radius = 10; source.p_mask = makeDisc(Nx, Ny, cx, cy, radius, 'true'); % define a time varying sinusoidal source source_freq = 0.25e6; % [Hz] source_mag = 0.5; % [Pa] source.p = source_mag * sin(2 * pi * source_freq * kgrid.t_array); % filter the source to remove high frequencies not supported by the grid source.p = filterTimeSeries(kgrid, medium, source.p); % create a display mask to display the transducer display_mask = source.p_mask; % create a sensor mask covering the entire computational domain using the % opposing corners of a rectangle sensor.mask = [1, 1, Nx, Ny].'; % set the record mode capture the final wave-field and the statistics at % each sensor point sensor.record = {'p_final', 'p_max', 'p_rms'}; % assign the input options input_args = {'DisplayMask', display_mask, 'PMLInside', false, 'PlotPML', false}; % run the simulation sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:}); % ========================================================================= % VISUALISATION % ========================================================================= % add the source mask onto the recorded wave-field sensor_data.p_final(source.p_mask ~= 0) = 1; sensor_data.p_max(source.p_mask ~= 0) = 1; sensor_data.p_rms(source.p_mask ~= 0) = 1; % plot the final wave-field figure; subplot(1, 3, 1); imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, sensor_data.p_final, [-1 1]); colormap(getColorMap); ylabel('x-position [mm]'); xlabel('y-position [mm]'); axis image; title('Final Wave Field'); % plot the maximum recorded pressure subplot(1, 3, 2); imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, sensor_data.p_max, [-1 1]); colormap(getColorMap); ylabel('x-position [mm]'); xlabel('y-position [mm]'); axis image; title('Maximum Pressure'); % plot the rms recorded pressure subplot(1, 3, 3); imagesc(kgrid.y_vec * 1e3, kgrid.x_vec * 1e3, sensor_data.p_rms, [-1 1]); colormap(getColorMap); ylabel('x-position [mm]'); xlabel('y-position [mm]'); axis image; title('RMS Pressure'); scaleFig(2, 1); Bradley Treeby on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-4891 Fri, 12 Dec 2014 12:02:18 +0000 Bradley Treeby 4891@http://www.k-wave.org/forum/ Hi gondolio, If you don't want to focus the field, just remove the line of code that calls the <code>focus</code> function, and instead assign the <code>input_signal</code> directly to <code>source.p</code>. This will then drive the source aperture with the same signal, without any beamforming delays. The parameters <code>medium.alpha_coeff</code> (a0) and <code>medium.alpha_power</code> (y) describe the power law acoustic attenuation in the medium, where the attenuation is of the form a = a0*f^y. The parameter <code>medium.BonA</code> is the nonlinearity parameter. For water, these values are given by <code>medium.alpha_coeff = 2.17e-3</code>, <code>medium.alpha_power = 2</code>, <code>medium.BonA = 4.96</code>. In water, at lower frequencies and amplitudes, the attenuation and nonlinearity will only have a small effect on the wave field. At higher frequencies and amplitudes, these effects become more important. Brad. gondolio on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-4876 Sat, 06 Dec 2014 03:58:58 +0000 gondolio 4876@http://www.k-wave.org/forum/ Hi Bradley! Thank you. That is what I ended up doing, however, I realized that the source I created is focused which I do not want. Currently, I am trying to create a source with a pulsed signal that has these parameters: <a href="http://imgur.com/SioqhAe" rel="nofollow">http://imgur.com/SioqhAe</a> Which results in a beam pattern that looks like this when the medium is water: <a href="http://imgur.com/3aRLw7s" rel="nofollow">http://imgur.com/3aRLw7s</a> I am attempting to accomplish this using the toneBurst function, but I am not having much luck. '% define properties of the input signal source_strength = 0.076e6; % [Pa] tone_burst_freq = 1.5e6; % [Hz] tone_burst_cycles = 32; % create the input signal using toneBurst input_signal = toneBurst(1/kgrid.dt, tone_burst_freq, tone_burst_cycles,'Plot',true,'Envelope','Rectangular','SignalLength',(1/kgrid.dt)*200e-6); % 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;' Do you have any suggestions? I apologize if this is a stupid question, I am very new to this material. Also, when defining the medium, what do medium.alpha_power and medium.BonA refer to and where could I find these values for water? Bradley Treeby on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-4872 Fri, 05 Dec 2014 17:59:25 +0000 Bradley Treeby 4872@http://www.k-wave.org/forum/ A first glance (without having run your code), it looks like your source and sensor are both defined in the x/y plane (facing in the z-dimension). If you want to see a beam profile, then you could try putting the sensor in the x/z plane, or defining your source to face in a different direction. Hope that helps, Brad. gondolio on "Source Positioning" http://www.k-wave.org/forum/topic/source-positioning#post-4863 Wed, 03 Dec 2014 06:40:58 +0000 gondolio 4863@http://www.k-wave.org/forum/ Hey guys, I'm super new to this toolbox and have barely any idea what I'm doing, so bare with me. I am trying to model how ultrasound waves from a single element,circular transducer propagate through skin, tissue, and bone. I originally tried using the makeTransducer function, however, as I've seen from other posts, there are better ways of modeling this. So I am currently using the method where you create a circular source. I think I am having trouble positioning this source so that it resembles the beam patterns created by transducer as in the example_us_beam_patterns file. This is what I want it to look like: <a href="http://imgur.com/NQxiXax" rel="nofollow">http://imgur.com/NQxiXax</a> This is the view that I am getting: <a href="http://imgur.com/13vNmEM" rel="nofollow">http://imgur.com/13vNmEM</a> Here is my code: <pre><code>% Simulating Ultrasound Beam Patterns Example % % This example shows how the nonlinear beam pattern from an ultrasound % transducer can be modelled. It builds on the Defining An Ultrasound % Transducer and Simulating Transducer Field Patterns examples. % % author: Bradley Treeby % date: 27th July 2011 % last update: 25th September 2012 % % This function is part of the k-Wave Toolbox (<a href="http://www.k-wave.org" rel="nofollow">http://www.k-wave.org</a>) % Copyright (C) 2009-2014 Bradley Treeby and Ben Cox % This file is part of k-Wave. k-Wave is free software: you can % redistribute it and/or modify it under the terms of the GNU Lesser % General Public License as published by the Free Software Foundation, % either version 3 of the License, or (at your option) any later version. % % k-Wave is distributed in the hope that it will be useful, but WITHOUT ANY % WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS % FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for % more details. % % You should have received a copy of the GNU Lesser General Public License % along with k-Wave. If not, see <http://www.gnu.org/licenses/>. clear all; % simulation settings DATA_CAST = 'single'; % set to 'single' or 'gpuArray-single' to speed up computations MASK_PLANE = 'xy'; % set to 'xy' or 'xz' to generate the beam pattern in different planes USE_STATISTICS = true; % set to true to compute the rms or peak beam patterns, set to false to compute the harmonic beam patterns % ========================================================================= % 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 = 128 - 2*PML_X_SIZE; % [grid points] Ny = 256 - 2*PML_Y_SIZE; % [grid points] Nz = 64 - 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; % [kg/m^3] medium.alpha_coeff = 0.75; % [dB/(MHz^y cm)] medium.alpha_power = 1.5; medium.BonA = 6; % create the time array t_end = 45e-6; % [s] kgrid.t_array = makeTime(kgrid, medium.sound_speed, [], 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; % define a circular source element source_radius = 15; % [grid points] source_shape = makeDisc(Nx, Ny, Nx/2, Ny/2, source_radius); % assign source element to source mask source_z_pos = 12; % [grid points] source.p_mask = zeros(Nx, Ny, Nz); source.p_mask(:, :, source_z_pos) = source_shape; % define source input source_freq = 2e6; % [Hz] num_cycles = 5; source_signal = toneBurst(1/kgrid.dt, source_freq, num_cycles); % create time delays to focus the transducer and assign to source focus_position = round([0,0,0]); source.p = focus(kgrid, source_signal, source.p_mask, focus_position, 1450); % ========================================================================= % DEFINE SENSOR MASK % ========================================================================= % define a sensor mask through the central plane sensor.mask = zeros(Nx, Ny, Nz); switch MASK_PLANE case 'xy' % define mask sensor.mask(:, :, Nz/2) = 1; % store y axis properties Nj = Ny; j_vec = kgrid.y_vec; j_label = 'y'; case 'xz' % define mask sensor.mask(:, Ny/2, :) = 1; % store z axis properties Nj = Nz; j_vec = kgrid.z_vec; j_label = 'z'; end % set the record mode such that only the rms and peak values are stored if USE_STATISTICS sensor.record = {'p_rms', 'p_max'}; end % ========================================================================= % RUN THE SIMULATION % ========================================================================= % set the input settings input_args = {'DisplayMask', source.p_mask, ... 'PMLInside', false, 'PlotPML', false, 'PMLSize', [PML_X_SIZE, PML_Y_SIZE, PML_Z_SIZE], ... 'DataCast', DATA_CAST, 'DataRecast', true, 'PlotScale', [-source_strength/2, source_strength/2]}; % stream the data to disk in blocks of 100 if storing the complete time % history if ~USE_STATISTICS input_args = [input_args {'StreamToDisk', 100}]; end % run the simulation sensor_data = kspaceFirstOrder3D(kgrid, medium, source, sensor, input_args{:}); % ========================================================================= % COMPUTE THE BEAM PATTERN USING SIMULATION STATISTICS % ========================================================================= if USE_STATISTICS % reshape the returned rms and max fields to their original position sensor_data.p_rms = reshape(sensor_data.p_rms, [Nx, Nj]); sensor_data.p_max = reshape(sensor_data.p_max, [Nx, Nj]); % plot the beam pattern using the pressure maximum figure; imagesc(j_vec*1e3, (kgrid.x_vec - min(kgrid.x_vec(:)))*1e3, sensor_data.p_max/1e6); xlabel([j_label '-position [mm]']); ylabel('x-position [mm]'); title('Total Beam Pattern Using Maximum Of Recorded Pressure'); colormap(jet(256)); c = colorbar; ylabel(c, 'Pressure [MPa]'); axis image; % plot the beam pattern using the pressure rms figure; imagesc(j_vec*1e3, (kgrid.x_vec - min(kgrid.x_vec(:)))*1e3, sensor_data.p_rms/1e6); xlabel([j_label '-position [mm]']); ylabel('x-position [mm]'); title('Total Beam Pattern Using RMS Of Recorded Pressure'); colormap(jet(256)); c = colorbar; ylabel(c, 'Pressure [MPa]'); axis image; % end the example return end % ========================================================================= % COMPUTE THE BEAM PATTERN FROM THE AMPLITUDE SPECTRUM % ========================================================================= % reshape the sensor data to its original position so that it can be % indexed as sensor_data(x, j, t) sensor_data = reshape(sensor_data, [Nx, Nj, kgrid.Nt]); % compute the amplitude spectrum [freq, amp_spect] = spect(sensor_data, 1/kgrid.dt, 'Dim', 3); % compute the index at which the source frequency and its harmonics occur [f1_value, f1_index] = findClosest(freq, tone_burst_freq); [f2_value, f2_index] = findClosest(freq, tone_burst_freq*2); % extract the amplitude at the source frequency and store beam_pattern_f1 = amp_spect(:, :, f1_index); % extract the amplitude at the second harmonic and store beam_pattern_f2 = amp_spect(:, :, f2_index); % extract the integral of the total amplitude spectrum beam_pattern_total = sum(amp_spect, 3); % plot the beam patterns figure; imagesc(j_vec*1e3, (kgrid.x_vec - min(kgrid.x_vec(:)))*1e3, beam_pattern_f1/1e6); xlabel([j_label '-position [mm]']); ylabel('x-position [mm]'); title('Beam Pattern At Source Fundamental'); colormap(jet(256)); c = colorbar; ylabel(c, 'Pressure [MPa]'); axis image; figure; imagesc(j_vec*1e3, (kgrid.x_vec - min(kgrid.x_vec(:)))*1e3, beam_pattern_f2/1e3); xlabel([j_label '-position [mm]']); ylabel('x-position [mm]'); title('Beam Pattern At Second Harmonic'); colormap(jet(256)); c = colorbar; ylabel(c, 'Pressure [kPa]'); axis image; figure; imagesc(j_vec*1e3, (kgrid.x_vec - min(kgrid.x_vec(:)))*1e3, beam_pattern_total/1e6); xlabel([j_label '-position [mm]']); ylabel('x-position [mm]'); title('Total Beam Pattern Using Integral Of Recorded Pressure'); colormap(jet(256)); c = colorbar; ylabel(c, 'Pressure [MPa]'); axis image; % ========================================================================= % PLOT DIRECTIVITY PATTERN AT FOCUS % ========================================================================= % compute the directivity at each of the harmonics directivity_f1 = squeeze(beam_pattern_f1(round(transducer.focus_distance/dx), :)); directivity_f2 = squeeze(beam_pattern_f2(round(transducer.focus_distance/dx), :)); % normalise directivity_f1 = directivity_f1./max(directivity_f1(:)); directivity_f2 = directivity_f2./max(directivity_f2(:)); % compute relative angles from transducer if strcmp(MASK_PLANE, 'xy') horz_axis = ((1:Ny) - Ny/2)*dy; else horz_axis = ((1:Nz) - Nz/2)*dz; end angles = 180*atan2(horz_axis, transducer.focus_distance)/pi; % plot the directivity figure; plot(angles, directivity_f1, 'k-', angles, directivity_f2, 'k--'); axis tight; set(gca, 'FontSize', 12); xlabel('Angle [deg]'); ylabel('Normalised Amplitude'); legend('Fundamental', 'Second Harmonic', 'Location', 'NorthWest');</code></pre> Sorry if this doesn't make any sense! Thank you.