-1

I would like to find the phase of a transfer function in python symbolically

I would like for example to find the phase of H_lp1 given below

from sympy import *
import cmath

import numpy as np

fr = Symbol('fr', real=True)

fo = Symbol('fo', real=True)
fn = Symbol('fn', real=True)
f = Symbol('f',real=True)
q = Symbol('q', real=True)
qn = Symbol('qn', real=True)

wr = 2*pi*fr
wo = 2*pi*fo
wn = 2*pi*fn
w = 2*pi*f
s = I*w


H_lp1 = wo/(s+wo)

Thanks

Adding a more complicated example

LPN - s domain

H_lpn2 = (s**2*(wo/wn)**2+wo**2)/(s**2+s*(wo/q)+wo**2)

(re_part,im_part) = H_lpn2.expand(complex=True).as_real_imag()

H_lpn2_ph_expanded = atan(im_part/re_part)

H_lpn2_ph = simplify(H_lpn2_ph_expanded)

produces atan(f*fo/(q*(f**2 - fo**2)))

This is incorrect as it dropped the term that contains the fn

The correct expression should be 

pi*sign(f**2-fn**2)/2 -atan(q*(f**2-fo**2)/(f*fo))
Lefti
  • 3
  • 5

2 Answers2

0

Normally is better to use Sympy own constants. I suppose that you expect fr and f to be real numbers, that would return a simpler expression. I would rewrite your example as:

from sympy import pi, Symbol, I, atan, im, re

# Defining the values as real numbers allows to simplify
# the resulting expression as im(fr) or re(fr)

fr = Symbol('fr', real=True)
f = Symbol('f', real=True)

wo = 2*pi*fr
w = 2*pi*f
s = I*w
H_lp1 = wo/(s+wo)
# The expresion has to be expanded first to have the I in the
# numerator
(re_part, im_part) = H_lp1.expand(complex=True).as_real_imag()
atan(im_part/ re_part)

That gives as an answer: 

-atan(f/fr)

Or alternative you can also: 

H_lp1 = H_lp1.expand(complex=True)
atan(im(H_lp1)/re(H_lp1))

With the same result

Juan Osorio
  • 378
  • 3
  • 15
  • This is good and works for the case above. When I try a more complicated example along with simplify it does not work. – Lefti Sep 01 '15 at 01:41
0

The result you get is correct:

atan(f*fo/(q*(f**2 - fo**2)))

the numerator of the expression is a real number so the angle is minus times the angle of the denominator and this is just:

-atan(w*wo/(q*(wo**2 - w**2))

equivalent to what you get from sympy!!.

Juan Osorio
  • 378
  • 3
  • 15