1

I have started optimising my code using SSE. Essentially it is a ray tracer that processes 4 rays at a time by storing the coordinates in __m128 data types x, y, z (the coordinates for the four rays are grouped by axis). However I have a branched statement which protects against divide by zero I can't seem to convert to SSE. In serial this is:

const float d = wZ == -1.0f ? 1.0f/( 1.0f-wZ) : 1.0f/(1.0f+wZ);

Where wZ is the z-coordinate and this calculation needs to be done for all four rays.

How could I translate this into SSE?

I have been experimenting using the SSE equals comparison as follows (now wz pertains to a __m128 data type containing the z values for each of the four rays):

_mm_cmpeq_ps(_mm_set1_ps(-1.0f) , wZ )

And then using this to identify cases where wZ[x] = -1.0, taking the absolute value of this case and then continue the calculation as normal.

However I have not had much success in this endeavour.

cubiclewar
  • 1,569
  • 3
  • 20
  • 37

1 Answers1

5

Here's a fairly straightforward solution which just implements the scalar code with SSE without any further optimisation. It can probably be made a little more efficient, e.g. by exploiting the fact that the result will be 0.5 when wZ = -1.0, or perhaps even by just doing the division regardless and then converting the INFs to 0.5 after the fact.

I've #ifdefd for SSE4 versus pre-SSE4, since SSE4 has a "blend" instruction which may be a little more efficient that the three pre-SSE4 instructions that are otherwise needed to mask and select values.

#include <emmintrin.h>
#ifdef __SSE4_1__
#include <smmintrin.h>
#endif

#include <stdio.h>

int main(void)
{
    const __m128 vk1 = _mm_set1_ps(1.0f);       // useful constants
    const __m128 vk0 = _mm_set1_ps(0.0f);

    __m128 wZ, d, d0, d1, vcmp;
#ifndef __SSE4_1__  // pre-SSE4 implementation
    __m128 d0_masked, d1_masked;
#endif

    wZ = _mm_set_ps(-1.0f, 0.0f, 1.0f, 2.0f);   // test inputs

    d0 = _mm_add_ps(vk1, wZ);                   // d0 = 1.0 - wZ
    d1 = _mm_sub_ps(vk1, wZ);                   // d1 = 1.0 + wZ
    vcmp = _mm_cmpneq_ps(d1, vk0);              // test for d1 != 0.0, i.e. wZ != -1.0
#ifdef __SSE4_1__   // SSE4 implementation
    d = _mm_blendv_ps(d0, d1, vcmp);
#else               // pre-SSE4 implementation
    d0_masked = _mm_andnot_ps(vcmp, d0);
    d1_masked = _mm_and_ps(vcmp, d1);
    d = _mm_or_ps(d0_masked, d1_masked);       // d = wZ == -1.0 ? 1.0 / (1.0 - wZ) : 1.0 / (1.0 + wZ)
#endif
   d = _mm_div_ps(vk1, d);

   printf("wZ = %vf\n", wZ);
   printf("d = %vf\n", d);

   return 0;
}
Paul R
  • 208,748
  • 37
  • 389
  • 560
  • exactly what I was after. There is a few operations I will need to read up on to fully understand the code but I generates the correct results. Out of curiosity can inf or nan (what ever 1/0 evaluates to) be easily identified in SSE and replaced? – cubiclewar Nov 04 '11 at 09:34
  • I haven't tried it but I *think* that you can exploit the fact that `_mm_cmpeq_ps(v, v)` will return false when `v` is `INF` or `NaN` - I may try for another solution using this method later if I get time... – Paul R Nov 04 '11 at 09:41
  • 1
    I tried @PaulR's suggestion to filter out `INF`/`NaN` using `_mm_cmpeq_ps(v, v)` as a bitmask and it seems to work fine. – Rotem Apr 30 '15 at 20:55
  • 1
    It works to identify INF or NaN. To distinguish between the two, you can compare the absolute value against _mm_set1_ps(__builtin_inff()) or on Windows, _mm_set1_ps(HUGE_VALF) – Jens Munk Jun 18 '15 at 13:55