Drawing tangent lines for each point in bezier curve

Question

I managed to draw a bezier curve like:

glColor3f(0,1,0);
glBegin(GL_LINE_STRIP);
for (int i = 3; i < nPt; i+=3) {
    glColor3f(0,0,0);
    for (float k = 0; k < NLINESEGMENT+1; k++) {
        float x = pow(1.0-k/NLINESEGMENT,3)*ptList[i-3].x +
            3*(k/NLINESEGMENT)*pow(1.0-k/NLINESEGMENT, 2) * ptList[i-2].x +
            3*(1.0-k/NLINESEGMENT)*pow(k/NLINESEGMENT, 2) * ptList[i-1].x +
            pow(k/NLINESEGMENT, 3)*ptList[i].x;
        float y = pow(1.0-k/NLINESEGMENT,3)*ptList[i-3].y +
            3*(k/NLINESEGMENT)*pow(1.0-k/NLINESEGMENT, 2) * ptList[i-2].y +
            3*(1.0-k/NLINESEGMENT)*pow(k/NLINESEGMENT, 2) * ptList[i-1].y +
            pow(k/NLINESEGMENT, 3)*ptList[i].y;
        glVertex2d(x,y);
    }
}
glEnd();

Now I want to add tangent arrows for each point, how can I do that? I am given a function that draws an arrow. So I believe I need to just rotate the reference frame and draw that arrow. But how do I compute the rotation? I think I need to differenciate the equations, but the question still remains, how do I use that?

UPDATE

As every 4th point is put, a curve is drawn.

I am supposed to achieve something like below

Full Source

UPDATE 2

Ok I made an attempt at drawing the tangents like:

glColor3f(0,1,0);
for (int i = 3; i < nPt; i+=3) {
    for (int n = 0; n < NOBJECTONCURVE; n++) {
        float t = (float)n/NOBJECTONCURVE;
        float x0 = points[i-3].x,
                x1 = points[i-2].x,
                x2 = points[i-1].x, 
                x3 = points[i].x;
        float y0 = points[i-3].y,
                y1 = points[i-2].y,
                y2 = points[i-1].y, 
                y3 = points[i].y;

        float x = pow(1.0-t, 3) * points[i-3].x +
            3 * t * pow(1.0 - t, 2) * points[i-2].x +
            3 * (1.0 - t) * pow(t, 2) * points[i-1].x +
            pow(t, 3)*points[i].x;
        float y = pow(1.0-t, 3) * points[i-3].y +
            3 * t * pow(1.0 - t, 2) * points[i-2].y +
            3 * (1.0 - t) * pow(t, 2) * points[i-1].y +
            pow(t, 3)*points[i].y;

        float dx = -3*(1-t)*x0 + 3*x1*((2*t)*(t-1)+pow((1-t),2)) + 3*x2*(2*t*(1-t)-pow(t,2)) + 3*pow(t,2)*x3;
        float dy = -3*(1-t)*y0 + 3*y1*((2*t)*(t-1)+pow((1-t),2)) + 3*y2*(2*t*(1-t)-pow(t,2)) + 3*pow(t,2)*y3;
        float angle = atan(dy/dx);

        glPushMatrix();
        glTranslatef(x, y, 0);
        glRotatef(angle * 180 / 3.14159265, 0, 0, 1);
        drawRightArrow();
        glPopMatrix();
    }
}

But as you can see the tangents appear to be incorrect especially in the middle of a bezier curve?

Answer 1

Because we don't want to interrupt the line strip, it is reasonable to make another loop that draws the arrows. In this loop, we can skip some steps, because we probably don't want arrows after each step. The arrows can be drawn as a 'GL_LINES' primitive.

For readability I'd recommend to define the parameter t in the inner loop as

float t = (float)k/NLINESEGMENT;

Now we need to calculate the derivative of the curve with respect to t at this point. This derivative can be calculated for each coordinate independently. It looks like this question is a homework, so I'll leave that to you

float dx = ... //derivative of x-component with respect to t
float dy = ... //derivative of y-component with respect to t

We will also need the curve point. Ideally you have saved it in the previous loop in an array or similar structure.

So we can draw the arrow's base line (where s is a custom scale factor):

glVertex2d(x, y);
glVertex2d(x + s * dx, y + s * dy);

We will also need the actual arrow. An orthogonal vector to the arrow direction is (-dy, dx). So we can just combine the direction and the orthogonal direction to get the arrow:

glVertex2d(x + s * dx, y + s * dy);
glVertex2d(x + 0.9 * s * dx - 0.1 * dy, y + 0.9 * s * dy + 0.1 * dx);
glVertex2d(x + s * dx, y + s * dy);
glVertex2d(x + 0.9 * s * dx + 0.1 * dy, y + 0.9 * s * dy - 0.1 * dx);

This will result in a 90° arrow. You can change the angle by adjusting the coefficients (0.9 and 0.1).

Update

Your derivative is closer to the actual solution, but still contains some mistakes:

Furthermore, when calculating the angle, use the atan2 function. This allows you to get angles greater than PI/2:

float angle = atan2(dy,dx);