Qucs-core  0.0.19
spfile.cpp
Go to the documentation of this file.
00001 /*
00002  * spfile.cpp - S-parameter file class implementation
00003  *
00004  * Copyright (C) 2004, 2005, 2006, 2008, 2009 Stefan Jahn <stefan@lkcc.org>
00005  *
00006  * This is free software; you can redistribute it and/or modify
00007  * it under the terms of the GNU General Public License as published by
00008  * the Free Software Foundation; either version 2, or (at your option)
00009  * any later version.
00010  *
00011  * This software is distributed in the hope that it will be useful,
00012  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00014  * GNU General Public License for more details.
00015  *
00016  * You should have received a copy of the GNU General Public License
00017  * along with this package; see the file COPYING.  If not, write to
00018  * the Free Software Foundation, Inc., 51 Franklin Street - Fifth Floor,
00019  * Boston, MA 02110-1301, USA.
00020  *
00021  * $Id$
00022  *
00023  */
00024 
00025 #if HAVE_CONFIG_H
00026 # include <config.h>
00027 #endif
00028 
00029 #include "component.h"
00030 #include "matvec.h"
00031 #include "dataset.h"
00032 #include "strlist.h"
00033 #include "poly.h"
00034 #include "spline.h"
00035 #include "interpolator.h"
00036 #include "spfile.h"
00037 
00038 using namespace qucs;
00039 
00040 // Constructor for S-parameter file vector.
00041 spfile_vector::spfile_vector () {
00042   v = f = 0;
00043   isreal = 1;
00044   inter = NULL;
00045   r = c = 0;
00046 }
00047 
00048 // Destructor for S-parameter file vector.
00049 spfile_vector::~spfile_vector () {
00050   delete inter;
00051 }
00052 
00053 // Passes vectors and their data types to the S-parameter file vector.
00054 void spfile_vector::prepare (qucs::vector * _v, qucs::vector * _f,
00055                              bool _isreal, int it, int dt) {
00056   v = _v;
00057   f = _f;
00058   isreal = _isreal;
00059   inter = new interpolator ();
00060   if (isreal) {
00061     inter->rvectors (v, f);
00062     inter->prepare (it, REPEAT_NO, dt | DATA_REAL);
00063   }
00064   else {
00065     inter->cvectors (v, f);
00066     inter->prepare (it, REPEAT_NO, dt | DATA_COMPLEX);
00067   }
00068 }
00069 
00070 // Returns interpolated data.
00071 nr_complex_t spfile_vector::interpolate (nr_double_t x) {
00072   if (isreal)
00073     return inter->rinterpolate (x);
00074   else
00075     return inter->cinterpolate (x);
00076 }
00077 
00078 // Constructor creates an empty and unnamed instance of the spfile class.
00079 spfile::spfile () : circuit () {
00080   data = NULL;
00081   sfreq = nfreq = NULL;
00082   spara = FMIN = SOPT = RN = NULL;
00083   interpolType = dataType = 0;
00084   type = CIR_SPFILE;
00085   setVariableSized (true);
00086 }
00087 
00088 // Destructor deletes spfile object from memory.
00089 spfile::~spfile () {
00090   delete[] spara;
00091   delete RN;
00092   delete FMIN;
00093   delete SOPT;
00094 #if DEBUG && 0
00095   if (data) {
00096     data->setFile ("spfile.dat");
00097     data->print ();
00098   }
00099 #endif
00100   delete data;
00101 }
00102 
00103 void spfile::calcSP (nr_double_t frequency) {
00104 
00105   // nothing to do if the given file type had errors
00106   if (spara == NULL || sfreq == NULL) return;
00107 
00108   // set interpolated S-parameters
00109   setMatrixS (expandSParaMatrix (getInterpolMatrixS (frequency)));
00110 }
00111 
00112 /* This function returns the S-parameter matrix of the circuit for the
00113    given frequency.  It uses interpolation for frequency points which
00114    are not part of the original touchstone file. */
00115 matrix spfile::getInterpolMatrixS (nr_double_t frequency) {
00116 
00117   // first interpolate the matrix values
00118   matrix s (getSize () - 1);
00119   for (int r = 0; r < getSize () - 1; r++) {
00120     for (int c = 0; c < getSize () - 1; c++) {
00121       int i = r * getSize () + c;
00122       s.set (r, c, spara[i].interpolate (frequency));
00123     }
00124   }
00125 
00126   // then convert them to S-parameters if necessary
00127   switch (paraType) {
00128   case 'Y':
00129     s = ytos (s);
00130     break;
00131   case 'Z':
00132     s = ztos (s);
00133     break;
00134   case 'H':
00135     s = htos (s);
00136     break;
00137   case 'G':
00138     s = gtos (s);
00139     break;
00140   }
00141   return s;
00142 }
00143 
00144 void spfile::calcNoiseSP (nr_double_t frequency) {
00145   // nothing to do if the given file type had errors
00146   if (spara == NULL || nfreq == NULL) return;
00147   setMatrixN (calcMatrixCs (frequency));
00148 }
00149 
00150 matrix spfile::calcMatrixCs (nr_double_t frequency) {
00151   // set interpolated noise correlation matrix
00152   nr_double_t r = real (RN->interpolate (frequency));
00153   nr_double_t f = real (FMIN->interpolate (frequency));
00154   nr_complex_t g = SOPT->interpolate (frequency);
00155   matrix s = getInterpolMatrixS (frequency);
00156   matrix n = correlationMatrix (f, g, r, s);
00157   matrix c = expandNoiseMatrix (n, expandSParaMatrix (s));
00158   return c;
00159 }
00160 
00161 /* This function expands the actual S-parameter file data stored
00162    within the touchstone file to have an additional reference one-port
00163    whose S-parameter is -1 (i.e. ground). */
00164 matrix spfile::expandSParaMatrix (matrix s) {
00165   assert (s.getCols () == s.getRows ());
00166   int r, c, ports = s.getCols () + 1;
00167   nr_double_t g = -1;
00168   nr_complex_t fr, ss, sr, sc, sa;
00169   matrix res (ports);
00170 
00171   // compute S'mm
00172   for (sa = 0, r = 0; r < ports - 1; r++)
00173     for (c = 0; c < ports - 1; c++) sa += s.get (r, c);
00174   ss = (2 - g - ports + sa) / (1 - ports * g - sa);
00175   res.set (ports - 1, ports - 1, ss);
00176   fr = (1.0 - g * ss) / (1.0 - g);
00177 
00178   // compute S'im
00179   for (r = 0; r < ports - 1; r++) {
00180     for (sc = 0, c = 0; c < ports - 1; c++) sc += s.get (r, c);
00181     res.set (r, ports - 1, fr * (1.0 - sc));
00182   }
00183 
00184   // compute S'mj
00185   for (c = 0; c < ports - 1; c++) {
00186     for (sr = 0, r = 0; r < ports - 1; r++) sr += s.get (r, c);
00187     res.set (ports - 1, c, fr * (1.0 - sr));
00188   }
00189 
00190   // compute S'ij
00191   for (r = 0; r < ports - 1; r++) {
00192     for (c = 0; c < ports - 1; c++) {
00193       fr = g * res (r, ports - 1) * res (ports - 1, c) / (1.0 - g * ss);
00194       res.set (r, c, s.get (r, c) - fr);
00195     }
00196   }
00197 
00198   return res;
00199 }
00200 
00201 /* The function is the counterpart of the above expandSParaMatrix()
00202    function.  It shrinks the S-parameter matrix by removing the
00203    reference port. */
00204 matrix spfile::shrinkSParaMatrix (matrix s) {
00205   assert (s.getCols () == s.getRows () && s.getCols () > 0);
00206   int r, c, ports = s.getCols ();
00207   nr_double_t g = -1;
00208   matrix res (ports - 1);
00209 
00210   // compute S'ij
00211   for (r = 0; r < ports - 1; r++) {
00212     for (c = 0; c < ports - 1; c++) {
00213       res.set (r, c, s (r, c) + g * s (r, ports - 1)  *
00214                s (ports - 1, c) / (1.0 - g * s (ports - 1, ports - 1)));
00215     }
00216   }
00217   return res;
00218 }
00219 
00220 /* This function expands the actual noise correlation matrix to have an
00221    additional reference one-port whose S-parameter is -1
00222    (i.e. ground).  The given S-parameter matrix is required to perform
00223    this transformation and is obtained using the expandSParaMatrix()
00224    function. */
00225 matrix spfile::expandNoiseMatrix (matrix n, matrix s) {
00226   assert (s.getCols () == s.getRows () && n.getCols () == n.getRows () &&
00227           n.getCols () == s.getCols () - 1);
00228   nr_double_t T = getPropertyDouble ("Temp");
00229   int r, c, ports = n.getCols () + 1;
00230   nr_double_t g = -1;
00231 
00232   // create K matrix
00233   matrix k (ports, ports - 1);
00234   for (r = 0; r < ports - 1; r++) {
00235     for (c = 0; c < ports - 1; c++) {
00236       if (r == c)
00237         k.set (r, c, 1.0 + g * (s.get (r, ports - 1) - 1.0));
00238       else
00239         k.set (r, c, g * s.get (r, ports - 1));
00240     }
00241   }
00242   for (c = 0; c < ports - 1; c++)
00243     k.set (ports - 1, c, g * s.get (ports - 1, ports - 1) - 1.0);
00244 
00245   // create D vector
00246   matrix d (ports, 1);
00247   for (r = 0; r < ports - 1; r++) d.set (r, 0, s.get (r, ports - 1));
00248   d.set (ports - 1, 0, s.get (ports - 1, ports - 1) - 1.0);
00249 
00250   // expand noise correlation matrix
00251   matrix res (ports);
00252   res = (k * n * adjoint (k) - celsius2kelvin (T) / T0 * fabs (1 - norm (g)) *
00253          d * adjoint (d)) * norm (1 / (1 - g));
00254   return res;
00255 }
00256 
00257 /* The function is the counterpart of the above expandNoiseMatrix()
00258    function.  It shrinks the noise correlation matrix by removing the
00259    reference port.  The given S-parameter matrix is required to perform
00260    this transformation and is obtained using the expandSParaMatrix()
00261    function. */
00262 matrix spfile::shrinkNoiseMatrix (matrix n, matrix s) {
00263   assert (s.getCols () == s.getRows () && n.getCols () == n.getRows () &&
00264           n.getCols () == s.getCols () && n.getCols () > 0);
00265   int r, ports = n.getCols ();
00266   nr_double_t g = -1;
00267   nr_double_t T = getPropertyDouble ("Temp");
00268 
00269   // create K' matrix
00270   matrix k (ports - 1, ports);
00271   for (r = 0; r < ports - 1; r++) k.set (r, r, 1);
00272   for (r = 0; r < ports - 1; r++)
00273     k.set (r, ports - 1, g * s.get (r, ports - 1) /
00274            (1.0 - g * s.get (ports - 1, ports - 1)));
00275 
00276   // create D' vector
00277   matrix d (ports - 1, 1);
00278   for (r = 0; r < ports - 1; r++) d.set (r, 0, s.get (r, ports - 1));
00279 
00280   // shrink noise correlation matrix
00281   matrix res (ports - 1);
00282   res = k * n * adjoint (k) + celsius2kelvin (T) / T0 * fabs (1.0 - norm (g)) /
00283     norm (1.0 - g * s.get (ports - 1, ports - 1)) * d * adjoint (d);
00284   return res;
00285 }
00286 
00287 void spfile::prepare (void) {
00288 
00289   // check type of data
00290   const char * const dtype = getPropertyString ("Data");
00291   if (!strcmp (dtype, "rectangular")) {
00292     // rectangular data
00293     dataType = DATA_RECTANGULAR;
00294   }
00295   else if (!strcmp (dtype, "polar")) {
00296     // polar data
00297     dataType = DATA_POLAR;
00298   }
00299 
00300   // check type of interpolator
00301   const char * const itype = getPropertyString ("Interpolator");
00302   if (!strcmp (itype, "linear")) {
00303     interpolType = INTERPOL_LINEAR;
00304   }
00305   else if (!strcmp (itype, "cubic")) {
00306     interpolType = INTERPOL_CUBIC;
00307   }
00308 
00309   // load S-parameter file
00310   const char * file = getPropertyString ("File");
00311   if (data == NULL) data = dataset::load_touchstone (file);
00312   if (data != NULL) {
00313     // determine the number of ports defined by that file
00314     int ports = (int) std::sqrt ((double) data->countVariables ());
00315     if (ports == getSize () - 1) {
00316       if (spara == NULL) {
00317         // find matrix vector entries in touchstone dataset
00318         createIndex ();
00319       }
00320       if (sfreq == NULL) {
00321         logprint (LOG_ERROR, "ERROR: file `%s' contains no `frequency' "
00322                   "vector\n", file);
00323       }
00324     }
00325     else {
00326       logprint (LOG_ERROR, "ERROR: file `%s' specifies a %d-port, `%s' "
00327                 "requires a %d-port\n", file, ports, getName (),
00328                 getSize () - 1);
00329     }
00330   }
00331 }
00332 
00333 void spfile::initSP (void) {
00334   // allocate S-parameter matrix
00335   allocMatrixS ();
00336   // initialize data
00337   prepare ();
00338 }
00339 
00340 /* The function creates an additional data vector for the given matrix
00341    entry and adds it to the dataset. */
00342 void spfile::createVector (int r, int c) {
00343   int i = r * getSize () + c;
00344   spara[i].r = r;
00345   spara[i].c = c;
00346   qucs::vector * v = new qucs::vector (matvec::createMatrixString ("S", r, c),
00347                                sfreq->getSize ());
00348   v->setDependencies (new strlist ());
00349   v->getDependencies()->add (sfreq->getName ());
00350   data->addVariable (v);
00351   spara[i].v = v;
00352 }
00353 
00354 /* This function goes through the dataset stored within the original
00355    touchstone file and looks for the S-parameter matrices and
00356    frequency vector.  It also tries to find the noise parameter
00357    data. */
00358 void spfile::createIndex (void) {
00359   qucs::vector * v; int s = getSize ();
00360   char * n;
00361   const char *name;
00362   int r, c, i;
00363 
00364   // go through list of dependency vectors and find frequency vectors
00365   for (v = data->getDependencies (); v != NULL; v = (::vector *) v->getNext ()) {
00366     if ((name = v->getName ()) != NULL) {
00367       if (!strcmp (name, "frequency")) sfreq = v;
00368       else if (!strcmp (name, "nfreq")) nfreq = v;
00369     }
00370   }
00371 
00372   // create vector index
00373   spara = new spfile_vector[s * s] ();
00374 
00375   // go through list of variable vectors and find matrix entries
00376   for (v = data->getVariables (); v != NULL; v = (::vector *) v->getNext ()) {
00377     if ((n = matvec::isMatrixVector (v->getName (), r, c)) != NULL) {
00378       // save matrix vector indices
00379       i = r * s + c;
00380       spara[i].r = r;
00381       spara[i].c = c;
00382       spara[i].prepare (v, sfreq, false, interpolType, dataType);
00383       paraType = n[0];  // save type of touchstone data
00384       free (n);
00385     }
00386     if ((name = v->getName ()) != NULL) {
00387       // find noise parameter vectors
00388       if (!strcmp (name, "Rn")) {
00389         RN = new spfile_vector ();
00390         RN->prepare (v, nfreq, true, interpolType, dataType);
00391       }
00392       else if (!strcmp (name, "Fmin")) {
00393         FMIN = new spfile_vector ();
00394         FMIN->prepare (v, nfreq, true, interpolType, dataType);
00395       }
00396       else if (!strcmp (name, "Sopt")) {
00397         SOPT = new spfile_vector ();
00398         SOPT->prepare (v, nfreq, false, interpolType, dataType);
00399       }
00400     }
00401   }
00402 }
00403 
00404 /* This function computes the noise correlation matrix of a twoport
00405    based upon the noise parameters and the given S-parameter
00406    matrix. */
00407 matrix spfile::correlationMatrix (nr_double_t Fmin, nr_complex_t Sopt,
00408                                   nr_double_t Rn, matrix s) {
00409   assert (s.getCols () == s.getRows () && s.getCols () == 2);
00410   matrix c (2);
00411   nr_complex_t Kx = 4 * Rn / z0 / norm (1.0 + Sopt);
00412   c.set (0, 0, (Fmin - 1) * (norm (s.get (0, 0)) - 1) +
00413          Kx * norm (1.0 - s.get (0, 0) * Sopt));
00414   c.set (1, 1, norm (s.get (1, 0)) * ((Fmin - 1) + Kx * norm (Sopt)));
00415   c.set (0, 1, s.get (0, 0) / s.get (1, 0) * c.get (1, 1) -
00416          conj (s.get (1, 0)) * conj (Sopt) * Kx);
00417   c.set (1, 0, conj (c.get (0, 1)));
00418   return c;
00419 }
00420 
00421 /* The function computes the noise figure and noise parameters for the
00422    given S-parameter and noise correlation matrices of a twoport. */
00423 nr_double_t spfile::noiseFigure (matrix s, matrix c, nr_double_t& Fmin,
00424                                  nr_complex_t& Sopt, nr_double_t& Rn) {
00425   assert (s.getCols () == s.getRows () && c.getCols () == c.getRows () &&
00426           s.getCols () == 2 && c.getCols () == 2);
00427   nr_complex_t n1, n2;
00428   n1 = c.get (0, 0) * norm (s.get (1, 0)) -
00429     2 * real (c.get (0, 1) * s.get (1, 0) * conj (s.get (0, 0))) +
00430     c.get (1, 1) * norm (s.get (0, 0));
00431   n2 = 2.0 * (c.get (1, 1) * s.get (0, 0) -
00432               c.get (0, 1) * s.get (1, 0)) / (c.get (1, 1) + n1);
00433 
00434   // optimal source reflection coefficient
00435   Sopt = 1 - norm (n2);
00436   if (real (Sopt) < 0.0)
00437     Sopt = (1.0 + std::sqrt (Sopt)) / n2;  // avoid a negative radicant
00438   else
00439     Sopt = (1.0 - std::sqrt (Sopt)) / n2;
00440 
00441   // minimum noise figure
00442   Fmin = real (1.0 + (c.get (1, 1) - n1 * norm (Sopt)) /
00443                norm (s.get (1, 0)) / (1.0 + norm (Sopt)));
00444 
00445   // equivalent noise resistance
00446   Rn = real ((c (0, 0) - 2.0 *
00447               real (c (0, 1) * conj ((1.0 + s (0, 0)) / s (1, 0))) +
00448               c (1, 1) * norm ((1.0 + s (0, 0)) / s (1, 0))) / 4.0);
00449   Rn = Rn * z0;
00450 
00451   // noise figure itself
00452   return real (1.0 + c.get (1, 1) / norm (s.get (1, 0)));
00453 }
00454 
00455 void spfile::initDC (void) {
00456   // get appropriate property value
00457   const char * const dc = getPropertyString ("duringDC");
00458 
00459   // a short during DC including the reference node
00460   if (!strcmp (dc, "shortall")) {
00461     int v, n, lastnode = getSize () - 1;
00462     setVoltageSources (lastnode);
00463     allocMatrixMNA ();
00464     // place zero voltage sources
00465     for (v = VSRC_1, n = NODE_1; n < lastnode; n++, v++) {
00466       voltageSource (v, n, lastnode);
00467     }
00468     return;
00469   }
00470   // a short during DC excluding the reference node
00471   if (!strcmp (dc, "short")) {
00472     int v, n, lastnode = getSize () - 2;
00473     setVoltageSources (lastnode);
00474     allocMatrixMNA ();
00475     // place zero voltage sources
00476     for (v = VSRC_1, n = NODE_1; n < lastnode; n++, v++) {
00477       voltageSource (v, n, lastnode);
00478     }
00479     return;
00480   }
00481   // an open during DC
00482   else if (!strcmp (dc, "open")) {
00483     setVoltageSources (0);
00484     allocMatrixMNA ();
00485     return;
00486   }
00487   // none specified, DC value of IDFT ?
00488   else {
00489     setVoltageSources (0);
00490     allocMatrixMNA ();
00491   }
00492 }
00493 
00494 void spfile::initAC (void) {
00495   setVoltageSources (0);
00496   allocMatrixMNA ();
00497   initSP ();
00498 }
00499 
00500 void spfile::calcAC (nr_double_t frequency) {
00501   // nothing to do if the given file type had errors
00502   if (spara == NULL || sfreq == NULL) return;
00503   // calculate interpolated S-parameters
00504   calcSP (frequency);
00505   // convert S-parameters to Y-parameters
00506   setMatrixY (stoy (getMatrixS ()));
00507 }
00508 
00509 void spfile::calcNoiseAC (nr_double_t frequency) {
00510   // nothing to do if the given file type had errors
00511   if (spara == NULL || nfreq == NULL) return;
00512   setMatrixN (cstocy (calcMatrixCs (frequency), getMatrixY () * z0) / z0);
00513 }
00514 
00515 void spfile::initTR (void) {
00516   initDC ();
00517 }
00518 
00519 // properties
00520 PROP_REQ [] = {
00521   { "File", PROP_STR, { PROP_NO_VAL, "spfile.snp" }, PROP_NO_RANGE },
00522   PROP_NO_PROP };
00523 PROP_OPT [] = {
00524   { "Data", PROP_STR, { PROP_NO_VAL, "polar" },
00525     PROP_RNG_STR2 ("rectangular", "polar") },
00526   { "Interpolator", PROP_STR, { PROP_NO_VAL, "linear" },
00527     PROP_RNG_STR2 ("linear", "cubic") },
00528   { "Temp", PROP_REAL, { 26.85, PROP_NO_STR }, PROP_MIN_VAL (K) },
00529   { "duringDC", PROP_STR, { PROP_NO_VAL, "open" },
00530     PROP_RNG_STR4 ("open", "short", "shortall", "unspecified") },
00531   PROP_NO_PROP };
00532 struct define_t spfile::cirdef =
00533   { "SPfile",
00534     PROP_NODES, PROP_COMPONENT, PROP_NO_SUBSTRATE, PROP_LINEAR, PROP_DEF };