differential equations - How to plot the derivative of NDSolve interpolation function?

So I got a question. Everytime I use a derivative of a function during plotting or ongoing evaluation, it seems Mathematica misunderstands me, and just derives a value (a constant) which always results in a 0. What could be the problem? I demonstrate this with the solution of the Blassius boundary layer equation:

sol1 = NDSolve[{
    D[f1[ξ], ξ] == g1[ξ], D[g1[ξ], ξ] == h1[ξ],
    D[h1[ξ], ξ] + f1[ξ]*h1[ξ] == 0, 
    f1[0] == 0, g1[0] == 1, g1[10] == 0},

   {f1, g1, h1},
   {ξ, 0, 10}];

f1rand = f1[10] /. sol1;
g1rand = g1[10] /. sol1;
h1rand = h1[10] /. sol1;

Plot[D[g1[ξ], ξ] /. sol1, {ξ, 0, 10}, PlotRange -> All]
Plot[h1[ξ] /. sol1, {ξ, 0, 10}, PlotRange -> All]

Any help would be greatly appreciated!

Answer

The correct way to do it is to use g1'[ξ]

 Plot[g1'[ξ] /. sol1, {ξ, 0, 10}, PlotRange -> All]

Plotting f1'[ξ], g1'[ξ] and h1'[ξ] in one plot,

Plot[Evaluate[{f1'[ξ], g1'[ξ], h1'[ξ]} /. sol1], {ξ,0, 10}, PlotRange -> All, 
PlotStyle -> {Red, Green, Blue},Frame -> True, PlotLegends -> {"f1'", "g1'", "h1'"}]