source: Sophya/trunk/SophyaPI/DemoPIApp/resonance.pic

################################################################
### SOPHYA class library - piapp interactive analysis tool 
### (C) LAL / Univ. Paris Sud, CNRS/IN2P3    1998-2012 
### (C) IRFU-SPP / CEA                       1998-2012 
################################################################
#  Cours S1SM - Oscillateur harmonique - Courbe de resonance 
#          Variation en fonction de la courbe d'amortissement
#                                  Reza Ansari - Janvier 2003 

#  ...........................................................
#  Script Calcul/trace frequence / gain de resonance
# Frequence OH non amorti $1 : Omega0
# Amortissement $2 : Lambda
# Force exterieur $3 : F0/m
# $4 : graphic_att
# $5 : if $5 < 0 --> Pas de trace, calcul omrlast et gainreslast uniquement
# 

defscript omegares ' Calcul frequence / gain de resonance '
  rpneval omres  $1 sq 2 $2 sq * - sqrt
  rpneval gainres  2. $2 * $1 sq $2 sq - sqrt * $3 x<>y / 
  echo "OmegaRes= $omres  GainRes= $gainres"
  if ( $5 > -1 ) then  
    addmarker $omres $gainres "$4"
    addline $omrlast $gainreslast  $omres $gainres 'black dottedline'
  endif
  set omrlast $omres
  set gainreslast $gainres
endscript 
#  ...........................................................

#  ...........................................................
#  Script trace de courbe de Gain = f (omega)
# Frequence OH non amorti $1 : Omega0
# Amortissement $2 : Lambda
# Force exterieur $3 : F0/m
# $4 $5: Limite Min,Max en Omega 
# $6: graphic_att

defscript plotgain ' Trace de courbe de resonance - Gain '
  func $3/sqrt(($1*$1-x*x)*($1*$1-x*x)+4*$2*$2*x*x) $4 $5 200 "$6" 
endscript 
#  ...........................................................

#  ...........................................................
#  Script trace de courbe de Phase = f (omega)
# Frequence OH non amorti $1 : Omega0
# Amortissement $2 : Lambda
# Force exterieur $3 : F0/m
# $4 $5: Limite Min,Max en Omega 
# $6: graphic_att

defscript plotphase
  func atan2(-2*$2*x,$1*$1-x*x) $4 $5 200 "$6" 
endscript 
#  ...........................................................


#  Choix des parametres   
om0 = 1.
lambda = 0.45
f0om = 3.
omegares $om0 $lambda $f0om "boxmarker5 $col" -1

#  Choix de la liste de couleurs et des valeurs d'amortissement
set colors ( red siennared orange yellow green blue blueviolet )
set lambdas ( 0.4 0.3   0.2   0.15  0.1  0.075 0.05 )
setaxesatt 'minorticks timesfont bolditalicfont fixedfontsize'

# echo "DBG222 ------ " 

#  On trace courbe de gain et courbe de phase l'une sur l'autre 
zone 1 2
gratt = 'notitle xylimits=0.,2.2,0.25,32. '
for i 0:7
  set col $colors[i]
  set gratt "$gratt $col"
  set lambda $lambdas[i]
  plotgain $om0 $lambda $f0om 0. 3. "$gratt"
  omegares $om0 $lambda $f0om "boxmarker5 $col" $i
  if ( $i < 1 ) then 
    setaxelabels 'Frequence' 'Gain' 'timesfont bigbolditalicfont'
  endif
  addctext 1.6 $gainres 'l/w' '' '0' "symbolfont bigfont horizright vertcenter $col"
  rap = $lambda/$om0  
  addtext 1.65 $gainres "= $rap"  "timesfont bigboldfont horizleft vertcenter $col"
  gratt = 'same '
end


gratt = 'notitle xylimits=0.,2.2,-3.2,0.2 '
for i 0:7
  set col $colors[i]
  set gratt "$gratt $col"
  set lambda $lambdas[i]
  plotphase $om0 $lambda $f0om 0. 3. "$gratt"
  if ( $i < 1 ) then 
    setaxelabels 'Frequence' 'Phase' 'timesfont bigbolditalicfont'
  endif
  yt = -0.2-0.3*$i 
  addctext 1.5 $yt 'l/w' '' '0' "symbolfont bigfont horizright vertcenter $col"
  rap = $lambda/$om0  
  addtext 1.55 $yt "= $rap"  "timesfont bigboldfont horizleft vertcenter $col"
  gratt = 'same '
end

w2eps resonance.eps