How to automatically have ContourLabels in middle of Contour lines if specific levels are not specified?

Question

I looked at all the questions/answers relating to placement of ContourLabels but could not find an answer to this. Some of the answers before are based on user giving before-hand, the Contours levels to use, using the option Contours -> {v1, v2, ....}.

But if one just want to use Contours->number, then it gets hard to figure how to solve this.

The labels that show up are not in the nice locations and hard to read. I'd like them to be in the middle of the contour lines. This is how Matlab actually does it by default.

Here is a MWE, and I also show the Matlab result.

Clear[x,y]
f = (11 - x - y)^2 + (1 + x + 10 y - x y)^2;
ContourPlot[f, {x, 0, 20}, {y, 0, 15}, 
 Contours -> 10, ContourShading -> None,
 ContourLabels -> Function[{x, y, z},
   Text[z, {x, y}]], AspectRatio -> Automatic, PlotRangePadding -> 2]

I am trying to get closer output to Matlab's

clear; close all;
x = 0:0.05:20;
y = 0:0.05:15;
[x,y] = meshgrid(x,y);
z = (11-x-y).^2 + (1+x+10*y-x.*y).^2;
[C,h] =contour(x,y,z,10);
clabel(C,h)

I also noticed that Matlab seems to have additional heuristics in placing labels. It does not label all the lines if the space is tight. Notice the top levels are not labeled in Matlab. Also the text is not written over the lines, but there is space in the line where the text goes. This makes it easier to read.

But my main question I have is how to put the labels in the middle of the level lines, where it is more pleasing to see than on the edges.

If someone knows of a package to use for this, I will try it. I looked at few here http://packagedata.net/ but did not see anything yet. Looked at sciDraw package also. I think Contour label placement in Mathematica should be improved more.

Version 10.3.1. windows.

Related questions

How to control the positions of contour labels?

Custom contour labels in ContourPlot

How do I add contour labels to contour plot?

Added:

FYI, I found the book Mathematica graphics which has an function called LabelContourLines[], which from the example shown, seems to do what I want. Page 348 in the book:

Here is how Mathematica 10.3.1. gives, for the same example in the book. Notice the difference :

ContourPlot[x y, {x, -2, 2}, {y, -2, 2}, ContourShading -> False, 
 Contours -> 10, ContourLabels -> All]

I do not have this book around myself, but will try to see if the author has the software on the net to use or I can look for it at the university library.

It looks like the above function is what I am looking for. But the book is old (1994), and do not know if book software will even work on current Mathematica even if I find it.

Answer 1

This is a good question which I had tried to work on several years ago, and I just revisited it to see if the code could be useful.

I had initially only worked on rotating the contour labels to be tangential with the contours, but in the code it turned out to be quite easy to add a line that also attempts to center the labels along the contour. My simple approach is to estimate the center of the contour by just choosing the midpoint in the Line specification defining that particular contour.

The code comes with comments, but it's quite lengthy:

Clear[rotateContourLabels];
rotateContourLabels::usage = 
  "The function rotateContourLabels accepts the output of ContourPlot \
or ListContourPlot, assuming they were made with the option \
ContourLabel-> All. The plot is passed via the required first \
argument. The function rotates the labels of all contours to be \
approximately parallel to the iso-lines. The optional LabelFunction \
can specify a custom label style in the form of a function f[#1, #2, \
#3], where {#1, #2} is the 2D location vector of the label, and #3 \
the value of the plotted function at that location. For examples see \
the documentation on ContourLabels. The default label style is given \
by the function Text[#3,{#1,#2}]&. A second option, Alignment, \
influences the placement of the labels along the contour. The default \
Alignment -> Automatic leaves the original placement intact, any \
other value will shift the labels to an estimated center position \
along the contour.";

