//
//     InterDigitalFilter
//
//
//  from a C++ program by Dale Heatherington, WA4DSY (1996);
//  the original BASIC  program was published in the January 1985 issue of
//  Ham Radio magazine on page 12.  The authors are Jerry Hinshaw, N6JH
//  and Shahrokh Monemzadeh.
//  Überarbeitet von OE7WPA, Oktober 2020

// Konstanten
const INCH = 25.4;   // inch --> mm
const M10 = 1.0/Math.log(10.0);
// global variables
var jkplot;

var G=[200], C=[200], RK=[200], AK=[200], A=[200], B=[200];

function log10( x) {
  //
  //   Dekadischer Logarithmus
  //
  if(x > 0.0)   return M10*Math.log(x);
  return -1001.0;
}

function IDF() {
  // Deklariere Methoden
  this.compute = compute;
  this.showPlot = showPlot;

  function FuncFNRJ(ta, b, c, d) {
    return ((b*c)-(ta*d)) / ((c*c) + (d*d));
  }

  function chebyResponse(s, norder, ripple) {
    var eps2 = Math.pow(10.0, ripple/10.0) - 1.0;
    var t2 = tcheby(norder, s); t2=t2*t2;
    return (1.0+eps2*t2);
  }
  
  function tcheby(norder, x) {
    switch (norder)  {
      case 0:  return 1.0;
      case 1:  return x;
//
//      Rekursive Vorgangsweise
//
      default:  return 2.0*x*tcheby(norder-1, x) - tcheby(norder-2, x);
    }
  }

  function cheb(ele, rip) {
    //
    //    Berechnet die normalisierten LC Filterwerte
    //
    var beta, gamma, c;
    var k;
    c = (2.0 * rip) / 17.37;
    beta = Math.log((Math.exp(c)+1.0) / (Math.exp(c) -1.0));
    gamma = 0.5 * (Math.exp(beta/(2.0*ele)) - Math.exp(-beta/(2.0*ele)));
    for(k=0; k<ele; k++)    {
      A[k] = Math.sin(0.5 * (2.0 * (k+1) - 1.0) * Math.PI / ele);
      B[k] = (gamma * gamma) + Math.pow(Math.sin((k+1)*Math.PI/ele) , 2.0);
    }
    G[0] = 2.0 * A[0] / gamma;
    for(k=1; k<ele; k++)   {
      G[k] = 4.0 * A[k-1] * A[k] / (B[k-1] * G[k-1]);
    }
    var nn = ele / 2;
    var nnn = (ele+1) / 2;
    var s = nnn - nn;
    if(s <= 0)   { 
      G[ele] = Math.pow(((Math.exp(beta/2)+1) / (Math.exp(beta/2)-1)) , 2);
    }
    G[ele] = 1.0;
  }

  function butter(ele) {
    //
    //    Berechnet die normalisierten LC Filterwerte
    //
    var k;
    for(k=0;k<ele;k++)   {
      G[k] = 2.0 * Math.sin(1.57079633 * (2.0*(k+1) - 1) / ele);
    }
    G[ele] =  1.0;
  }

//-----------------------------------------------------------
  function compute() {
    var k, ip, conter=0;
    var bw3gc, bwrgc, bw3, b, ca;
    //
    //    Berechne grundsätzliche Eigenschaften die Kopplung zwischen den Elementen
    //
    var F1 = frequency - 0.0005 * bandwidth;
    var F2 = frequency + 0.0005 * bandwidth;
    if (ripple == 0.0)    {
      bw3gc = F2 - F1;
      bwrgc = 0.0;
      bw3   = 1.0;
      butter(elements);
    } else {
      b = 1.0 / Math.sqrt(Math.pow(10.0, 0.1*ripple)-1);
      ca = Math.log(b + Math.sqrt(b * b - 1)) / elements;
      bw3 = (Math.exp(ca) + Math.exp(-ca)) / 2.0;
      cheb(elements, ripple);
      bwrgc = F2 - F1;
      bw3gc = bwrgc * bw3;
    }

    var w = 2.0 * (F2 - F1) / (F2 + F1);
    var qf = frequency / bw3gc;
    var nfm = elements-1;
    var qwvl = 11.8028 / (4.0 * frequency);
    for(k = 0; k<nfm; k++)      {
      AK[k] = 1.0 / (bw3 * Math.sqrt(G[k] * G[k+1]));
      RK[k] = AK[k] / qf;
    }
    var AK0 = G[0]*bw3;
    AK[elements-1] = G[0] * bw3;  
    AK[elements] = 0.0;
    var qs = G[0] * bw3 * qf;
    var canh = (Math.exp(2.0*Math.PI*endspacing/height)-1.0) / (Math.exp(2.0*Math.PI*endspacing/height)+1.0);
    var zm = 59.9585 * Math.log(4.0*height / (Math.PI*diameter));
    var ze = 59.9585 * Math.log(canh * height * 4.0 / (Math.PI * diameter));
    var rkm = RK[0] * Math.sqrt(zm/ze);
    var z = Math.PI * diameter / (2.0 * height);
    var coth = (Math.exp(z) + 1.0) / (Math.exp(z) - 1.0);
    var y = Math.PI * rkm / 4.0;
    var t = Math.pow(coth, y);
    C[0] = height/Math.PI * Math.log((t+1.0) / (t-1.0));
    for(k=2; k<=elements-2; k++)      {
      y = Math.PI * RK[k-1] / 4.0;
      t = Math.pow(coth, y);
      C[k-1] = height/Math.PI * Math.log((t+1.0) / (t-1.0));
    }
    C[elements-2] = C[0];
    var x = Math.sqrt(Math.PI * rload / (4.0 * ze * qs));
    var aq = 2.0 * qwvl * Math.atan(x / Math.sqrt(1.0 - x * x)) / Math.PI;
    var qunloaded = 2200.0 * height * Math.sqrt(frequency);
    var sumg = 0.0;
    for(k=0; k<elements; k++) {
      sumg += G[k];
    }
    var bloss = 4.34 * frequency * sumg / (qunloaded * (F2 - F1));
    var delay = sumg / (2.0 * Math.PI * (F2 - F1));


    var output = document.getElementById("Ergebnis");
    txt = "<b>"+elements+" Element Interdigital Bandpass Filter</b><br><br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Bandpass Welligkeit = "+fformat(ripple, 1)+" dB<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Center Frequenz = "+fformat(frequency*1000.0, 3)+" MHz<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; 3 dB Bandbreite = "+fformat(bw3gc*1000.0, 3)+" MHz<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Fractional Bandbreite = "+fformat(w, 4)+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Ripple Bandbreite = "+fformat(bwrgc * 1000.0, 2)+" MHz<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Cutoff Frequencies = "+fformat(F1*1000.0, 5)+" and "+fformat(F2*1000.0, 5)+" MHz<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Filter Q = "+fformat(qf*1.0, 5)+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Gesch&auml;tztes unbelastetes Qu = "+fformat(qunloaded, 3)+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Insertion Loss basierend auf Qu = "+fformat(bloss, 3)+" dB<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Verz&ouml;gerung der Center-Frequenz = "+fformat(delay, 3)+" ns<br>";
    //
    //        Berechnung der Länge der inneren Elemente
    //
    var wo = 2.0 * Math.PI * frequency * 1.0e+9;
    var f = diameter / height;
    var cf = (-0.0000422 + 0.0857397 * f + 0.0067853 * f*f 
               -9.092165e-2 * f*f*f + 0.169088 * f*f*f*f) * Math.PI * height * INCH/10.0;
    var ww   = wo * 1.0e-12;
    var b2   = Math.PI * aq / ( 2.0 * qwvl);
    var gg   = 1.0 / rload;
    var bb   = - 1.0/ (ze * Math.tan(b2));
    var el1  = 0.8 * qwvl;

    var ang  = el1 * Math.PI / (2.0 * qwvl);
    var b1   = ang - b2;
    var yl   = -1.0 / (zm * Math.tan(ang));
    var cp   = ww * (cf + 0.17655 * diameter * diameter / (qwvl - el1));
    var y1   = cp + yl;
    var el2  = 0.87 * qwvl;

    ang  = el2 * Math.PI / ( 2.0* qwvl);
    var b4   = ang - b2;
    yl   = -1.0 / (zm * Math.tan(ang));
    var cd   = ww * (cf + 0.17655 * diameter * diameter / (qwvl - el2));
    var y2   = cd + yl;
    var el3  = 0.95 * qwvl;

    ang  = el3 * Math.PI / (2.0 * qwvl);
    var  b5  = ang - b2;
    yl  = -1.0 / ( zm * Math.tan(ang));
    var  cq  = ww * (cf + 0.17655 * diameter * diameter / (qwvl - el3));
    var  y3  = cq + yl;
    var elem = y3 * y2 * el1 / ((y1 - y2) * (y1-y3)) + y1 * y3 * el2
                    / ((y2-y1)*(y2-y3)) + y1 * y2 * el3 / ((y3-y1) * (y3-y2));

  //
  //        Berechnung der Länge der End-Elemente
  //
    var tann = Math.tan(b1);
    yl  = FuncFNRJ(gg, bb+tann/ze, 1-ze*bb*tann, ze*gg*tann);
    y1  = cp + yl;
    tann = Math.tan(b4);
    yl = FuncFNRJ(gg, bb+tann/ze, 1-ze*bb*tann, ze*gg*tann);
    y2  = cd + yl;
    tann = Math.tan(b5);
    yl = FuncFNRJ(gg, bb+tann/ze, 1-ze*bb*tann, ze*gg*tann);
    y3  = cq + yl;
    var eleq = y3 * y2 * el1 / ((y1 - y2) * (y1 - y3)) + y1 * y3 * el2
             / ((y2 - y1) * (y2 - y3)) + y1 * y2 * el3 / ((y3 - y1) * (y3 - y2));

    txt +=  "<br>Viertel-Wellenl&auml;nge: "+fformat(qwvl*INCH, 5)+" mm<br>";
    txt +=  "L&auml;nge der inneren Elemente: "+fformat(elem*INCH, 3)+" mm<br>";
    txt +=  "L&auml;nge der End-Elemente: "+fformat(eleq*INCH, 3)+" mm<br>";
    txt +=  "H&ouml;he der Elemente über der Ground Plane: "+fformat(height*INCH/2.0, 3)+" mm<br>";
    txt +=  "Element-Durchmesser: "+fformat(diameter*INCH, 3)+" mm<br>";
    txt +=  "Endplatten befinden sich "+fformat(endspacing*INCH, 3)+" mm vom Zentrum der End-Elemente<br>";
    txt +=  "Gamma-Match: Punkt befindet sich "+fformat(aq*INCH, 3)+" mm vom kurzgeschlossenen Ende<br>";
    txt +=  "Line Impedanzen: End-Element: "+fformat(ze, 1)+", andere Elemente: "+fformat(zm, 1)+"<br>";
    txt +=  "Load Impedanz: "+fformat(rload, 1)+" Ohm<br>";

    var Tab = 10;
    txt +=  "<br>";
    txt +=  "Dimensionen: (mm)<br>";
    txt +=  "<table>";
    txt +=  "<tr><td>Element No.</td><td>&nbsp;&nbsp;Abstand von </td><td>&nbsp;&nbsp;Abstand zu </td><td>&nbsp;G[k]</td><td>&nbsp;Q/Coup</td></tr>";
    txt +=  "<tr><td>           </td><td>Endplatte</td>     <td>n&auml;chstes Element</td><td> </td><td></td></tr>";
    txt +=  "<tr><td>Ende</td><td>0.0</td><td>"+fformat(endspacing*INCH*1.0,4)+"</td><td></td><td></td></tr>";
    var dom = endspacing;
    C[nfm] = endspacing;
    for(var i=0;i<=nfm;i++)  {
      txt +=  "<tr><td>"+(i+1)+"</td><td>"+fformat(dom*INCH*1.0,4)+"</td><td>"+fformat(C[i]*INCH*1.0,4)+"</td><td>"
                      +fformat(G[i],4)+"</td><td>&nbsp;&nbsp;"+fformat(AK[i],4)+"</td></tr>";
      dom += C[i];
    }
    txt +=  "<tr><td>Ende</td><td>"+fformat(dom*INCH*1.0,4)+"</td><td></td><td></td><td></tr>";
    txt +=  "</table><br>";
    txt +=  "Innen-Abmessungen der Kavit&auml;t [mm]: "+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; Breite:    "+fformat(qwvl*INCH*1.0, 4)+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; L&auml;nge:     "+fformat(dom*INCH*1.0, 4)+"<br>";
    txt +=  "&nbsp;&nbsp;&nbsp; H&ouml;he: "+fformat(height*INCH*1.0, 4)+"<br>";

    output.innerHTML = txt;
  }
  function showPlot()  {
    jkplot.startCurve(0,-100,1,"red");
    var h1 = 0.0, h2 = 0.0, h3 = 0.0, f = 0.0;
    var nf = 200;
    var f1 = frequency - frange;
    var f2 = frequency + frange;
    var df = (f2-f1)/(nf);
    var deltaf;
    var F1 = frequency - 0.0005 * bandwidth;
    var F2 = frequency + 0.0005 * bandwidth;
    var w = 2.0 * (F2 - F1) / (F2 + F1);
    for(var i=0; i<nf; i++) {
      f = f1 + i*df;
      deltaf = Math.abs((2.0 * (f-frequency))/(w*frequency));
      if(ripple == 0.0)   {
        h1 = 10.0 * log10(1.0 + Math.pow(deltaf,2.0*elements));
      } else {
        h1 = 10.0 * log10(chebyResponse(deltaf, elements, ripple));
      }
      jkplot.addCurve(f, -h1);
    }      
  }
}


