for (double theta=0.0, thetaIndex=0; theta<PI; thetaIndex++, theta+=PI/numAngleCells) { rho = rhoMax + x*cos(theta) + y*sin(theta); accum[thetaIndex][(int)rho]++; }