http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm Support for the k-Wave MATLAB toolbox en-US Tue, 12 May 2026 23:07:20 +0000 http://bbpress.org/?v=1.0.2 q http://www.k-wave.org/forum/search.php http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8896 Thu, 17 Aug 2023 02:52:35 +0000 zengqw2021 8896@http://www.k-wave.org/forum/ <p>Dear anyone,<br /> My sensor has a center frequency and finite bandwidth. So, should we need to add the center frequency and bandwidth of the transducer in the forward simulation as follows:</p> <p>sensor.frequency_response = [2e6 70] % center frequency and percentage bandwidth<br /> % run the simulation to get the PA signal<br /> sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:}); </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8855 Fri, 30 Jun 2023 15:45:56 +0000 zengqw2021 8855@http://www.k-wave.org/forum/ <p>Dear Bradley Treeby,<br /> As you said, I'm trying to replicate the ultrasonic signal reception process with sensors of finite bandwidth. So, should I consider the center frequency and bandwidth of the transducer in the forward simulation rather than in the back projection reconstruction？If so, how could I set the transducer frequency and bandwidth?<br /> I found a possibly relevant attribute in the kspaceFirstOrder2D function named"sensor.frequency_response", which says 'two element array specifying the center frequency and percentage bandwidth of a frequency domain Gaussian filter applied to the sensor_data'. Therefore, I think I can simulate the process of receiving ultrasonic signal by ultrasonic transducer with specific center frequency and bandwidth as follows：</p> <p> 。。。<br /> sensor.frequency_response = [2e6 70] % center frequency and percentage bandwidth<br /> % run the simulation to get the PA signal<br /> sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor, input_args{:});</p> <p> In addition, If I take into account the center frequency and bandwidth of the transducer in the forward simulation, do I also need to take into account the back projection process? </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8803 Wed, 07 Jun 2023 19:17:01 +0000 Bradley Treeby 8803@http://www.k-wave.org/forum/ <p>It depends entirely on what processes you’re trying to replicate. For example, if your sensor has a finite bandwidth, you might want to include it in your simulation. </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8799 Wed, 07 Jun 2023 08:52:35 +0000 zengqw2021 8799@http://www.k-wave.org/forum/ <p>Dear,<br /> I found some differences between compensating for acoustic attenuation before or after filtering.<br /> So, which way should I choose.</p> <p>*************************************************************</p> <p>Comp_data = attenComp(pa_data, Ts, speed, alpha_0, y0);</p> <p>% filter the sensor data using a Gaussian filter<br /> center_freq = fc; % [Hz]<br /> bandwidth = 70; % [%]<br /> sensor_data_gaussian1 = gaussianFilter(pa_data, freqs, center_freq, bandwidth);<br /> sensor_data_gaussian2 = gaussianFilter(Comp_data, freqs, center_freq, bandwidth);</p> <p>% add noise to the recorded sensor data<br /> signal_to_noise_ratio = 40; % [dB]<br /> sensor_data1 = addNoise(sensor_data_gaussian1, signal_to_noise_ratio, 'peak');<br /> sensor_data2 = addNoise(sensor_data_gaussian2, signal_to_noise_ratio, 'peak');</p> <p>%----Attenuation compensation after filtering<br /> Comp_data2 = attenComp(sensor_data1, Ts, speed, alpha_0, y0); </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8798 Wed, 07 Jun 2023 08:52:33 +0000 zengqw2021 8798@http://www.k-wave.org/forum/ <br /> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8777 Mon, 29 May 2023 17:07:35 +0000 zengqw2021 8777@http://www.k-wave.org/forum/ <p>Dear Bradley Treeby,<br /> The back-projection I coded up is listed as follows:</p> <p>function Rec_img = UBP_2D(kgrid,medium,sensor,sensor_data)</p> <p>% INPUTS:<br /> % kgrid,<br /> % kgrid.Nx、kgrid.Ny,kgrid.dx、kgrid.dy,kgrid.x,kgrid.y<br /> % kgrid.t_array,kgrid.dt<br /> % medium,<br /> % medium.speed<br /> % sensor,<br /> % sensor.mask<br /> % sensor.scan_radius<br /> % sensor_data，<br /> % Nd×Nt，<br /> %<br /> % Xu M, Wang LV. Universal back-projection algorithm<br /> %-----------doi:10.1103/PhysRevE.71.016706--------------</p> <p>NumElem = sensor.Nd ; %<br /> sensor_radius = sensor.scan_radius ; % </p> <p>Nx = kgrid.Nx;<br /> Ny = kgrid.Ny;</p> <p>% dx = kgrid.dx;<br /> % dy = kgrid.dy;</p> <p>SpeedofSound = medium.sound_speed; % [m/s]</p> <p>t = kgrid.t_array;<br /> dt = kgrid.dt;</p> <p>sensor_location_x = sensor.mask(1,:); % 2*len [x1,y1;x2,y2.....]<br /> sensor_location_y = sensor. Mask(2,:);</p> <p>p0_recon1 = zeros(Nx,Ny);<br /> p0_recon_location_x = kgrid.x;<br /> p0_recon_location_y = kgrid.y;</p> <p>% arc = 2 * pi * sensor_radius / NumElem;</p> <p>arc = sensor.trans_diameter/sensor.scan_radius;</p> <p>sum_w_omega = 0;<br /> for i = 1:NumElem<br /> distance_x = p0_recon_location_x - sensor_location_x(i);<br /> distance_y = p0_recon_location_y - sensor_location_y(i);<br /> distance_xy = sqrt(distance_y.^2 + distance_x.^2);</p> <p> distance_xy_time = distance_xy./SpeedofSound;</p> <p> distance_xy_index = floor(distance_xy_time ./dt);<br /> distance_xy_index(distance_xy_index&lt;1) = 1;</p> <p> unit_normal_vector = -[sensor_location_x(i),sensor_location_y(i)]./ norm([sensor_location_x(i),sensor_location_y(i)]);<br /> unit_dp_vector_x = distance_x./distance_xy;<br /> unit_dp_vector_y = distance_y./distance_xy;<br /> w_omega = arc ./ distance_xy.^2 .* (unit_dp_vector_x .* sensor_location_x(i) + ...<br /> unit_dp_vector_y .* sensor_location_y(i)) ./ -norm([sensor_location_x(i),sensor_location_y(i)]);</p> <p> p_i = sensor_data(i,:);</p> <p> %% b(r,t) = 2*p(r,t) - 2*t*dp(r,t)/d(t)<br /> bp_i = 2*p_i - 2*t .* diff([0,p_i])./dt;</p> <p> p0_i = w_omega .* bp_i(distance_xy_index);</p> <p> sum_w_omega = sum_w_omega + w_omega;</p> <p> p0_recon1 = p0_recon1 + p0_i;<br /> % break<br /> end<br /> p0_recon12 = p0_recon1 ./sum_w_omega;</p> <p>Rec_img = p0_recon12; </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8773 Tue, 23 May 2023 06:27:09 +0000 Bradley Treeby 8773@http://www.k-wave.org/forum/ <p>Have you looked at <a href="http://www.k-wave.org/documentation/attenComp.php"><code>attencomp</code></a>? I'm not sure which back-projection you've coded up, but might be worth checking if there is a scaling factor also. </p> http://www.k-wave.org/forum/topic/attenuation-compensation-in-the-back-projection-algorithm#post-8772 Mon, 22 May 2023 13:50:53 +0000 zengqw2021 8772@http://www.k-wave.org/forum/ <p>Dear all,<br /> In the forward simulation, I have set the medium properties as follows and an initial pressure distribution is created. The received ultrasound signal is used to reconstruct the image by the back projection algorithm. However, I found that the reconstructed sound pressure was much smaller than the initial value. Can I make some attenuation compensation to make up for this difference？I hope someone can give me some help.<br /> Best.</p> <p>% define the properties of the propagation medium<br /> medium.density = 1.29 ; % [kg/m^3]<br /> medium.sound_speed = 330 ; % [m/s]<br /> medium.alpha_coeff = 1.59; % [dB/(MHz^y cm)]<br /> medium.alpha_power = 2.0;</p> <p>%%=============================================<br /> % create initial pressure distribution<br /> Cv = 0.72 ; % Specific heat coefficient of constant volume<br /> Heat_beta = 1/273.15; % the coefficient of the thermal expansion[1/K]</p> <p>%----Construct a Gaussian pressure source<br /> [X,Y] = meshgrid(kgrid.x_vec,kgrid.y_vec);</p> <p>disc_H1 = 30; % [kJ/m3]<br /> disc_x_pos1 = 0; % actual length [m]<br /> disc_y_pos1 = 0; % actual length [m]<br /> FWHM = 120e-6; %<br /> disc_radius1 = FWHM*sqrt(log(2))/2; % actual length [m]<br /> Gaus_1 = exp(-((X-disc_x_pos1).^2+(Y-disc_y_pos1).^2)/disc_radius1^2);<br /> disc_1 = disc_H1*Heat_beta*medium.sound_speed^2/Cv*Gaus_1; </p>