Производительность OpenCL

#pragma OPENCL EXTENSION cl_khr_fp64 : enable
__constant double rzem = 6371.1;
__constant double eps   = 1e-9;
__constant sampler_t sampler = CLK_NORMALIZED_COORDS_FALSE | CLK_ADDRESS_CLAMP_TO_EDGE | CLK_FILTER_NEAREST;
__constant uchar wt_water_bit = 128;
__constant uchar wt_cold_bit = 2;
  
__kernel void keGroundWaterCount    (   __global double2 * pgws,
                                                                        __global double2 * da,
                                                                        __global double * steps,
                                                                        __global int * isteps,
                                                                        const int iter,
                                                                        const int iterOffs,
                                                                        const double2 station,
                                                                        const int4 map_param,
                                                                        __read_only image3d_t map
                                                                    )
{
    int ipoint = get_global_id(0);
     
    int water_hot = 0;
    int water_cold = 0;
    double ar = radians( da[ipoint].y );
    uchar v0 = 0;
     
    double2 coord,d,norm;
    double dist,c3,sinar,sinc3,PiSubSt,dd,ss;
    int2 tile,m,p;
    int z;
     
    sinar = sin( ar );
    PiSubSt = radians( 90.0 - station.y );
    double cosPiSubStcosar = cos(PiSubSt)*cos( ar );
    double sinPiSubSt = sin(PiSubSt);
     
    double sdist = steps[ipoint];
    int i = iter*iterOffs;
    int limit = min( isteps[ipoint], i + iterOffs );
     
    while( i < limit )
    {
        dist = i * sdist;
        c3 = dist / rzem;
        if( c3 < eps )
        {
            coord = station;
        }else
        {  
            sinc3 = sin( c3 );
            dd = atan( sinar /( sinPiSubSt/(sinc3/cos(c3)) - cosPiSubStcosar ));
            ss = sinc3*sinar / sin(dd);
            if( ss >= 1.0 )
            {
                ss = 90.0;
            }else
            {
                if( ss <= -1.0 )
                {
                    ss = 270.0;
                }else
                {
                    ss = degrees( atan( ss / sqrt( 1.0 - ss*ss )));
                }
            }
            coord = (double2)( station.x + degrees(dd), 90.0 - ss );
            if( coord.x < -180.0 )
                coord.x += 360.0;
            if( coord.x >= 180.0 )
                coord.x -= 360.0;
        }
         
        d = (double2)( coord.x + 180.0, 90.0 - coord.y );
     
        tile = convert_int2( d );
        norm = d - convert_double2( tile );
  
        m = tile - map_param.xy;   
        z = m.y * map_param.z + m.x;
      
        p = convert_int2( norm * (double2)( get_image_width( map ), get_image_height( map ) ));
         
        v0 = read_imagei( map, sampler, (int4)( p.x, p.y, z, 0 ) ).x;  
         
        if( v0 & wt_water_bit )
        {
            if( v0 & wt_cold_bit )
            {
                water_cold++;
            }else
            {
                water_hot++;
            }
        }
        i++;
    }
    pgws[ipoint].x += water_cold;
    pgws[ipoint].y += water_hot;
}