Refine `Plot` around y=1 only

Question

I have a complicated function $f$ which turns an input (let's say a single number, to simplify) into either 1, -1, -2 or -3, each value corresponding to different categories.

I can plot the categories as a function of the input, and that's a typical results:

What I am interested in is $f^{-1}(\{1\})$. For example, for the plot above it would be something like $[-11.1,-10.9]\cup [-9.7,7.6]\cup[-7.4,-0.5]\cup\dots$, but much more accurately (typically, with five digits).

One possibility would be to increase PlotPoints and/or MaxRecursions, but this will also refine the parts in which I am not intersted, i.e. the transitions between -1 and -2, etc.

---

Here is a simplified example; f is a very simple function given in the end, but assume it is not known explicitly.

plot = Plot[f[x], {x, -5, 5}, Exclusions -> None, PlotPoints -> 10, AspectRatio -> 1/5];

gives

The points can be extracted from the plot and the transitions obtained with BobHanlon's code:

transitions = Cases[plot, Line[pts_] :> pts, Infinity][[1]] //.
                  {s___, {x1_, y1_}, {x2_, y1_}, e___} :> {s, {x1, y1}, e}

(* {{-5., -2.}, {-3.13482, 1.}, {0.0167145, -1.}, {3.14997, -2.}} *)

Since I'm interested in the $x$ values of image $1$, I simply select:

pos = Position[transitions, _?(Last[#] == 1 &), {1}];
Take[transitions[[All, 1]], Flatten@{pos, pos + 1}]

(* {-3.13482, 0.0167145} *)

Which gives me an estimation of $f^{-1}(\{1\})$ of $[-3.13482, 0.0167145]$. In the present case, f is build up so that the exact transitions are known: the solution is -Pi,0. My goal is to get more accurate values, in an efficient manner.

One possibility is to increase PlotPoints and MaxRecusions: for example, same code with PlotPoints -> 2000, MaxRecursion -> 5:

 plot = Plot[f[x], {x, -5, 5}, Exclusions -> None, AspectRatio -> 1/5, 
   ImageSize -> 700, PlotPoints -> 2000, MaxRecursion -> 5];
transitions = 
 Cases[plot, Line[pts_] :> pts, 
    Infinity][[1]] //. {s___, {x1_, y1_}, {x2_, y1_}, 
    e___} :> {s, {x1, y1}, e}
pos = Position[transitions, _?(Last[#] == 1 &), {1}];
Take[transitions[[All, 1]], Flatten@{pos, pos + 1}]

(* {-3.14151, 0.0000912312} *)

It's much better... but I'm looking for a more efficient way to do this, in particular without refining the transitions which are not of interest. In other words, focus on making very fine transitions only when the functions reaches or leaves the value 1.

I am also interested in approaches which are not plot-based; I just want to get the intervals of image 1 as accurate as possible.

---

Function f used in the illustration code:

f[x_] = -UnitStep[x - Pi] + 3 UnitStep[x + Pi] - 2 UnitStep[x] - 2;

Answer 1

Plot has an suboption "ControlValue" which specify the max angle between successive points (in Radians)

f[x_] := -UnitStep[x - Pi] + 3 UnitStep[x + Pi] - 2 UnitStep[x] - 2;
    gr = Plot[f[x], {x, -5, 5},
      Method -> {"Refinement" -> {"ControlValue" -> 2 Pi/100}}, 
      Mesh -> All, MaxRecursion -> 15, PlotPoints -> 10, 
      MeshStyle -> PointSize[0.015]]

intervals = Cases[Normal[gr], Point[___], {1, Infinity}] //
      First /@ # & //
     Sort //
    Split[#, Abs[#1[[2]] - #2[[2]]] < 0.1 &] & //
   #[[All, {1, -1}, 1]] &;

intervals // Grid

error around -Pi :

-3.1415587451696543` - -3.141626562009932`

0.0000678168

There are a lot of related informations here and here

Answer 2

f[x_] = -UnitStep[x - 2] + 3 UnitStep[2 x + 4] - 2 UnitStep[x] - 2;

transitions = 
 Cases[Plot[f[x], {x, -5, 5}, Exclusions -> None], Line[pts_] :> pts, 
    Infinity][[1]] //.
  {s___, {x1_, y1_}, {x2_, y1_}, e___} :>
   {s, {x1, y1}, e}

(*  {{-5., -2.}, {-1.99783, 1.}, {0.000769648, -1.}, {2.00061, -2.}}  *)

ListStepPlot[transitions]

EDIT:

Transitions occur at discontinuities in function's derivative. Select transition pairs for which function has desired value between those transitions.

intervals[func_, val_: 1] :=
 Module[{x, transitions},
  transitions = 
   SortBy[x /. Solve[#, x][[1]] & /@ 
     Cases[func'[x], {Indeterminate, rel_} :> rel, Infinity], N[#] &];
  Select[
   Partition[transitions, 2, 1],
   func[Mean[#]] == val &]]

f[x_] = -UnitStep[x - 2] + 3 UnitStep[2 x + 4] - 2 UnitStep[x] - 2;

intervals[f]

(*  {{-2, 0}}  *)

f2[x_] = -UnitStep[x - Pi] + 3 UnitStep[x + Pi] - 2 UnitStep[x] - 2;

intervals[f2]

(*  {{-π, 0}}  *)