For a contour plot use,

ListContourPlot[Flatten[Table[
    s = NDSolveValue[{D[xtraj[t], t] == 
        (Sinh[2 xtraj[t] (t - 2)])/(Cosh[2 xtraj[t] (t - 2)]), 
        xtraj[0] == n}, xtraj, {t, 0, 4}];
    Table[{s[t], t, n}, {t, 0, 4, .2}], {n, -3.2, 3.2, .2}], 1], 
 Contours -> Range[-4, 4], ImageSize -> Large, PlotRange -> {{-3, 3}, Automatic}]

The code for a density plot is similar.

For a contour plot use,

ListContourPlot[Flatten[Table[
    s = NDSolveValue[{D[xtraj[t], t] == 
        (Sinh[2 xtraj[t] (t - 2)])/(Cosh[2 xtraj[t] (t - 2)]), 
        xtraj[0] == n}, xtraj, {t, 0, 4}];
    Table[{s[t], t, n}, {t, 0, 4, .2}], {n, -3.2, 3.2, .2}], 1], 
 Contours -> Range[-4, 4], ImageSize -> Large, PlotRange -> {{-3, 3}, Automatic},
 AspectRatio -> 2/3, ColorFunction -> "Rainbow"]

The code for a density plot is similar.

ListDensityPlot[Flatten[Table[
    s = NDSolveValue[{D[xtraj[t], t] == 
        (Sinh[2 xtraj[t] (t - 2)])/(Cosh[2 xtraj[t] (t - 2)]), 
        xtraj[0] == n}, xtraj, {t, 0, 4}];
    Table[{s[t], t, n}, {t, 0, 4, .2}], {n, -3, 3, .2}], 1], 
ImageSize -> Large, AspectRatio -> 2/3, ColorFunction -> "Rainbow"]