k-Wave User Forum » User Favorites: DLam

k-Wave User Forum » User Favorites: DLam dlam Support for the k-Wave MATLAB toolbox en-US Tue, 12 May 2026 22:59:45 +0000 http://bbpress.org/?v=1.0.2 <![CDATA[Search]]> q http://www.k-wave.org/forum/search.php Bradley Treeby on "Strange increase in starting magnitude when validating attenuation" http://www.k-wave.org/forum/topic/strange-increase-in-starting-magnitude-when-validating-attenuation#post-5805 Wed, 04 Jan 2017 14:39:02 +0000 Bradley Treeby 5805@http://www.k-wave.org/forum/ Hi DLam, Yes, that's correct regarding <code>source.p</code>. If you spectrally decompose your signal, do you get the attenuation you expect? It's possible that using the peak magnitude is not a good measure of attenuation when the medium is dispersive, as the wave will change shape. One way to test this would be to set <code>medium.alpha_mode = no_absorption</code> and just examine the effects of dispersion. For an inhomogeneous medium, you should scale the source based on the local value of the sound speed. Brad DLam on "Strange increase in starting magnitude when validating attenuation" http://www.k-wave.org/forum/topic/strange-increase-in-starting-magnitude-when-validating-attenuation#post-5790 Thu, 15 Dec 2016 16:16:55 +0000 DLam 5790@http://www.k-wave.org/forum/ Brad, Thanks! Sorry for the late reply. Just to clarify, do you mean multiplying source.p by the scaling_factor? I have tried this (0.1 mm grid resolution, f = 1MHz) , and it fixes the simulated source magnitude at the source; but as the wave propagates, the simulated peak magnitude again becomes larger than the theoretical power law. I have removed the BonA parameter, so I am testing this on linear propagation equation. Is this unavoidable? Furthermore, for an inhomogenous medium the medium.sound_speed will vary, thus varying the sound_speed_dispersion with spatial position. Is it possible for the code to account for this as well? Bradley Treeby on "Strange increase in starting magnitude when validating attenuation" http://www.k-wave.org/forum/topic/strange-increase-in-starting-magnitude-when-validating-attenuation#post-5766 Fri, 25 Nov 2016 13:43:42 +0000 Bradley Treeby 5766@http://www.k-wave.org/forum/ Hi DLam, Thanks for the interesting post. After some digging, the reason is due to the dependence of the source scaling factors on sound speed (see Eqs. 2.18 and 2.19 in the k-Wave manual). In an absorbing medium, the sound speed depends on frequency (this is required for causality), and so the source scaling should use the dispersive sound speed. However, k-Wave currently uses the value of <code>medium.sound_speed</code> as shown in the manual. This means the amplitude of the source depends on the absorption values, as you have observed. If you only have a single frequency and the source points all have the same background value, you can correct for this by introducing an additional scaling factor, and multiplying your source by this: <pre><code>sound_speed_dispersive = powerLawKramersKronig(2*pi*source_freq, 0, medium.sound_speed, db2neper(medium.alpha_coeff, medium.alpha_power), medium.alpha_power); scaling_factor = medium.sound_speed ./ sound_speed_dispersive;</code></pre> where it's assumed that <code>medium.sound_speed</code> is a scalar. Alternatively, if sound speed dispersion is not important, you can switch this off by setting: <code>medium.alpha_mode = 'no_dispersion'</code>. Note, a few points about your simulation. The comparison that you're doing is only valid in 1D, i.e., for a plane wave. You can fix this either by moving your simulation to 1D, or by making your source the whole width of the domain with the PML inside, and turning off the PML in the lateral direction by setting: <code>'PMLAlpha', [2, 0]</code> Brad. DLam on "Strange increase in starting magnitude when validating attenuation" http://www.k-wave.org/forum/topic/strange-increase-in-starting-magnitude-when-validating-attenuation#post-5744 Wed, 16 Nov 2016 13:09:18 +0000 DLam 5744@http://www.k-wave.org/forum/ Sorry, I think I have solved my own problem after perusing the 2010 paper on modelling absorption. It is because my alpha_power is too close to 1. Is there a way around this? DLam on "Strange increase in starting magnitude when validating attenuation" http://www.k-wave.org/forum/topic/strange-increase-in-starting-magnitude-when-validating-attenuation#post-5742 Tue, 15 Nov 2016 19:40:49 +0000 DLam 5742@http://www.k-wave.org/forum/ Hello, I am wondering if anyone can give me some help on an odd issue I've been having? I was modelling the ultrasound wave propagating through a homogeneous, linear, attenuating medium in 2D. When I plotted the results, I spotted that the peak of the maximum pressure of the source wave is higher than the source magnitude that I have set by about 10-40%, from attenuation_coeff of 8 to 16 respectively. I am sure that there is no aliasing going on (dx = dy = 0.25e-3, frequency = 0.25 MHz). This is using a continuously-powered sinusoidal line source that I am sure is placed outside of the Perfectly-Matched-Layer. I also moved the line source to the centre of my simulation block, to check that the effect is not caused by minor reflections off the discrete matching layer elements, but this effect is still there. I am wondering if this is an artefact related to the implementation of the fractional Laplacian, or whether I have forgotten something else? Any help is appreciated. This is my code: % k-wave simulation % B. E. Treeby and B. T. Cox, "k-Wave: MATLAB toolbox for the simulation and reconstruction of photoacoustic wave-fields," J. Biomed. Opt., vol. 15, no. 2, p. 021314, 2010. % B. E. Treeby, J. Jaros, A. P. Rendell, and B. T. Cox, "Modeling nonlinear ultrasound propagation in heterogeneous media with power law absorption using a k-space pseudospectral method," J. Acoust. Soc. Am., vol. 131, no. 6, pp. 4324-4336, 2012. % addpath('C:\Users\DLam\Desktop\k-Wave') % in case of non-admin access clear all; close all; filename = 'test_homogenous_medium.mat'; % filename to save to. If don't want to save, use []. mode = 'flat' %'flat' sensor.record = {'p_final', 'p_max', 'p_rms', 'I_avg'}; % this is what to record. If out of memory, delete 'p', 'I', 'I_avg' % create the computational grid Nx = 256; % number of grid points in the x (down-row) direction Ny = 256; % number of grid points in the y (towards right column) direction dx = 0.25*1e-3; % grid point spacing in the x direction [m] dy = 0.25*1e-3; % grid point spacing in the y direction [m] t_end = 60e-6; % end time of sim offset = 22; % [grid points] distance between transducer source and top edge of simulation block %boundary_gp = Nx/2; % [grid points] for flat Tx only. This is the row of the interface. si = [Nx, Ny]; kgrid = makeGrid(Nx, dx, Ny, dy); % CARE: makeGrid implements a 20 grid-point (10 in 3D) thick "border" around the % simulation block that helps with absorbing the incoming wave % Source settings: create pressure distribution source_magnitude = 1e06; % [Pa] source_freq = 0.25e06; % [Hz] %Flat Tx if strcmp(mode,'flat') == 1 starting_gp = 22; source.p_mask = zeros(si); source.p_mask(starting_gp,22:end-22) = 1; source.p_mask = logical(source.p_mask); end % MEDIUM: define the properties of the propagation medium medium.alpha_power = 1.1; % for absorption (exponent). Can only be scalar for k-wave. For tissue. Source: F.A. Duck 1990, "Physical Properties of Tissue". % Soft Tissue medium.sound_speed = 1540*ones(Nx,Ny); % [m/s] tissue medium.density = 1050*ones(Nx,Ny); % [kg/m^3] % tissue medium.alpha_coeff = 8*ones(Nx,Ny);% [dB / (cm.MHz)] fictional number % create time array. This is needed for a continuously-powered transducer. [kgrid.t_array] = makeTime(kgrid, medium.sound_speed, 0.3, t_end); % default CFL criterion is 0.3 source.p = source_magnitude.*sin(2*pi*source_freq*kgrid.t_array); %According to k-wave documentation, this applies "a finite impulse response %(FIR) filter designed using the Kaiser windowing method" to filter away %high harmonics. BUT this also changes the source magnitude. %source.p = filterTimeSeries(kgrid, medium, source.p); % SIMULATION % define the rectangular sensor region by specifying the location of % opposing corners. Leave this the way it is now to get all information % from the sim block. rect1_x_start = 1; rect1_y_start = 1; rect1_x_end = Nx; rect1_y_end = Ny; sensor.mask = [rect1_x_start; rect1_y_start; rect1_x_end; rect1_y_end]; % create a display mask to display the transducer display_mask = source.p_mask; % 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{:}); intensity_grid = sqrt((sensor_data.Ix_avg.^2+sensor_data.Iy_avg.^2)); % test_plot_attenuation input = sensor_data.p_final; MaxP = max(max(input)) sensor_data_mean = mean(input(:,round(Ny/3):round(2*Ny/3),end),2); % reduce effect of reflections from PML at side. y_raw = sensor_data_mean; x = [1:1:Nx].*dx*100; %[cm] alpha_dB = (medium.alpha_coeff(Nx/2,Ny/2)).*((source_freq/1e6).^medium.alpha_power); % [dB/(cm. MHz^power)] alpha_neper = alpha_dB/(8.685889638); % [Np / (cm.MHz^power)] %Theory y = source_magnitude*exp(-alpha_neper*(x-(starting_gp*dx*100))); h1 = plot(x,y_raw); y1 = get(gca,'ylim'); xl = get(gca,'xlim'); hold on h2 = line([starting_gp*dx*100 starting_gp*dx*100],y1,'Color','r'); % boundary h4 = line(xl, [source_magnitude source_magnitude], 'Color','m'); % source magnitude line h3 = plot(x,y,'g'); legend([h1 h2 h3 h4],'Final Pressure as a function of distance', 'Transducer Face Position', 'Theoretical Attenuation','Source Magnitude'); ylabel('Pressure (Pa)') xlabel('Distance (cm)') title('Attenuation of pressure with distance using linear k-wave first order') hold off