Cubic Interpolation with the official formula fails

Question

I am trying to implement the Cubic Interpolation method using the next formula when a=-0.5 as usual.

My Linear Interpolation and Nearest Neighbor interpolation is working great but for some reason the Cubic interpolation fails with white pixels and turn them sometimes to turquoise color and sometimes messing around with another colors.

for example using rotation: (NOTE: please look carefully on the right image and you will notice the problems)

Another Example with much more black pixels. It almost seems to work perfectly but look on the dog's tongue. (strong white pixels turn to turquoise again)

you can see that my implementation of the Linear Interpolation is working great:

Since the actual rotation worked, I think I have a small mistake in the code that I did not notice, or maybe it's a numeric error or a double / float error. It is important to note that I read the image normally and store the destination image as follows:

cv::Mat img = cv::imread("../dogfails.jpeg");
cv::Mat rotatedImageCubic(img.rows,img.cols,CV_8UC3);

Clarifications:

1. Inside my cubic interpolation function, srcPoint (newX and newY) is the "landing point" from the inverse transformation.

2. In my inverse transformations I am not using matrix multiplication with the pixels, right now I am just using the formulas for rotation. It might be important for the "numerical errors". For example:

   rotatedX = x * cos(angle * toRadian) + y * sin(angle * toRadian);
   rotatedY = x * (-sin(angle * toRadian)) + y * cos(angle * toRadian);
   

Here is my code for the Cubic Interpolation

double cubicEquationSolver(double d,double a) {
    d = abs(d);
    if( 0.0 <= d  && d <= 1.0) {
        double score = (a + 2.0) * pow(d, 3.0) - ((a + 3.0) * pow(d, 2.0)) + 1.0;
        return score;
    }
    else if(1 < d && d <= 2) {
        double score = a * pow(d, 3.0) - 5.0*a * pow(d, 2.0) + 8.0*a * d - 4.0*a;
        return score;
    }
    else
        return 0.0;
}

void Cubic_Interpolation_Helper(const cv::Mat& src, cv::Mat& dst, const cv::Point2d& srcPoint, cv::Point2i& dstPixel) {
    double newX = srcPoint.x;
    double newY = srcPoint.y;
    double dx = abs(newX - round(newX));
    double dy = abs(newY - round(newY));

    double sumCubicBValue = 0;
    double sumCubicGValue = 0;
    double sumCubicRValue = 0;
    double sumCubicGrayValue = 0;
    double uX = 0;
    double uY = 0;

    if (floor(newX) - 1  < 0 || floor(newX) + 2  > src.cols - 1 || floor(newY) < 0 || floor(newY)  > src.rows - 1) {
        if (dst.channels() > 1)
            dst.at<cv::Vec3b>(dstPixel) = cv::Vec3b(0, 0,0);
        else
            dst.at<uchar>(dstPixel) = 0;
    }
    else {
        for (int cNeighbor = -1; cNeighbor <= 2; cNeighbor++) {
            for (int rNeighbor = -1; rNeighbor <= 2; rNeighbor++) {
                uX = cubicEquationSolver(rNeighbor + dx, -0.5);
                uY = cubicEquationSolver(cNeighbor + dy, -0.5);
                if (src.channels() > 1) {
                    sumCubicBValue = sumCubicBValue + (double) src.at<cv::Vec3b>(
                            cv::Point2i(round(newX) + rNeighbor, cNeighbor + round(newY)))[0] * uX * uY;
                    sumCubicGValue = sumCubicGValue + (double) src.at<cv::Vec3b>(
                            cv::Point2i(round(newX) + rNeighbor, cNeighbor + round(newY)))[1] * uX * uY;
                    sumCubicRValue = sumCubicRValue + (double) src.at<cv::Vec3b>(
                            cv::Point2i(round(newX) + rNeighbor, cNeighbor + round(newY)))[2] * uX * uY;
                } else {
                    sumCubicGrayValue = sumCubicGrayValue + (double) src.at<uchar>(
                            cv::Point2i(round(newX) + rNeighbor, cNeighbor + round(newY))) * uX * uY;
                }
            }
        }
        if (dst.channels() > 1)
            dst.at<cv::Vec3b>(dstPixel) = cv::Vec3b((int) round(sumCubicBValue), (int) round(sumCubicGValue),
                                                    (int) round(sumCubicRValue));
        else
            dst.at<uchar>(dstPixel) = sumCubicGrayValue;
    }

I hope someone here will be able to help me, Thanks!