Demo entry 6625531

back projection

   

Submitted by anonymous on Jun 21, 2017 at 10:42
Language: C. Code size: 3.4 kB.

#define BP_NPIX_X 855
#define BP_NPIX_Y 956
#define N_RANGE_UPSAMPLED 1067
#define N_PULSES 2179

extern double sqrt(double _);
extern double sin(double _);
extern double cos(double _);

inline void norm(double *n, double x, double y, double z) {
    *n = sqrt(x*x + y*y + z*z);
}

// Input: arg
//   Output: *sine = sin(arg); *cosine = cos(arg) 
inline void sin_cos(float *sine, float *cosine, double arg) {
    *sine = (float) sin(arg);
    *cosine = (float) cos(arg);
}

void bin_sample(float * sample_r, float * sample_i, 
                float data_r[N_RANGE_UPSAMPLED],
                float data_i[N_RANGE_UPSAMPLED],
                const double * bin) {  
    const double b = *bin;
    if (b >= 0) {
	if (b < N_RANGE_UPSAMPLED-1) {
	    /* interpolation range is [bin_floor, bin_floor+1] */
	    const int bin_floor = (int) b;
	    /* interpolation weight */
	    const float w = (float) (b - (double) bin_floor);
	    /* linearly interpolate to obtain a sample at bin */
	    *sample_r = (1.0f-w)*data_r[bin_floor] + w*data_r[bin_floor+1];
	    *sample_i = (1.0f-w)*data_i[bin_floor] + w*data_i[bin_floor+1];
	    /* compute the complex exponential for the matched filter */
	    return;
	}
    }
    *sample_r = 0.0f;
    *sample_i = 0.0f;
}

void sar_backprojection2(
    float image_r[BP_NPIX_Y][BP_NPIX_X],
    float image_i[BP_NPIX_Y][BP_NPIX_X],
    float (* const data_r)[N_RANGE_UPSAMPLED],
    float (* const data_i)[N_RANGE_UPSAMPLED],
    const double platpos_x[N_PULSES],
    const double platpos_y[N_PULSES],
    const double platpos_z[N_PULSES],
    double * ku,
    double * R0,
    double * dR,
    double * dxdy,
    double * z0)
{
    const double dR_inv = 1.0/(*dR);
    const double dxdy_2 = (*dxdy) * 0x1.0p-1;
    const double ku_2 = (*ku) * 2.0;
    int p, ix, iy;

    for (iy = 0; iy < BP_NPIX_Y; ++iy)
    {
        const double py = ((double)(-BP_NPIX_Y + 1 + 2 * iy)) * dxdy_2;
        for (ix = 0; ix < BP_NPIX_X; ++ix)
        {
          const double px = ((double)(-BP_NPIX_X + 1 + 2 * ix)) * dxdy_2;

	  double contrib_r = 0.0;
	  double contrib_i = 0.0;

	  for (p = 0; p < N_PULSES; ++p)
            {
              float sample_r, sample_i, matched_filter_r, matched_filter_i, 
                prod_r, prod_i;
              double R;
              /* calculate the range R from the platform to this pixel */
              const double xdiff = platpos_x[p] - px;
              const double ydiff = platpos_y[p] - py;
              const double zdiff = platpos_z[p] - (*z0);
              double bin;
              norm(&R, xdiff, ydiff, zdiff);

              /* convert to a range bin index */
              bin = (R-(*R0))*dR_inv;
              
              bin_sample(&sample_r, &sample_i, data_r[p], data_i[p], &bin);
	      sin_cos(&matched_filter_i, &matched_filter_r, ku_2*R);

              /* scale the interpolated sample by the matched filter */
              prod_r = 
                sample_r*matched_filter_r - sample_i * matched_filter_i;
              prod_i = 
                sample_r*matched_filter_i + sample_i * matched_filter_r;
              /* accumulate this pulse's contribution into the pixel */
	      contrib_r += prod_r;
	      contrib_i += prod_i;
            }
	  image_r[iy][ix] = (float) contrib_r;
	  image_i[iy][ix] = (float) contrib_i;

        }
    }
}

This snippet took 0.01 seconds to highlight.

Back to the Entry List or Home.

Delete this entry (admin only).