k-Wave

1. jtimo
   

   Member
   Joined: Jun '14
   Posts: 4

   Hi,

   I'm using K-wave to run simulations varying different parameters in hopes of optimizing the image reconstruction and I was trying to perform a reconstruction with the source having a pulse frequency instead of just single emission.

   

   This is all I could get. With the kspaceFirstOrder2D I had to delete the p-final part of the sensor structure so that sensor would be just a double instead of a structure, and I'm guessing that's part of the problem. I don't know if there is a way to reconstruct an image with a pulse speed, but if there is any help would be appreciated.

   ux = disc_magnitude*(disc_1 + disc_2);
   
   source.u_mask = zeros(Nx, Ny);
   source.u_mask = disc_magnitude*(disc_1 + disc_2);
   
   % define a time varying sinusoidal velocity source in the x-direction
   source_freq = 0.125e6;       % [Hz]
   source_mag = 2/(medium.sound_speed*medium.density);
   source.ux = -source_mag*sin(2*pi*source_freq*kgrid.t_array);
   
   % filter the source to remove high frequencies not supported by the grid
   source.ux = filterTimeSeries(kgrid, medium, source.ux);
   
   % define a single sensor point
   sensor.mask = zeros(Nx, Ny);
   sensor.mask(Nx/4, Ny/2) = 1;
   
   % define the acoustic parameters to record
   sensor.record = {'p', 'p_final'};

   This is a portion of the code I used. I tried combining the 2D modeling with directivity and the time varying source in a homogeneous medium examples.

   Posted 9 years ago #

2. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,640

   Hi jtimo,

   I'm not sure I follow what you are trying to do. Are you trying to reconstruct a photoacoustic image given a continuous wave (CW) source, rather than an impulsive source? What is it that is being pulsed?

   Brad.

   Posted 9 years ago #

3. jtimo
   

   Member
   Joined: Jun '14
   Posts: 4

   Hi Brad,

   Sorry I should have elaborated more.

   I had made 2 discs using the makeDisc functions and set them to ux. Then based off of the code in the example for the time varying source problems (the homogeneous medium monopole specifically) I tried to set the discs to have multiple pulses, but it made a very odd occurrence as you can see.

   the for the time reversal I had to use only 'p' from the sensor.record because an error occurred when attempting to use the sensor structure saying it only accepted doubles.
   The 'p_final' is probably why it's not working (maybe among other things), but I was just wondering if there was some other way to do it.

   Thanks!

   Posted 9 years ago #

4. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,640

   Hi jtimo,

   From your code snippet, your source is set to produce a continuous wave (CW) sinusoid at a particular frequency, rather than multiple "pulses". Is that what you intended? Incidentally, is there any reason you chose to use a velocity source (which produces a dipole field) rather than a pressure source?

   In the case of a CW source, it's not possible to directly use time reversal reconstruction as in the k-Wave examples, as there is no longer a single t_0 that corresponds to the arrival of the optical pulse.

   Perhaps you can you say a bit more about the experimental setup that you are trying to replicate?

   Brad.

   Posted 9 years ago #

5. jtimo
   

   Member
   Joined: Jun '14
   Posts: 4

   Sorry, as you can tell I'm new to all of this, so thank you so much for your help!

   What I want to set up is have the source emit several pulses to see its effects on reconstruction since we will have the ability to do so with the equipment we have.

   Posted 9 years ago #

6. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,640

   Hi jtimo,

   Here's one example of how you might create the forward data for multiple laser pulses. An important thing to note is that for a given initial pressure distribution, assigning this to source.p0 is the same as assigning this to two successive time points in source.p after scaling by 1/(4*CFL). This is used in this example to create an initial pressure distribution using source.p instead of source.p0. Using the forward data, you can then investigate how this might affect the reconstruction.

   Good luck with your simulations,

   Brad.

   % Example of using k-Wave to generate photoacoustic forward data using
   % multiple laser pulses.
   %
   % author: Bradley Treeby
   % date: 27th June 2014
   
   clear all;
   
   % create the computational grid
   Nx = 128;           % number of grid points in the x (row) direction
   Ny = 128;           % number of grid points in the y (column) direction
   dx = 0.1e-3;        % grid point spacing in the x direction [m]
   dy = 0.1e-3;        % grid point spacing in the y direction [m]
   kgrid = makeGrid(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 time array
   t_end = 20e-6;
   CFL = 0.3;
   kgrid.t_array = makeTime(kgrid, medium.sound_speed, CFL, t_end);
   
   % create initial pressure distribution
   disc_magnitude = 5; % [Pa]
   disc_x_pos = 50;    % [grid points]
   disc_y_pos = 50;    % [grid points]
   disc_radius = 8;    % [grid points]
   p0 = disc_magnitude*makeDisc(Nx, Ny, disc_x_pos, disc_y_pos, disc_radius);
   
   % smooth the initial pressure distribution
   p0 = smooth(kgrid, p0, true);
   
   % reorder p0 to be a 1D vector following MATLAB's linear matrix indexing
   p0 = reshape(p0, Nx*Ny, 1);
   
   % create mask for time varying pressure source
   source.p_mask = ones(Nx, Ny);
   
   % define pulsed source settings
   pulse_sep  = 3e-6;  % [s]
   num_pulses = 3;
   
   % convert to units of time points
   pulse_sep_tp  = round(pulse_sep/kgrid.dt); % [time points]
   
   % create empty time varying pressure source
   source.p = zeros(Nx*Ny, pulse_sep_tp*3 + 1);
   
   % assign initial pressure as pulses
   for j = 1:num_pulses
       source.p(:, 1 + (j-1)*pulse_sep_tp) = p0 / (4*CFL);
       source.p(:, 2 + (j-1)*pulse_sep_tp) = p0 / (4*CFL);
   end
   
   % define a centered circular sensor
   sensor_radius = 4e-3;   % [m]
   num_sensor_points = 50;
   sensor.mask = makeCartCircle(sensor_radius, num_sensor_points);
   
   % run the simulation
   sensor_data = kspaceFirstOrder2D(kgrid, medium, source, sensor);
   
   % plot the simulated sensor data
   figure;
   imagesc(sensor_data, [-1, 1]);
   colormap(getColorMap);
   ylabel('Sensor Position');
   xlabel('Time Step');
   colorbar;

   Posted 9 years ago #

7. jtimo
   

   Member
   Joined: Jun '14
   Posts: 4

   This was amazing! Thank you so much, everything is running smooth now!

   I just have one more quick question if you don't mind.

   Is it possible to set a sensor array of transducers?

   I want to model a specific model (512 transducers in a hemisphere with certain transducer propterties)

   Posted 9 years ago #

8. Bradley Treeby
   

   Developer
   Joined: Mar '10
   Posts: 1,640

   There is lots of information on how to define the sensor in the k-Wave User Manual that should help you get started. If you get stuck on specific issues, let us know.

   Brad.

   Posted 9 years ago #

RSS feed for this topic

Reply

You must log in to post.