commit 98d76e06f5c2161d305808ac75b0ca38ef8526c3
parent a33a3e9722d8156792bb4cdd38bdd0522455f847
Author: Egor Achkasov <eaachkasov@edu.hse.ru>
Date: Thu, 16 Jan 2025 10:35:49 +0100
Rm sampling_freq, num_samples; Fix LoS and scattering bug
Diffstat:
| M | compute_paths.c | | | 165 | +++++++++++++++++++++++++++++++++++-------------------------------------------- |
| M | compute_paths.h | | | 66 | +++++++++++++++++++++++++++++++++++++++--------------------------- |
| M | compute_paths_pybind11.cpp | | | 77 | ++++++++++++++++++++++++++++++++++++++++++++++++----------------------------- |
| M | test.c | | | 55 | +++++++++++++++++++++++++++++++------------------------ |
| M | test.py | | | 42 | +++++++++++++++++++++++++++++------------- |
5 files changed, 219 insertions(+), 186 deletions(-)
diff --git a/compute_paths.c b/compute_paths.c
@@ -378,29 +378,33 @@ void refl_coefs(
/* ==== MAIN FUNCTION ==== */
void compute_paths(
- IN const char *mesh_filepath, /* path to the mesh file */
- IN const float *rx_positions, /* shape (num_rx, 3) */
- IN const float *tx_positions, /* shape (num_tx, 3) */
- IN const float *rx_velocities, /* shape (num_rx, 3) */
- IN const float *tx_velocities, /* shape (num_tx, 3) */
- IN float carrier_frequency, /* > 0.0 (IN GHz!) */
- IN float sampling_frequency, /* > 0.0 (IN Hz!) */
- IN size_t num_rx, /* number of receivers */
- IN size_t num_tx, /* number of transmitters */
- IN size_t num_paths, /* number of paths */
- IN size_t num_bounces, /* number of bounces */
- IN size_t num_samples, /* number of samples */
- OUT size_t *num_paths_out, /* number of paths computed */
- OUT float *a_te_re, /* output array real parts of TE gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_te_im, /* output array imaginary parts of TE gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_tm_re, /* output array real parts of TM gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_tm_im, /* output array imaginary parts of TM gains (num_rx, num_tx, num_paths_out) */
- OUT float *tau /* output array of delays (num_rx, num_tx, num_paths_out) */
+ IN const char *mesh_filepath, /* path to the mesh file */
+ IN const float *rx_positions, /* shape (num_rx, 3) */
+ IN const float *tx_positions, /* shape (num_tx, 3) */
+ IN const float *rx_velocities, /* shape (num_rx, 3) */
+ IN const float *tx_velocities, /* shape (num_tx, 3) */
+ IN float carrier_frequency, /* > 0.0 (IN GHz!) */
+ IN size_t num_rx, /* number of receivers */
+ IN size_t num_tx, /* number of transmitters */
+ IN size_t num_paths, /* number of paths */
+ IN size_t num_bounces, /* number of bounces */
+ /* LoS */
+ OUT float *a_te_re_los, /* output array real parts of TE gains (num_rx, num_tx) */
+ OUT float *a_te_im_los, /* output array imaginary parts of TE gains (num_rx, num_tx) */
+ OUT float *a_tm_re_los, /* output array real parts of TM gains (num_rx, num_tx) */
+ OUT float *a_tm_im_los, /* output array imaginary parts of TM gains (num_rx, num_tx) */
+ OUT float *tau_los, /* output array of delays (num_rx, num_tx) */
+ /* Scatter */
+ OUT float *a_te_re_scat, /* output array real parts of TE gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_te_im_scat, /* output array imaginary parts of TE gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_tm_re_scat, /* output array real parts of TM gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_tm_im_scat, /* output array imaginary parts of TM gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *tau_scat /* output array of delays (num_bounces, num_rx, num_tx, num_paths) */
)
{
/* Init loop indices and offset variables */
size_t i, j;
- size_t off;
+ size_t off, off_scat;
/* Load the scene */
Mesh mesh = load_mesh_ply(mesh_filepath, carrier_frequency);
@@ -435,41 +439,30 @@ void compute_paths(
free(ray_directions);
/* Initialize variables needed for the RT */
- /* Number of output paths is unknown at this point */
- *num_paths_out = 0;
- /* As num_paths_out is not known in advance,
- * we will realloc the output arrays during propagation.
- * We will start with a capacity of num_paths_out_realloc_cap
- * and reallocate the arrays when needed.
- */
- const size_t num_paths_out_realloc_cap = 1024;
- size_t num_cur_allocs = 1;
- /* Init the output gains and delays */
- /* shape (num_rx, num_tx, num_paths_out) */
- a_te_re = (float*)malloc(num_rx * num_tx * num_paths_out_realloc_cap * sizeof(float));
- a_te_im = (float*)malloc(num_rx * num_tx * num_paths_out_realloc_cap * sizeof(float));
- a_tm_re = (float*)malloc(num_rx * num_tx * num_paths_out_realloc_cap * sizeof(float));
- a_tm_im = (float*)malloc(num_rx * num_tx * num_paths_out_realloc_cap * sizeof(float));
- tau = (float*)malloc(num_rx * num_tx * num_paths_out_realloc_cap * sizeof(float));
- /* Bouncing rays gains and delays */
+
+ /* Bouncing (reflecting) rays gains and delays */
/* shape (num_tx, num_paths) */
- float *a_te_re_t = (float*)malloc(num_tx * num_paths * sizeof(float));
- float *a_te_im_t = (float*)calloc(num_tx * num_paths, sizeof(float));
- float *a_tm_re_t = (float*)malloc(num_tx * num_paths * sizeof(float));
- float *a_tm_im_t = (float*)calloc(num_tx * num_paths, sizeof(float));
+ float *a_te_re_refl = (float*)malloc(num_tx * num_paths * sizeof(float));
+ float *a_te_im_refl = (float*)calloc(num_tx * num_paths, sizeof(float));
+ float *a_tm_re_refl = (float*)malloc(num_tx * num_paths * sizeof(float));
+ float *a_tm_im_refl = (float*)calloc(num_tx * num_paths, sizeof(float));
for (size_t i = 0; i < num_tx * num_paths; ++i)
- a_te_re_t[i] = a_tm_re_t[i] = 1.f;
+ a_te_re_refl[i] = a_tm_re_refl[i] = 1.f;
float *tau_t = (float*)calloc(num_tx * num_paths, sizeof(float));
+
/* Active rays bitmask. 1 if the ray is still traced, 0 if it has left the scene */
uint8_t *active = (uint8_t*)malloc((num_tx * num_paths / 8 + 1) * sizeof(uint8_t));
for (size_t i = 0; i < num_tx * num_paths / 8 + 1; ++i)
active[i] = 0xff;
+ size_t num_active = num_tx * num_paths;
+
/* Temporary variables */
float t, r_te_re, r_te_im, r_tm_re, r_tm_im;
float a_te_re_new, a_te_im_new, a_tm_re_new, a_tm_im_new;
int32_t ind;
Ray r, *r_ptr;
Vec3 *h = (Vec3*)malloc(sizeof(Vec3));
+
/* Friis free space loss multiplier.
Multiply this by a distance and take the second power to get a free space loss. */
float free_space_loss_multiplier = 2.f * PI * carrier_frequency * 1e9 / SPEED_OF_LIGHT;
@@ -478,40 +471,35 @@ void compute_paths(
/* Calculate LoS paths directly from tx to rx */
for (i = 0; i != num_rx; ++i)
for (j = 0; j != num_tx; ++j) {
+ off = i * num_tx + j;
t = -1.f;
ind = -1;
r.o = (Vec3){tx_positions[j * 3], tx_positions[j * 3 + 1], tx_positions[j * 3 + 2]};
r.d = (Vec3){rx_positions[i * 3], rx_positions[i * 3 + 1], rx_positions[i * 3 + 2]};
r.d = vec3_sub(&r.d, &r.o);
moeller_trumbore(&r, &mesh, &t, &ind, &theta);
- if (ind == -1) continue;
+ if (ind != -1 && t <= 1.f) {
+ /* An obstacle between the tx and the rx has been hit */
+ a_te_re_los[off] = a_te_im_los[off] = a_tm_re_los[off] = a_tm_im_los[off] = tau_los[off] = 0.f;
+ continue;
+ }
/* No triangle has been hit, the path is LoS */
- /* Calculate teh distance */
+ /* Calculate the distance */
t = sqrtf(vec3_dot(&r.d, &r.d));
- /* Calculate the delay */
- tau[*num_paths_out] = t / SPEED_OF_LIGHT;
/* Calculate the free space loss */
free_space_loss = free_space_loss_multiplier * t;
free_space_loss = free_space_loss * free_space_loss;
/* Set the gains for this LoS path */
- a_te_re[*num_paths_out] = a_tm_re[*num_paths_out] = 1.f / free_space_loss;
- a_te_im[*num_paths_out] = a_tm_im[*num_paths_out] = 0.f;
- /* Increment the number of paths */
- ++(*num_paths_out);
- /* Realloc the output arrays if needed */
- if (*num_paths_out % num_paths_out_realloc_cap == 0) {
- ++num_cur_allocs;
- a_te_re = (float*)realloc(a_te_re, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_te_im = (float*)realloc(a_te_im, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_tm_re = (float*)realloc(a_tm_re, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_tm_im = (float*)realloc(a_tm_im, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- tau = (float*)realloc(tau, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- }
+ a_te_re_los[off] = a_tm_re_los[off] = 1.f / free_space_loss;
+ a_te_im_los[off] = a_tm_im_los[off] = 0.f;
+ /* Calculate the delay */
+ tau_los[off] = t / SPEED_OF_LIGHT;
}
/* Bounce the rays. On each bounce:
* - Add path length / SPEED_OF_LIGHT to tau
* - Update gains using eqs. (31a)-(31b) from ITU-R P.2040-3
+ * - Spawn refraction rays with gains as per eqs. (31c)-(31d) from ITU-R P.2040-3
*/
for (i = 0; i < num_bounces; ++i) {
for (off = 0; off != num_tx * num_paths; ++off) {
@@ -527,9 +515,7 @@ void compute_paths(
moeller_trumbore(r_ptr, &mesh, &t, &ind, &theta);
if (ind == -1) { /* Ray hit nothing */
active[off / 8] &= ~(1 << (off % 8));
- a_te_im_t[off] = a_te_re_t[off] = 0.f;
- a_tm_im_t[off] = a_tm_re_t[off] = 0.f;
- tau_t[off] = -1.f;
+ --num_active;
continue;
}
/* Calculate the reflection coefficients
@@ -546,14 +532,14 @@ void compute_paths(
r_tm_re /= free_space_loss;
r_tm_im /= free_space_loss;
/* Update the gains */
- a_te_re_new = a_te_re_t[off] * r_te_re - a_te_im_t[off] * r_te_im;
- a_te_im_new = a_te_re_t[off] * r_te_im + a_te_im_t[off] * r_te_re;
- a_tm_re_new = a_tm_re_t[off] * r_tm_re - a_tm_im_t[off] * r_tm_im;
- a_tm_im_new = a_tm_re_t[off] * r_tm_im + a_tm_im_t[off] * r_tm_re;
- a_te_re_t[off] = a_te_re_new;
- a_te_im_t[off] = a_te_im_new;
- a_tm_re_t[off] = a_tm_re_new;
- a_tm_im_t[off] = a_tm_im_new;
+ a_te_re_new = a_te_re_refl[off] * r_te_re - a_te_im_refl[off] * r_te_im;
+ a_te_im_new = a_te_re_refl[off] * r_te_im + a_te_im_refl[off] * r_te_re;
+ a_tm_re_new = a_tm_re_refl[off] * r_tm_re - a_tm_im_refl[off] * r_tm_im;
+ a_tm_im_new = a_tm_re_refl[off] * r_tm_im + a_tm_im_refl[off] * r_tm_re;
+ a_te_re_refl[off] = a_te_re_new;
+ a_te_im_refl[off] = a_te_im_new;
+ a_tm_re_refl[off] = a_tm_re_new;
+ a_tm_im_refl[off] = a_tm_im_new;
/* Update the delay */
tau_t[off] += t / SPEED_OF_LIGHT;
/* Advance the ray to the hit point */
@@ -562,47 +548,40 @@ void compute_paths(
/* Scatter a path from the hit point towards the rx */
r.o = r_ptr->o;
for (j = 0; j != num_rx; ++j) {
+ off_scat = i * num_rx * num_tx * num_paths + j * num_tx * num_paths + off;
r.d = (Vec3){rx_positions[j * 3], rx_positions[j * 3 + 1], rx_positions[j * 3 + 2]};
r.d = vec3_sub(&r.d, &r.o);
ind = -1;
moeller_trumbore(&r, &mesh, &t, &ind, &theta);
- if (ind == -1) continue;
+ if (ind != -1 && t <= 1.f) {
+ /* An obstacle between the hit point and the rx has been hit */
+ a_te_re_scat[off_scat] = a_te_im_scat[off_scat] = a_tm_re_scat[off_scat]
+ = a_tm_im_scat[off_scat] = tau_scat[off_scat] = 0.f;
+ continue;
+ }
/* No obstacle has been hit */
/* Calculate the distance */
t = sqrtf(vec3_dot(&r.d, &r.d));
/* Calculate the delay */
- tau[*num_paths_out] = tau_t[off] + t / SPEED_OF_LIGHT;
+ tau_scat[off_scat] = tau_t[off] + t / SPEED_OF_LIGHT;
/* Calculate the free space loss */
free_space_loss = free_space_loss_multiplier * t;
free_space_loss = free_space_loss * free_space_loss;
/* Set the gains */
- a_te_re[*num_paths_out] = a_te_re_t[off] / free_space_loss;
- a_te_im[*num_paths_out] = a_te_im_t[off] / free_space_loss;
- a_tm_re[*num_paths_out] = a_tm_re_t[off] / free_space_loss;
- a_tm_im[*num_paths_out] = a_tm_im_t[off] / free_space_loss;
- /* Increment the number of paths */
- ++(*num_paths_out);
- /* Realloc the output arrays if needed */
- if (*num_paths_out % num_paths_out_realloc_cap == 0) {
- ++num_cur_allocs;
- a_te_re = (float*)realloc(a_te_re, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_te_im = (float*)realloc(a_te_im, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_tm_re = (float*)realloc(a_tm_re, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- a_tm_im = (float*)realloc(a_tm_im, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- tau = (float*)realloc(tau, num_rx * num_tx * num_cur_allocs * num_paths_out_realloc_cap * sizeof(float));
- }
+ a_te_re_scat[off_scat] = a_te_re_refl[off] / free_space_loss;
+ a_te_im_scat[off_scat] = a_te_im_refl[off] / free_space_loss;
+ a_tm_re_scat[off_scat] = a_tm_re_refl[off] / free_space_loss;
+ a_tm_im_scat[off_scat] = a_tm_im_refl[off] / free_space_loss;
}
}
}
- /* Realloc the output arrays to the actual number of paths */
- a_te_re = (float*)realloc(a_te_re, num_rx * num_tx * *num_paths_out * sizeof(float));
- a_te_im = (float*)realloc(a_te_im, num_rx * num_tx * *num_paths_out * sizeof(float));
- a_tm_re = (float*)realloc(a_tm_re, num_rx * num_tx * *num_paths_out * sizeof(float));
- a_tm_im = (float*)realloc(a_tm_im, num_rx * num_tx * *num_paths_out * sizeof(float));
- tau = (float*)realloc(tau, num_rx * num_tx * *num_paths_out * sizeof(float));
-
/* Free */
+ free(a_te_re_refl);
+ free(a_te_im_refl);
+ free(a_tm_re_refl);
+ free(a_tm_im_refl);
+ free(tau_t);
free(h);
free(active);
free(rays);
diff --git a/compute_paths.h b/compute_paths.h
@@ -16,38 +16,50 @@
* \param rx_velocities receiver velocities, shape (num_rx, 3)
* \param tx_velocities transmitter velocities, shape (num_tx, 3)
* \param carrier_frequency carrier frequency in GHz. Must be > 0.0
- * \param sampling_frequency sampling frequency in Hz. Must be > 0.0
* \param num_rx number of receivers. Must be > 0
* \param num_tx number of transmitters. Must be > 0
* \param num_paths number of paths to compute. Must be > 0
* \param num_bounces number of bounces to compute. Must be > 0
- * \param num_samples number of samples in the original signal. Must be > 0
- * \param num_paths_out output number of paths computed. That includes LoS, scatter, and bounces that never left the scene.
- * \param a_te_re output array of real parts of transverse electric gains, shape (num_rx, num_tx, num_paths)
- * \param a_te_im output array of imaginary parts of transverse electric gains, shape (num_rx, num_tx, num_paths)
- * \param a_tm_re output array of real parts of transverse magnetic gains, shape (num_rx, num_tx, num_paths)
- * \param a_tm_im output array of imaginary parts of transverse magnetic gains, shape (num_rx, num_tx, num_paths)
- * \param tau output array of delays in seconds, shape (num_rx, num_tx, num_paths)
+ *
+ * Outputs:
+ *
+ * LoS:
+ * \param a_te_re_los output array of real parts of transverse electric gains for LoS, shape (num_rx, num_tx)
+ * \param a_te_im_los output array of imaginary parts of transverse electric gains for LoS, shape (num_rx, num_tx)
+ * \param a_tm_re_los output array of real parts of transverse magnetic gains for LoS, shape (num_rx, num_tx)
+ * \param a_tm_im_los output array of imaginary parts of transverse magnetic gains for LoS, shape (num_rx, num_tx)
+ * \param tau_los output array of delays for LoS in seconds, shape (num_rx, num_tx)
+ *
+ * Scatter:
+ * \param a_te_re_scat output array of real parts of transverse electric gains for scatter, shape (num_bounces, num_rx, num_tx, num_paths)
+ * \param a_te_im_scat output array of imaginary parts of transverse electric gains for scatter, shape (num_bounces, num_rx, num_tx, num_paths)
+ * \param a_tm_re_scat output array of real parts of transverse magnetic gains for scatter, shape (num_bounces, num_rx, num_tx, num_paths)
+ * \param a_tm_im_scat output array of imaginary parts of transverse magnetic gains for scatter, shape (num_bounces, num_rx, num_tx, num_paths)
+ * \param tau_scat output array of delays for scatter in seconds, shape (num_bounces, num_rx, num_tx, num_paths)
*/
void compute_paths(
- IN const char *mesh_filepath, /* path to the mesh file */
- IN const float *rx_positions, /* shape (num_rx, 3) */
- IN const float *tx_positions, /* shape (num_tx, 3) */
- IN const float *rx_velocities, /* shape (num_rx, 3) */
- IN const float *tx_velocities, /* shape (num_tx, 3) */
- IN float carrier_frequency, /* > 0.0 (IN GHz!) */
- IN float sampling_frequency, /* > 0.0 (IN Hz!) */
- IN size_t num_rx, /* number of receivers */
- IN size_t num_tx, /* number of transmitters */
- IN size_t num_paths, /* number of paths */
- IN size_t num_bounces, /* number of bounces */
- IN size_t num_samples, /* number of samples */
- OUT size_t *num_paths_out, /* number of paths computed */
- OUT float *a_te_re, /* output array real parts of TE gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_te_im, /* output array imaginary parts of TE gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_tm_re, /* output array real parts of TM gains (num_rx, num_tx, num_paths_out) */
- OUT float *a_tm_im, /* output array imaginary parts of TM gains (num_rx, num_tx, num_paths_out) */
- OUT float *tau /* output array of delays (num_rx, num_tx, num_paths_out) */
+ IN const char *mesh_filepath, /* path to the mesh file */
+ IN const float *rx_positions, /* shape (num_rx, 3) */
+ IN const float *tx_positions, /* shape (num_tx, 3) */
+ IN const float *rx_velocities, /* shape (num_rx, 3) */
+ IN const float *tx_velocities, /* shape (num_tx, 3) */
+ IN float carrier_frequency, /* > 0.0 (IN GHz!) */
+ IN size_t num_rx, /* number of receivers */
+ IN size_t num_tx, /* number of transmitters */
+ IN size_t num_paths, /* number of paths */
+ IN size_t num_bounces, /* number of bounces */
+ /* LoS */
+ OUT float *a_te_re_los, /* output array real parts of TE gains (num_rx, num_tx) */
+ OUT float *a_te_im_los, /* output array imaginary parts of TE gains (num_rx, num_tx) */
+ OUT float *a_tm_re_los, /* output array real parts of TM gains (num_rx, num_tx) */
+ OUT float *a_tm_im_los, /* output array imaginary parts of TM gains (num_rx, num_tx) */
+ OUT float *tau_los, /* output array of delays (num_rx, num_tx) */
+ /* Scatter */
+ OUT float *a_te_re_scat, /* output array real parts of TE gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_te_im_scat, /* output array imaginary parts of TE gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_tm_re_scat, /* output array real parts of TM gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *a_tm_im_scat, /* output array imaginary parts of TM gains (num_bounces, num_rx, num_tx, num_paths) */
+ OUT float *tau_scat /* output array of delays (num_bounces, num_rx, num_tx, num_paths) */
);
-#endif /* COMPUTE_PATHS_H */
+#endif /* COMPUTE_PATHS_H */
diff --git a/compute_paths_pybind11.cpp b/compute_paths_pybind11.cpp
@@ -12,9 +12,13 @@ extern "C" {
namespace py = pybind11;
std::tuple<
- py::array_t<float>, py::array_t<float>, // a_te_re, a_te_im
- py::array_t<float>, py::array_t<float>, // a_tm_re, a_tm_im
- py::array_t<float> > // tau
+ py::array_t<float>, py::array_t<float>, // a_te_re_los, a_te_im_los
+ py::array_t<float>, py::array_t<float>, // a_tm_re_los, a_tm_im_los
+ py::array_t<float>, // tau_los
+ py::array_t<float>, py::array_t<float>, // a_te_re_scat, a_te_im_scat
+ py::array_t<float>, py::array_t<float>, // a_tm_re_scat, a_tm_im_scat
+ py::array_t<float> // tau_scat
+>
compute_paths_wrapper(
const std::string &mesh_filepath,
py::array_t<float> rx_positions,
@@ -22,12 +26,10 @@ compute_paths_wrapper(
py::array_t<float> rx_velocities,
py::array_t<float> tx_velocities,
float carrier_frequency,
- float sampling_frequency,
int num_rx,
int num_tx,
int num_paths,
- int num_bounces,
- int num_samples
+ int num_bounces
) {
// Prepare input arrays (this is a basic implementation, check shapes and memory layout)
py::buffer_info rx_pos_info = rx_positions.request();
@@ -36,8 +38,16 @@ compute_paths_wrapper(
py::buffer_info tx_vel_info = tx_velocities.request();
// Output
- size_t num_paths_out;
- float *a_te_re, *a_te_im, *a_tm_re, *a_tm_im, *tau;
+ float *a_te_re_los = new float[num_rx * num_tx];
+ float *a_te_im_los = new float[num_rx * num_tx];
+ float *a_tm_re_los = new float[num_rx * num_tx];
+ float *a_tm_im_los = new float[num_rx * num_tx];
+ float *tau_los = new float[num_rx * num_tx];
+ float *a_te_re_scat = new float[num_bounces * num_rx * num_tx * num_paths];
+ float *a_te_im_scat = new float[num_bounces * num_rx * num_tx * num_paths];
+ float *a_tm_re_scat = new float[num_bounces * num_rx * num_tx * num_paths];
+ float *a_tm_im_scat = new float[num_bounces * num_rx * num_tx * num_paths];
+ float *tau_scat = new float[num_bounces * num_rx * num_tx * num_paths];
// Call the C function
compute_paths(
@@ -47,38 +57,49 @@ compute_paths_wrapper(
(const float*)rx_vel_info.ptr, // Rx velocities
(const float*)tx_vel_info.ptr, // Tx velocities
carrier_frequency, // Carrier frequency in GHz
- sampling_frequency, // Sampling frequency in Hz
(size_t)num_rx,
(size_t)num_tx,
(size_t)num_paths,
(size_t)num_bounces,
- (size_t)num_samples,
- &num_paths_out,
- a_te_re, a_te_im, a_tm_re, a_tm_im, // Gains
- tau // Delays
+ a_te_re_los, a_te_im_los, a_tm_re_los, a_tm_im_los, tau_los, // LoS outputs
+ a_te_re_scat, a_te_im_scat, a_tm_re_scat, a_tm_im_scat, tau_scat // Scatter outputs
);
// Convert output arrays into numpy arrays for easy use in Python
// TODO prevent copying data
// TODO construct np.complex64 array for a
- py::array_t<float> a_te_re_array = py::array_t<float>({num_rx, num_tx, (int)num_paths_out}, a_te_re);
- py::array_t<float> a_te_im_array = py::array_t<float>({num_rx, num_tx, (int)num_paths_out}, a_te_im);
- py::array_t<float> a_tm_re_array = py::array_t<float>({num_rx, num_tx, (int)num_paths_out}, a_tm_re);
- py::array_t<float> a_tm_im_array = py::array_t<float>({num_rx, num_tx, (int)num_paths_out}, a_tm_im);
- py::array_t<float> tau_array = py::array_t<float>({num_rx, num_tx, (int)num_paths_out}, tau);
+ py::array_t<float> a_te_re_los_array = py::array_t<float>({num_rx, num_tx}, a_te_re_los);
+ py::array_t<float> a_te_im_los_array = py::array_t<float>({num_rx, num_tx}, a_te_im_los);
+ py::array_t<float> a_tm_re_los_array = py::array_t<float>({num_rx, num_tx}, a_tm_re_los);
+ py::array_t<float> a_tm_im_los_array = py::array_t<float>({num_rx, num_tx}, a_tm_im_los);
+ py::array_t<float> tau_los_array = py::array_t<float>({num_rx, num_tx}, tau_los);
+ py::array_t<float> a_te_re_scat_array = py::array_t<float>({num_bounces, num_rx, num_tx, num_paths}, a_te_re_scat);
+ py::array_t<float> a_te_im_scat_array = py::array_t<float>({num_bounces, num_rx, num_tx, num_paths}, a_te_im_scat);
+ py::array_t<float> a_tm_re_scat_array = py::array_t<float>({num_bounces, num_rx, num_tx, num_paths}, a_tm_re_scat);
+ py::array_t<float> a_tm_im_scat_array = py::array_t<float>({num_bounces, num_rx, num_tx, num_paths}, a_tm_im_scat);
+ py::array_t<float> tau_scat_array = py::array_t<float>({num_bounces, num_rx, num_tx, num_paths}, tau_scat);
// Deallocate arrays
- delete[] a_te_re;
- delete[] a_te_im;
- delete[] a_tm_re;
- delete[] a_tm_im;
- delete[] tau;
+ delete[] a_te_re_los;
+ delete[] a_te_im_los;
+ delete[] a_tm_re_los;
+ delete[] a_tm_im_los;
+ delete[] tau_los;
+ delete[] a_te_re_scat;
+ delete[] a_te_im_scat;
+ delete[] a_tm_re_scat;
+ delete[] a_tm_im_scat;
+ delete[] tau_scat;
// Return the results as a tuple (gains, delays)
return std::make_tuple(
- a_te_re_array, a_te_im_array,
- a_tm_re_array, a_tm_im_array,
- tau_array);
+ a_te_re_los_array, a_te_im_los_array,
+ a_tm_re_los_array, a_tm_im_los_array,
+ tau_los_array,
+ a_te_re_scat_array, a_te_im_scat_array,
+ a_tm_re_scat_array, a_tm_im_scat_array,
+ tau_scat_array
+ );
}
PYBIND11_MODULE(rt, m) {
@@ -89,12 +110,10 @@ PYBIND11_MODULE(rt, m) {
py::arg("rx_velocities"),
py::arg("tx_velocities"),
py::arg("carrier_frequency"),
- py::arg("sampling_frequency"),
py::arg("num_rx"),
py::arg("num_tx"),
py::arg("num_paths"),
- py::arg("num_bounces"),
- py::arg("num_samples"));
+ py::arg("num_bounces"));
// Scene filepaths
m.def("get_scene_fp_box",
diff --git a/test.c b/test.c
@@ -15,14 +15,20 @@ int main(int argc, char **argv)
float rx_velocities[3] = {0.0, 0.0, 0.0};
float tx_velocities[3] = {0.0, 0.0, 0.0};
float carrier_frequency = 3.0; /* 3 GHz */
- float sampling_frequency = 3*1e9; /* 3 GHz */
size_t num_rx = 1;
size_t num_tx = 1;
size_t num_paths = 10000;
size_t num_bounces = 3;
- size_t num_samples = 150;
- size_t num_paths_out;
- float *a_te_re, *a_te_im, *a_tm_re, *a_tm_im, *tau;
+ float *a_te_re_los = (float*)malloc(num_rx * num_tx * sizeof(float));
+ float *a_te_im_los = (float*)malloc(num_rx * num_tx * sizeof(float));
+ float *a_tm_re_los = (float*)malloc(num_rx * num_tx * sizeof(float));
+ float *a_tm_im_los = (float*)malloc(num_rx * num_tx * sizeof(float));
+ float *tau_los = (float*)malloc(num_rx * num_tx * sizeof(float));
+ float *a_te_re_scat = (float*)malloc(num_bounces * num_rx * num_tx * num_paths * sizeof(float));
+ float *a_te_im_scat = (float*)malloc(num_bounces * num_rx * num_tx * num_paths * sizeof(float));
+ float *a_tm_re_scat = (float*)malloc(num_bounces * num_rx * num_tx * num_paths * sizeof(float));
+ float *a_tm_im_scat = (float*)malloc(num_bounces * num_rx * num_tx * num_paths * sizeof(float));
+ float *tau_scat = (float*)malloc(num_bounces * num_rx * num_tx * num_paths * sizeof(float));
compute_paths(
argv[1],
rx_positions,
@@ -30,26 +36,27 @@ int main(int argc, char **argv)
rx_velocities,
tx_velocities,
carrier_frequency,
- sampling_frequency,
- num_rx, num_tx, num_paths, num_bounces, num_samples,
- &num_paths_out,
- a_te_re, a_te_im, a_tm_re, a_tm_im,
- tau);
+ num_rx, num_tx, num_paths, num_bounces,
+ a_te_re_los, a_te_im_los, a_tm_re_los, a_tm_im_los, tau_los,
+ a_te_re_scat, a_te_im_scat, a_tm_re_scat, a_tm_im_scat, tau_scat
+ );
/* Save the results */
- FILE *f = fopen("a_te_re.bin", "wb");
- fwrite(a_te_re, sizeof(float), num_rx * num_tx * num_paths_out, f);
- fclose(f);
- f = fopen("a_te_im.bin", "wb");
- fwrite(a_te_im, sizeof(float), num_rx * num_tx * num_paths_out, f);
- fclose(f);
- f = fopen("a_tm_re.bin", "wb");
- fwrite(a_tm_re, sizeof(float), num_rx * num_tx * num_paths_out, f);
- fclose(f);
- f = fopen("a_tm_im.bin", "wb");
- fwrite(a_tm_im, sizeof(float), num_rx * num_tx * num_paths_out, f);
- fclose(f);
- f = fopen("tau.bin", "wb");
- fwrite(tau, sizeof(float), num_rx * num_tx * num_paths_out, f);
- fclose(f);
+
+ FILE *f;
+ #define WRITE_BIN(name, data, size) \
+ f = fopen(name, "wb"); \
+ fwrite(data, sizeof(float), size, f); \
+ fclose(f);
+
+ WRITE_BIN("a_te_re_los.bin", a_te_re_los, num_rx * num_tx);
+ WRITE_BIN("a_te_im_los.bin", a_te_im_los, num_rx * num_tx);
+ WRITE_BIN("a_tm_re_los.bin", a_tm_re_los, num_rx * num_tx);
+ WRITE_BIN("a_tm_im_los.bin", a_tm_im_los, num_rx * num_tx);
+ WRITE_BIN("tau_los.bin", tau_los, num_rx * num_tx);
+ WRITE_BIN("a_te_re_scat.bin", a_te_re_scat, num_bounces * num_rx * num_tx * num_paths);
+ WRITE_BIN("a_te_im_scat.bin", a_te_im_scat, num_bounces * num_rx * num_tx * num_paths);
+ WRITE_BIN("a_tm_re_scat.bin", a_tm_re_scat, num_bounces * num_rx * num_tx * num_paths);
+ WRITE_BIN("a_tm_im_scat.bin", a_tm_im_scat, num_bounces * num_rx * num_tx * num_paths);
+ WRITE_BIN("tau_scat.bin", tau_scat, num_bounces * num_rx * num_tx * num_paths);
}
diff --git a/test.py b/test.py
@@ -1,4 +1,5 @@
import numpy as np
+import matplotlib.pyplot as plt
from rt import compute_paths
# Define inputs
@@ -8,33 +9,48 @@ tx_positions = np.array([[0., 0., 2.5]], dtype=np.float64)
rx_velocities = np.array([[0., 0., 0.]], dtype=np.float64)
tx_velocities = np.array([[0., 0., 0.]], dtype=np.float64)
carrier_frequency = 3.0
-sampling_rate = 3*1e9
num_rx = 1
num_tx = 1
num_paths = 10000
num_bounces = 3
-num_samples = 150
# Call compute_paths
-a_te_re, a_te_im, a_tm_re, a_tm_im, tau = compute_paths(
+a_te_re_los, a_te_im_los, a_tm_re_los, a_tm_im_los, tau_los, a_te_re_scat, a_te_im_scat, a_tm_re_scat, a_tm_im_scat, tau_scat = compute_paths(
mesh_filepath,
rx_positions,
tx_positions,
rx_velocities,
tx_velocities,
carrier_frequency,
- sampling_rate,
num_rx,
num_tx,
num_paths,
- num_bounces,
- num_samples
+ num_bounces
)
-print(f"TE Gains: {a_te_re + 1.j*a_te_im}")
-print(f"TM Gains: {a_tm_re + 1.j*a_tm_im}")
-print(f"Delays: {tau}")
-print(f"a_te_re min, max: {np.min(a_te_re)}, {np.max(a_te_re)}")
-print(f"a_te_im min, max: {np.min(a_te_im)}, {np.max(a_te_im)}")
-print(f"a_tm_re min, max: {np.min(a_tm_re)}, {np.max(a_tm_re)}")
-print(f"a_tm_im min, max: {np.min(a_tm_im)}, {np.max(a_tm_im)}")
+print("#####################")
+print("LoS:")
+print("#####################")
+print(f"shape(a_te_re_los): {a_te_re_los.shape}")
+print(f"TE Gains: {a_te_re_los + 1.j*a_te_im_los}")
+print(f"TM Gains: {a_tm_re_los + 1.j*a_tm_im_los}")
+print(f"Delays: {tau_los}")
+print(f"a_te_re_los min, max: {np.min(a_te_re_los)}, {np.max(a_te_re_los)}")
+print(f"a_te_im_los min, max: {np.min(a_te_im_los)}, {np.max(a_te_im_los)}")
+print(f"a_tm_re_los min, max: {np.min(a_tm_re_los)}, {np.max(a_tm_re_los)}")
+print(f"a_tm_im_los min, max: {np.min(a_tm_im_los)}, {np.max(a_tm_im_los)}")
+
+print("\n#####################")
+print("Scatter:")
+print("#####################")
+print(f"shape(a_te_re_scat): {a_te_re_scat.shape}")
+print(f"TE Gains: {a_te_re_scat + 1.j*a_te_im_scat}")
+print(f"TM Gains: {a_tm_re_scat + 1.j*a_tm_im_scat}")
+print(f"Delays: {tau_scat}")
+print(f"a_te_re_scat min, max: {np.min(a_te_re_scat)}, {np.max(a_te_re_scat)}")
+print(f"a_te_im_scat min, max: {np.min(a_te_im_scat)}, {np.max(a_te_im_scat)}")
+print(f"a_tm_re_scat min, max: {np.min(a_tm_re_scat)}, {np.max(a_tm_re_scat)}")
+print(f"a_tm_im_scat min, max: {np.min(a_tm_im_scat)}, {np.max(a_tm_im_scat)}")
+
+assert a_te_re_los.shape == a_te_im_los.shape == a_tm_re_los.shape == a_tm_im_los.shape == tau_los.shape
+assert a_te_re_scat.shape == a_te_im_scat.shape == a_tm_re_scat.shape == a_tm_im_scat.shape == tau_scat.shape