Options[rotateContourLabels] = {LabelingFunction -> (Text[#3, {#1, \
#2}] &), Alignment -> Automatic};

rotateContourLabels[plot_, OptionsPattern[]] := 
 Module[{gc = plot[[1]] , valueList, pointList, rotatePoint, 
   labelCoordinatesIndexed, directionVector},
  (******* Created by Jens U. Nöckel, August 17, 2009, 
  last modified May 16, 2017  *******)
  Catch[
   labelCoordinatesIndexed = 
    Quiet@Check[gc[[-1, -1, All, 2]], Throw[plot]];
   valueList = Quiet@Check[gc[[-1, -1, All, 1]], Throw[plot]];
   If[Last[Dimensions[labelCoordinatesIndexed]]==2, 
     pointList = labelCoordinatesIndexed, 
     pointList = Quiet@Check[gc[[1, labelCoordinatesIndexed]],
      Throw[{"The plot must have the option ContourLabels->All !!!", 
        plot}]]];
   (** gc is the GraphicComplex, gc[[1]] a list of 2D points. 
   The list of contour labels is contained in plot[[1,-1,-1]] = 
   gc[[-1,-1]] and we want to replace this list. **)

   (* Local function definition: *)

   rotatePoint[labelPoint_] := 
    Module[{pt, extractSecant , allContourLines, contourPoints2D, 
      contourPointsIndexed , closestContourIndexed, 
      closestPointIndexed},
     pt = labelPoint;
     allContourLines = 
      Map[First,(** 
       The First extracts the list of line points from the Line \
expressions collected in the next line: **)

       Cases[Flatten@
         MapAll[If[SameQ[Head[#], Line], #, Apply[List, #]] &, 
          gc[[-1]]], _Line]];
     (** Line can be nested inside other non-
     list expressions such as Tooltip, 
     so the above flattens all levels of the expression except for \
subexpressions with head Line, 
     before it proceeds to collect the Lines using Cases. **) 

     {contourPoints2D, contourPointsIndexed} =
      (** Given a point pt, 
      find out for each contour what is its closest point to pt. 
      Return this information both as 2D points and as a list of \
indices counted within each contour. **)
      Transpose[
       Map[
        First@Nearest[
           # -> Transpose[{#, Range[Length[#]]}],
           pt
           ] (** End First@Nearest **)&,
        Map[gc[[1, #]] &, allContourLines]
        ](* End Map *)
       ] (* End Transpose *);

     {closestContourIndexed, closestPointIndexed} =
      (** Given a point pt, find closest contour. 
      Give its index in the list of contours, 
      and the index of the closest point within this contour. **)

        First[Nearest[
        contourPoints2D -> 
         Transpose[{Range[Length[contourPoints2D]], 
           contourPointsIndexed}],
        pt
        ]];(* End First@Nearest *)

     If[OptionValue[Alignment] =!= Automatic,
      (* Alignment refers to the placement of labels on contours. 
      If not Automatic, 
      try to approximately center the labels according to contour \
length, measured by points in their line: *)
      closestPointIndexed = 
       Floor[Length[allContourLines[[closestContourIndexed]]]/2]; 
      pt = Part[gc[[1]], 
        allContourLines[[closestContourIndexed, 
          closestPointIndexed]]]];
     extractSecant = (** Given the closest point to pt, 
      record it and its next neighbor on the same contour line. 
      The two points are returned as real 2D points, 
      determined from their indices. **)
      {
       allContourLines[[closestContourIndexed, closestPointIndexed]] (* 
       The first point *),
       Quiet@
        Check[allContourLines[[closestContourIndexed, 
           closestPointIndexed + 1]], 
         allContourLines[[closestContourIndexed, 
           closestPointIndexed - 1]]] (* 
       The second sequential point - 
       here I want to avoid falling outside the range of the contour \
point list. *)
       };

     (* The rotation aims to make this vector parallel to the \
contour. We approximate the contour direction by taking two points on \
the contour and forming their difference: *)

     directionVector = Apply[Subtract,
       (** The two required points are found from the list gc[[1]] : **)

              Part[gc[[1]],
        (* 
        The Part of gc[[1]] we're looking for is a set of two points \
that define the secant closest to the contour label: *)

        extractSecant] 
       ];

     (** Here starts the body of rotatePoint: **)
     Composition[
      (* The composition sandwiches a rotation around the origin \
between a translation to the origin and its inverse: *)

      TranslationTransform[pt],
      Quiet@Check[RotationTransform[{
          {1, 0},  (* 
          The reference direction for the rotation is the horizontal: \
*)
          (Sign[directionVector[[1]]] directionVector)
          (* 
          To avoid rotations outside [-Pi,Pi] the above \
statement places the direction vector into the first or fourth \
quadrant. *) 
          }] (* End RotationTransform *)
        , Identity (* 
        No rotation happens if for some reason we can't get a valid \
direction vector *)],
      TranslationTransform[-labelPoint]
      ] (* End Composition *) 
     ] ;(** End function rotatePoint **)

   (**** Main function body: *****)

   ReplacePart[plot, {1, -1, -1} -> 
     (** the replacement is a geometric transformation of the \
original list gc[[-1,-1]]: **)

     MapThread[GeometricTransformation,
      {MapThread[OptionValue[LabelingFunction], 
        Append[Transpose[pointList], valueList]],
       (** 
       What follows is the list of transformations for each contour \
label: **)
       Map[rotatePoint, pointList]}]]
   ]]

This defines a function rotateContourLabels which takes an existing plot and post-processes it. As always with post-processing, I have to make some assumptions about how the input was prepared:

In the ContourPlot or ListContourPlot that serves as the input, you must set either the option ContourLabels -> True or ContourLabels -> All.

With this assumption, you can then use the function as follows:

a = ContourPlot[Im[(x + I y)^(1/2)], {x, -1, 1}, {y, -1, 1},
   Contours -> 20, ContourLabels -> True]

This is the original output in Mathematica 10.2.

Now the usage example:

rotateContourLabels[a, Alignment -> Center, 
 LabelingFunction -> (Text[#3, {#1, #2}, 
     Background -> Directive[Opacity[.5], White]] &)]

I think this result is an improvement (given that it's fully automatic), although of course the centering is not perfect.

Edit

For completeness, here is the OP's minimum working example:

Clear[x, y]
f = (11 - x - y)^2 + (1 + x + 10 y - x y)^2;
b = ContourPlot[f, {x, 0, 20}, {y, 0, 15}, Contours -> 7, 
   ContourShading -> None, ContourLabels -> True, 
   AspectRatio -> Automatic, PlotRangePadding -> 2];

rotateContourLabels[b, Alignment -> Center, 
 LabelingFunction -> Function[{x, y, z}, Text[z, {x, y}]]]

Note that I added the label function as an option to rotateContourLabels instead of the original plot (which I called b). I also decreased the number of contours to 7 so the labels are more readable. The placement of the labels is obviously much better than the original. The background I chose in the first example would probably also make this plot look a little better.

Finally, here is my automatic version of the book example in the question:

c = ContourPlot[x y, {x, -2, 2}, {y, -2, 2}, ContourShading -> False, 
   Contours -> 10, ContourLabels -> All];

rotateContourLabels[c, Alignment -> Center, 
 LabelingFunction -> (Text[#3, {#1, #2}, Background -> White] &)]

One could probably also add an option to place the labels at the ends of the contours instead of the center, but I haven't done that yet (because it usually looks ugly, even though Mathematica seems to like doing that by default...).

Answer 2

The software is located at http://library.wolfram.com/infocenter/Books/3753/ as a .tar and a .zip file.

Update

How about a shameless lifting of one of the references you mention and adding a Rotate function?

ContourPlot[x y, {x, -2, 2}, {y, -2, 2}, ContourShading -> False, 
 Contours -> 10, 
 ContourLabels -> (Rotate[Text[" " <> ToString[#3] <> " ", {#1, #2}, 
      Background -> White], -ArcTan[#2/#1]] &)]

For other curves one would need to determine the tangent to the curve at point {#1, #2}.

Answer 3

I was unable to run the function written by Jens in Mathematica 12.0 cause of multiple errors caused by an altered (?) ContourPlot's realisation. So I have rewrote rotateContourLabels function to adopt it to a newer Mathatica version. For my own purposes I use the following realisation:

Clear[rotateContourLabels];

rotateContourLabels::usage = 
"The function rotateContourLabels accepts the output of ContourPlot \
or ListContourPlot, assuming they were made with the option \
ContourLabel-> All. The plot is passed via the required first \
argument. The function rotates the labels of all contours to be \
approximately parallel to the iso-lines. The optional LabelFunction \
can specify a custom label style in the form of a function f[#1, #2, \
#3], where {#1, #2} is the 2D location vector of the label, and #3 \
the value of the plotted function at that location. For examples see \
the documentation on ContourLabels. The default label style is given \
by the function Text[#3,{#1,#2}]&. A second option, Alignment, \
influences the placement of the labels along the contour. The default \
Alignment -> Automatic leaves the original placement intact, any \
other value will shift the labels to an estimated center position \
along the contour.";

Options[rotateContourLabels] = {
    LabelingFunction -> (Text[#3, {#1, #2}]&), 
    Alignment        -> Automatic
};

rotateContourLabels[plot_, OptionsPattern[]] := With[{
    gcIdx = Sequence[1,1], (* index of a GraphicsComplex inside a Graphics *)
    ppIdx = Sequence[1  ], (* index of a list of 2d-points in the GC *)
    ccIdx = Sequence[2,2], (* index of a list of contour lines in the GC *)
    llIdx = Sequence[2,3], (* index of a list of labels in the GC *)
    case  = Identity
}, Block[{gc = plot[[gcIdx]], 
    pointList, contourList, labelList, contourIdxs, contourPoss, labelValues, labelIdxs, 
    labelPoints, getClosestLine, getClosestPoints, makeTransforms, origins, transforms,
    plotRange, plotRatio, plotScale, computeLabel, rangeRatio
}, Catch[
    (*-- extract basic data containters --*)
    pointList   = Quiet@Check[gc[[ppIdx]], plot]; (* {{x,y} ...} *)
    contourList = Quiet@Check[gc[[ccIdx]], plot]; (* {{<style ...>, Line[List[..]]} ...} *)
    labelList   = Quiet@Check[gc[[llIdx]], plot]; (* {Inset[val, idx] ...} *)
    plotRange   = First /@ Differences /@ (PlotRange /. Options[plot]);
    rangeRatio  = plotRange[[2]] / plotRange[[1]];
    plotRatio   = AspectRatio /. Options[plot];
    
    (*-- normalize extracted data --*)
    contourIdxs = Cases[contourList, Line[idxs_] :> idxs, \[Infinity]]; (* {{point_idx ...} ...} *)
    contourPoss = Map[pointList[[#]]&, contourIdxs, {2}];
    labelValues = labelList[[;;, 1]];
    labelIdxs   = labelList[[;;, 2]];
    labelPoints = pointList[[labelIdxs]];
    If [Length@contourPoss == 0, 
        Return@plot
    ];
    
(*----> helper functions --*)
    (** Lablel rotation angle **) 
    getLabelAngle = Function[{baseVector, targetVector}, With[{
        (* Map targetVector to [-\[Pi], \[Pi]] angle *)
        base = baseVector   Sign[baseVector  [[1]]],
        dest = targetVector Sign[targetVector[[1]]]
    }, With[{
        \[Alpha] = VectorAngle[base, dest] Sign[(dest-base)[[2]]]
    },  
        (**   \[Alpha] \[LongDash] unscaled angle = arctan[y_0 / x_0]
         *    \[Beta] \[LongDash] scaled angle   = arctan[y_1 / x_1] = arctan[y_1 / x_0]
         *    y_1 = ar_target x_1 = ar_target x_0
         *    y_0 = ar_range  x_0
        ** => \[Beta] = arctan[ar_target/ar_range tan[\[Alpha]]] **)
        ArcTan[plotRatio/rangeRatio Tan[\[Alpha]]]
    ]]];
    (** find positions of the closest contour and the line's closest point **)
    getClosestLine = Function[{idxOrPos}, Block[{}, 
        Switch[Length @ idxOrPos
            , case @ 0, Block[{idx = idxOrPos, checks, positions}, 
                getClosestLine[pointList[[idx]]]
            ]
            , case @ 2, Block[{pos = idxOrPos, 
                closestPoint, closestPoints, closestLinePos, closestPointPos
            },
                closestPoints   = Map[First@Nearest[#1, pos]&, contourPoss];
                closestPoint    = First@Nearest[closestPoints, pos];
                closestLinePos  = Position[closestPoints, closestPoint][[1,1]];
                closestPointPos = Position[contourPoss[[closestLinePos]], closestPoint][[1,1]];
                {closestLinePos, closestPointPos}
            ]
            , _, Throw[{"unsupported index ot position type ", idxOrPos}]
        ]
    ]];
    (** the function finds points linked to a specified point **)
    getClosestPoints = Function[{linePos, pointPos}, Block[{idxs,
        lineLen = Length@contourIdxs[[linePos]]
    },
        idxs = With[{alignment = OptionValue[Alignment]}, Switch[alignment
            , Automatic, {-1, 0, 1} + pointPos
            , Center   , {-1, 0, 1} + \[LeftFloor]lineLen/2\[RightFloor]
            , _, Throw[{"Unsupported Alignment type ", OptionValue[Alignment]}]
        ]];
        idxs = idxs /. a_?(# < 1 || # > lineLen&) -> None;
        idxs
    ]];
    (** the function calclates a geometry transformation and a new label origin **)
    makeTransforms = Function[{pos}, Block[{
        idxOr2DPos, linePos, pointPos, basePointPoss, originPoint,
        basePointIdxs, secantVector, pointCoords, basePoints, transform
    },
        (* >> get closest points to a target label anchor *)
        {linePos, pointPos} = getClosestLine[labelIdxs[[pos]]];
        basePointPoss = getClosestPoints[linePos, pointPos];
        basePointIdxs = If[# =!= None, contourIdxs[[linePos, #]], None]& /@ basePointPoss;
        basePoints    = If[# =!= None, pointList[[#]]]& /@ basePointIdxs;
        originPoint   = basePoints[[2]];
        basePoints    = DeleteCases[basePoints, None];
        Assert[Length@basePoints >= 2];
        
        (* >> compute an right-direction-vector *)
        secantVector = basePoints[[-1]] - basePoints[[1]];
        
        (* >> composite a rotation operator *)
        With[{
            angle  = getLabelAngle[{1, 0}, secantVector],
            origin = originPoint
        },
            transform = Rotate[#, angle]&;
        ];
        {originPoint, transform}
    ]];
    computeLabel = Function[{x, y, labelValue, transform}, Block[{label},
        label = Replace[OptionValue[LabelingFunction][x, y, labelValue]
            , { 
                Text[v_, p_List, o_List, r___] :> Text[v, r],
                Text[v_, p_List, r___] :> Text[v, r],
                Text[v_, r___] :> Text[v, r]
            }
        ];
        Text[transform @ label, {x, y}]
    ]];
(*----< helper functions --*)
    
    (*-- main: transform original plot's labels --*)
    {origins, transforms} = Map[makeTransforms, Range[Length@labelPoints]]\[Transpose];
    ReplacePart[plot, {gcIdx, llIdx} -> MapThread[computeLabel
        , Join[origins\[Transpose], {labelValues, transforms}]
    ]]
]]]

As a result it's possible to get the same result on the presented examples like

c = ContourPlot[x y, {x, -2, 2}, {y, -2, 2}
    , ContourShading -> False
    , Contours -> 10
    , ContourLabels -> All
    , Background -> White
];

rotateContourLabels[c
    , Alignment -> Center
    , LabelingFunction -> (Text[#3, {#1, #2}, Background -> White] &)
]

and

a = ContourPlot[Im[(x + I y)^(1/2)], {x, -1, 1}, {y, -1, 1}
    , Contours -> 20
    , ContourLabels -> True
    , Background->White
];
rotateContourLabels[a
    , Alignment -> Center
    , LabelingFunction -> (Text[#3, {#1, #2}
        , Background -> Directive[Opacity[.5], White]
    ] &)
]

Edited:

Also, Mathemarica tries to adjust plot's aspect ratio to a specified value (Option: AspectRatio). The adjustment shears all plot contents, including all rotated lables (Text[..] filds are not sheared). To prevent labels from transformations I've started to wrap labels into a Text[, {position}] command, where can be Rotate[Text["lebel"], angle]. The artifact was clearly visible in plots where max[x_y_range] / min[x_y_range] >> 1

With[{
    equation = Subscript[\[Sigma], sb] (Subscript[\[Epsilon], w] Subscript[T, w]^4 - Subscript[\[Epsilon], 0] Subscript[T, 0]^4) + h (Subscript[T, w] - Subscript[T, 0]) /. {
        Subscript[\[Sigma], sb] -> 5.670373 10^-8,
        Subscript[T, 0]  -> 273,
        Subscript[\[Epsilon], 0]  -> 0.95,
        Subscript[\[Epsilon], w]  -> 0.35
    },
    vars = Sequence[{h, 2, 10}, {Subscript[T, w], 273, 273 + 100}]
},
    d = ContourPlot[equation, vars
        , ContourLabels  -> All
        , Contours       -> 10
        , Background     -> White
    ];
    d = rotateContourLabels[d
        , Alignment -> Center
        , LabelingFunction -> (Text[" " <> ToString@#3 <> " ", {#1, #2}
            , Background -> Directive[Opacity[.8], White]
        ]&)
    ]
]

Answer 4

JimB's answer is great, I wish I'd found it years ago. That said, sometimes I have data where the contours are very irregular, and I don't want a break in the contour line. Below is a simple version of the solution I use. This also gives control on label placement, one of the original concerns.

Considering the graphic used above;

cPlot = ContourPlot[x y, {x, -2, 2},{y, -2, 2}, ContourShading->False, Contours->10]

Then utilizing the Graphics option in Arrowheads

Arrowheads[{spec1, spec2, ...}]

Where each spec is {size, position, graphic}

Then use Replace (/.) to swap Line heads for Arrow, and prepend the contents of Tooltip with the desired Arrowheads, taking the Tooltip text to be the value used inside.

cPlot /. {Line@x___ :> Arrow@x} /. {Tooltip[{y___}, z_] :> Tooltip[{Arrowheads[{{Automatic, 0.5, Graphics@Text[ToString[z] <> "\n"]}}], y}, z]}

The 0.5 inside of Arrowheads places the label in the middle of each Line. Additional labels can be added by doing something like the following;

Arrowheads[{{s1, pos1, graphic1}, {s2, pos2, graphic2}, ...}]