In MATLAB there's a function called phase that is just like our Listable friend Arg, but that, when evaluated with a list, it tries to keep the result "continuous", allowing itself to return values outside the [-Pi, Pi] range. So, for example, sampling a complex exponential

vec = E^(I 2 π #/10) & /@ Range[100];
resMMA = Arg[vec] // N;

The red line is Mathematica's result with Arg, and the blue line is MATLAB's with phase

enter image fdfdhere

So, I want to implement a function that does a good job at allowing you to see phase properties such as "linearity", more easily, like phase, in Mathematica. Ideally, it should also work on some symbolic functions and not only lists, but what's more important to me is for it to work as good as possible on lists.

My attempt at the "list version" is very simple, but seems to be working reasonably, and fast. I also tested it with a function that's positive and negative, like vec2 = E^(I 2 Pi/16. #) Sinc[2 Pi #/20.] & /@ Range[100]; Here's the code

phase = Compile[{{l, _Complex, 1}}, 
    Function[{prev, new}, # + Round[prev - #, 2 Pi] &@Arg@new], 
    Arg@First@l, Rest@l], 
   CompilationTarget -> "C", RuntimeOptions -> "Speed"];

Any improvements, would be appreciated. For example, it would probably be better if it could look at more than just the nearest neighbouring point, to be more robust to an occasional noise. I'd also appreciate ideas on how to achieve something analogous for a continuous well-behaved, not-created-by-mathematicians-just-to-prove-people-wrong symbolic function?

  • $\begingroup$ Indeed, very useful. I wonder how to implement similar thing for Plot3D. It'd be needed for example when solving a Schrodinger equation and plotting Arg of wave-function versus time and space. $\endgroup$ – Vitaliy Kaurov May 19 '12 at 9:16
  • 8
    $\begingroup$ For completeness, here is a straightforward Mathematica translation of the code used within MATLAB: phase[vec_?VectorQ] := Module[{n = Length[vec], ph = Arg[vec], df}, df = Differences[ph]; Do[ph -= 2 Pi Sign[df[[k]]] PadRight[ConstantArray[0, k], n, 1], {k, Flatten[Position[Abs[df], _?(# > Pi &)]]}]; ph] $\endgroup$ – J. M.'s technical difficulties May 19 '12 at 16:11
  • 2
    $\begingroup$ forums.wolfram.com/mathgroup/archive/1998/May/msg00107.html $\endgroup$ – Daniel Lichtblau May 20 '12 at 16:26

For a continuous function you could do something like this:

SetAttributes[argPlot, HoldAll];
Options[argPlot] = Options[Plot];

argPlot[exp_, {x_, x0_, x1_}, opt : OptionsPattern[argPlot]] :=
 Module[{pts, pl},
  pl = Plot[Arg[exp], {x, x0, x1}, PlotRange -> All, 
    PlotPoints -> OptionValue[PlotPoints]];
  pts = SortBy[Cases[pl, Line[pts_] :> pts, Infinity], #[[1, 1]] &];
  pts = Reap[Fold[Module[{ptsn},
         ptsn = #2;
         ptsn[[All, 2]] -= Round[ptsn[[1, 2]] - #1, 2 Pi];
         ptsn[[-1, 2]]] &, 0, pts];][[2, 1]];      
  ListLinePlot[Flatten[pts, 1], opt]]

argPlot[3 + 2 Exp[3 I a] + Exp[(1 - I a^2)], {a, -8, 8}]

Mathematica graphics

Compared to an ordinary plot of the Arg[3 + 2 Exp[3 I a] + Exp[(1 - I a^2)]]

Plot[Arg[3 + 2 Exp[3 I a] + Exp[(1 - I a^2)]], {a, -8, 8}, PlotRange -> All]

Mathematica graphics

|improve this answer|||||
  • $\begingroup$ Nice test function. $\endgroup$ – lirtosiast Aug 10 '19 at 7:46

This solution only considers direct differences, but is a bit shorter than the original (it doesn't contain compilation, though, and I didn't test the speed):

phase[l: {_?NumericQ ..}] :=
  Module[{args = Arg[l]},
         args+Prepend[2Pi Accumulate@-IntegerPart@Differences[args/Pi],0]]

The following solution definitely isn't the fastest, but generally should give the best possible result:

phase[l: {_?NumericQ ..}]:=
  With[{vars = Table[Unique[], {Length@l}]},
       Arg@l+2 Pi Accumulate@vars /.
         Last@NMinimize[{Total@Abs@Differences[Arg@l+2\[Pi] Accumulate@vars],
                         (Alternatives@@vars) \[Element] Integers &&
|improve this answer|||||
  • $\begingroup$ Just a side note, the following is the compiled version of the first function: Clear[phase]; phase = Compile[{{l, _Complex, 1}}, Module[{args = Arg[l]}, args + Prepend[2 Pi Accumulate@-IntegerPart@Differences[args/Pi], 0]], RuntimeOptions -> "EvaluateSymbolically" -> False]; $\endgroup$ – xzczd Jun 11 '18 at 19:01

We can transfer the phase function from Matlab to Mathematica:

phase[vec_List] := 
  Module[{phi, df, len, i},
    phi = Arg @ vec;
    df = Differences @ phi;
    len = Length @ phi;
    i = Flatten @ Position[df, x_ /; Abs[x] > 3.5];
    Do[phi = phi - (2 Pi*Sign[df[[j]]]*UnitStep[# - (j + 1)] & /@ Range[len]),
      {j, i}]; 

data = Table[1/((3/2 Exp[I a])^4 - 1), {a, -8, 8, 0.01}];
data // phase // ListLinePlot

The output is

enter image description here


data // Arg // ListLinePlot


enter image description here

|improve this answer|||||

Here is a simple answer for an ordinary list. The OP asked about:

vec = E^(I 2 \[Pi] #/10) & /@ Range[100];
resMMA = Arg[vec] // N;

Just a simple few functions is enough:

unwraped =   resMMA + Accumulate[2*Pi*Unitize[Threshold[Prepend[Differences[resMMA], 0.], 1]]];
ListLinePlot[{resMMA, unwraped}]

enter image description here

Sometimes you might need to subtract instead of add, if your saw-tooth is upside down.

|improve this answer|||||
  • 1
    $\begingroup$ I'd suggest to use Threshold instead of Chop because the former preserves packing of lists. And using 0. instead of 0. $\endgroup$ – Henrik Schumacher Aug 9 '19 at 18:12
  • 2
    $\begingroup$ @HenrikSchumacher Thanks! I read about this, and agree. Changes made $\endgroup$ – axsvl77 Aug 10 '19 at 0:24
  • 1
    $\begingroup$ This code doesn't work on this dataset (which I import by resMMA = Flatten@Import["/tmp/angles.csv"];). $\endgroup$ – Ruslan Jan 21 at 10:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.