function projectTorus()

in packages/miew/src/gfx/geometries/ContactSurface.js [454:548]


  function projectTorus(a, b) {
    const aIdx = itemSize * a;
    const bIdx = itemSize * b;
    const xa = posRad[aIdx];
    const ya = posRad[aIdx + 1];
    const za = posRad[aIdx + 2];
    const r1 = posRad[aIdx + 3];
    let dx = mid.x = posRad[bIdx] - xa;
    let dy = mid.y = posRad[bIdx + 1] - ya;
    let dz = mid.z = posRad[bIdx + 2] - za;
    const innR2 = posRad[bIdx + 3];
    let d2 = dx * dx + dy * dy + dz * dz;

    // This check now redundant as already done in AVHash.withinRadii
    // if( d2 > (( r1 + r2 ) * ( r1 + r2 )) ){ return; }

    const d = Math.sqrt(d2);

    // Find angle between a->b vector and the circle
    // of their intersection by cosine rule
    const cosA = (r1 * r1 + d * d - innR2 * innR2) / (2.0 * r1 * d);

    // distance along a->b at intersection
    const dmp = r1 * cosA;

    mid.normalize();

    // Create normal to line
    normalToLine(n1, mid);
    n1.normalize();

    // Cross together for second normal vector
    n2.crossVectors(mid, n1);
    n2.normalize();

    // r is radius of circle of intersection
    const rInt = Math.sqrt(r1 * r1 - dmp * dmp);

    n1.multiplyScalar(rInt);
    n2.multiplyScalar(rInt);
    mid.multiplyScalar(dmp);

    mid.x += xa;
    mid.y += ya;
    mid.z += za;

    lastClip = -1;

    const ng = ngTorus;

    for (let innI = 0; innI < probePositions; innI++) {
      const cost = cosTable[innI];
      const sint = sinTable[innI];

      const px = mid.x + cost * n1.x + sint * n2.x;
      const py = mid.y + cost * n1.y + sint * n2.y;
      const pz = mid.z + cost * n1.z + sint * n2.z;

      if (obscured(px, py, pz, a, b) === -1) {
        // As above, iterate over our grid...
        // px, py, pz in grid coords
        const iax = Math.floor(scaleFactor * (px - min[0]));
        const iay = Math.floor(scaleFactor * (py - min[1]));
        const iaz = Math.floor(scaleFactor * (pz - min[2]));

        const minx = Math.max(0, iax - ng);
        const miny = Math.max(0, iay - ng);
        const minz = Math.max(0, iaz - ng);

        const maxx = Math.min(dim[0], iax + ng + 2);
        const maxy = Math.min(dim[1], iay + ng + 2);
        const maxz = Math.min(dim[2], iaz + ng + 2);

        for (let iz = minz; iz < maxz; iz++) {
          dz = pz - gridz[iz];
          const zOffset = dim[1] * dim[0] * iz;
          for (let iy = miny; iy < maxy; iy++) {
            dy = py - gridy[iy];
            const dzy2 = dz * dz + dy * dy;
            const zyOffset = zOffset + dim[0] * iy;
            for (let ix = minx; ix < maxx; ix++) {
              dx = px - gridx[ix];
              d2 = dzy2 + dx * dx;
              const idx = ix + zyOffset;
              const current = grid[idx];

              if (current > 0.0 && d2 < (current * current)) {
                grid[idx] = Math.sqrt(d2);
              }
            }
          }
        }
      }
    }
  }