/tomo/pyhst

#define HST_NEWTEX4                             //!< Assymetric: process 4 projections and 64 bins during each step.

 

//#define HST_SQUARE_PPT                                //!< NO. Little effect on performance (slighly reduces register usage, but doesn't affect computations as the loops are unrolled anyway)

//#define HST_FLOAT_LOOPS                               //!< NO. Use float loops to optimize register usage (mostly bad)

    const int bidx = get_group_id(0) * get_local_size(0);

    const int bidy = get_group_id(1) * get_local_size(1);

 

    const int xy_minor = tidx % 4;

    const int y_minor = xy_minor / 2;

    const int x_minor = xy_minor % 2;

    const int xy_major = tidy % 4;

    const int y_major = xy_major / 2;

    const int x_major = xy_major % 2;

 

#if defined(HST_NEWTEX4)

# ifdef SQUARE_PPT

    const int blocky = block/2;

    const int blockx = block%2;

 

    const int sidy = blocky * 4 + y_major * 2 + y_minor;

    const int sidx = blockx * 4 + x_major * 2 + x_minor;

#else

    const int ridy = tidy;

    const int ridx = tidx;

 

    const int sidy = y_major * 2 + y_minor;

    const int sidx = x_major * 2 + x_minor;

#endif

 

    const float x = bidx + sidx + apos_off_x;

    const float axis = c_all[0].z;

#endif

 

#ifdef HST_CACHE_CONST

    const int tid = 16 * tidy + tidx;

    __local float2 w_trig[CACHE_BLOCK];

    __local float  w_ofst[CACHE_BLOCK];

# endif

 

    for (int pblock = 0; pblock < num_proj; pblock += CACHE_BLOCK) {

        const int pblock_end = min(CACHE_BLOCK, num_proj - pblock);

        if (tid < CACHE_BLOCK) {

            const int g_proj = pblock + tid;

            w_trig[tid] = (float2){ c_all[g_proj].x, c_all[g_proj].y };

            w_ofst[tid] = c_all[g_proj].z;

#endif

        }

        barrier(CLK_LOCAL_MEM_FENCE);

 

        for (int proje = proj; proje < pblock_end; proje += PROJ_STEP) {

            const float projf = pblock + proje + 0.5f;

# ifdef HST_NO_OFFSETS

            const float4 all = (float4){w_trig[proje].x, w_trig[proje].y, axis, 0};

# else

            const float4 all = (float4){w_trig[proje].x, w_trig[proje].y, w_ofst[proje], 0};

# endif

#else

# if defined(HST_FLOAT_LOOPS)

        for (float projf =  proj + 0.5f; projf < num_proj; projf += PROJ_STEP) {