var elements = 3;              //Anzahl der Elemente
var ripple = 0;                //Ripple des Cheby-Filters
var frequency=2405.0, bandwidth=30.0, rload=50.0;   //Center-Frequenz in GHz, BW in MHz, Load Resistance in Ohm
var height=30.0, diameter=8.0, endspacing=15.0;  // Dicke, Element-Durchmesser, End-Abstand
var frange=250.0;                // Frequency Span inm Plot

function readFields() {
  elements   = Number(document.getElementById("n").value);
  ripple     = Number(document.getElementById("ripple").value);
  frequency  = Number(document.getElementById("freq").value) /1000.0;
  bandwidth  = Number(document.getElementById("bw").value);
  var F1 = frequency - 0.0005 * bandwidth;
  var F2 = frequency + 0.0005 * bandwidth;
  var w = 2.0 * (F2 - F1) / (F2 + F1);
  if(w < 0.01 || w > 0.1)  {
    document.getElementById("bw").style.background="yellow";   //Warnung bei unrealistischen Werten der fraktionalen Bandbreite!
  } else {
    document.getElementById("bw").style.background="white";
  }
  rload      = Number(document.getElementById("z").value);
  height     = Number(document.getElementById("thick").value) /INCH;
  diameter   = Number(document.getElementById("rod").value) /INCH;
  endspacing = Number(document.getElementById("rodend").value) /INCH;
  frange     = Number(document.getElementById("fspan").value) /1000.0;
}
function handleComputeButton() {
  readFields();
  jkplot.setXrange(Math.max(0.1, frequency-frange), frequency+frange);
  jkplot.setYrange(-80.0, 5.0);
  jkplot.clear(); 
  filter.compute();
  filter.showPlot(); 
}
function checkKey(e) {
  var k;
  if(window.event) {
    k = window.event.keyCode;
  } else {
    k = e.keyCode;
  } 
  if(k==13) {
    readFields();
    jkplot.setXrange(Math.max(0.1, frequency-frange), frequency+frange);   
    jkplot.clear();
    filter.showPlot(); 
    filter.compute(); 
  }
}
    

function setupButtons() {
  var b = document.getElementById("Berechnen");
  if(b.addEventListener) {
    b.addEventListener("click", handleComputeButton, false); 
  } else if(b.attachEvent) {
    b.addEventListener("onclick", handleComputeButton); 
  }
  document.onkeydown = checkKey;
}

function fformat(x, k) {
  var xx = Math.max(1.0e-40,Math.abs(x));
  var y = Math.pow(10.0, k-Math.floor(log10(xx)));
  return Math.floor(x*y+0.5)/y;
}  


window.onload = function() {
  setupButtons();
  jkplot = new JKPlot();
  jkplot.setup(1000,false,1.405,3.405,"X: ","Frequenz [GHz]",
               600,false,-140,5,"Y: ","S1,2 [dB]");
  jkplot.clear();
  filter = new IDF();
  readFields();
  filter.showPlot(); 
  filter.compute(); 
}