zqrsl.c

DTMK软件开发包,此为开源软件,是一款很好的医学图像开发资源.

    goto L250;
/*<    40 continue >*/
L40:

/*        set up to compute qy or qty. */

/*<          if (cqy) call zcopy(n,y,1,qy,1) >*/
    if (cqy) {
        zcopy_(n, &y[1], &c__1, &qy[1], &c__1);
    }
/*<          if (cqty) call zcopy(n,y,1,qty,1) >*/
    if (cqty) {
        zcopy_(n, &y[1], &c__1, &qty[1], &c__1);
    }
/*<          if (.not.cqy) go to 70 >*/
    if (! cqy) {
        goto L70;
    }

/*           compute qy. */

/*<             do 60 jj = 1, ju >*/
    i__1 = ju;
    for (jj = 1; jj <= i__1; ++jj) {
/*<                j = ju - jj + 1 >*/
        j = ju - jj + 1;
/*<                if (cabs1(qraux(j)) .eq. 0.0d0) go to 50 >*/
        i__2 = j;
        i__3 = j;
        z__1.r = qraux[i__3].r * 0. - qraux[i__3].i * -1., z__1.i = qraux[
                i__3].i * 0. + qraux[i__3].r * -1.;
        if ((d__1 = qraux[i__2].r, abs(d__1)) + (d__2 = z__1.r, abs(d__2)) ==
                0.) {
            goto L50;
        }
/*<                   temp = x(j,j) >*/
        i__2 = j + j * x_dim1;
        temp.r = x[i__2].r, temp.i = x[i__2].i;
/*<                   x(j,j) = qraux(j) >*/
        i__2 = j + j * x_dim1;
        i__3 = j;
        x[i__2].r = qraux[i__3].r, x[i__2].i = qraux[i__3].i;
/*<                   t = -zdotc(n-j+1,x(j,j),1,qy(j),1)/x(j,j) >*/
        i__2 = *n - j + 1;
        zdotc_(&z__3, &i__2, &x[j + j * x_dim1], &c__1, &qy[j], &c__1);
        z__2.r = -z__3.r, z__2.i = -z__3.i;
        z_div(&z__1, &z__2, &x[j + j * x_dim1]);
        t.r = z__1.r, t.i = z__1.i;
/*<                   call zaxpy(n-j+1,t,x(j,j),1,qy(j),1) >*/
        i__2 = *n - j + 1;
        zaxpy_(&i__2, &t, &x[j + j * x_dim1], &c__1, &qy[j], &c__1);
/*<                   x(j,j) = temp >*/
        i__2 = j + j * x_dim1;
        x[i__2].r = temp.r, x[i__2].i = temp.i;
/*<    50          continue >*/
L50:
/*<    60       continue >*/
/* L60: */
        ;
    }
/*<    70    continue >*/
L70:
/*<          if (.not.cqty) go to 100 >*/
    if (! cqty) {
        goto L100;
    }

/*           compute ctrans(q)*y. */

/*<             do 90 j = 1, ju >*/
    i__1 = ju;
    for (j = 1; j <= i__1; ++j) {
/*<                if (cabs1(qraux(j)) .eq. 0.0d0) go to 80 >*/
        i__2 = j;
        i__3 = j;
        z__1.r = qraux[i__3].r * 0. - qraux[i__3].i * -1., z__1.i = qraux[
                i__3].i * 0. + qraux[i__3].r * -1.;
        if ((d__1 = qraux[i__2].r, abs(d__1)) + (d__2 = z__1.r, abs(d__2)) ==
                0.) {
            goto L80;
        }
/*<                   temp = x(j,j) >*/
        i__2 = j + j * x_dim1;
        temp.r = x[i__2].r, temp.i = x[i__2].i;
/*<                   x(j,j) = qraux(j) >*/
        i__2 = j + j * x_dim1;
        i__3 = j;
        x[i__2].r = qraux[i__3].r, x[i__2].i = qraux[i__3].i;
/*<                   t = -zdotc(n-j+1,x(j,j),1,qty(j),1)/x(j,j) >*/
        i__2 = *n - j + 1;
        zdotc_(&z__3, &i__2, &x[j + j * x_dim1], &c__1, &qty[j], &c__1);
        z__2.r = -z__3.r, z__2.i = -z__3.i;
        z_div(&z__1, &z__2, &x[j + j * x_dim1]);
        t.r = z__1.r, t.i = z__1.i;
/*<                   call zaxpy(n-j+1,t,x(j,j),1,qty(j),1) >*/
        i__2 = *n - j + 1;
        zaxpy_(&i__2, &t, &x[j + j * x_dim1], &c__1, &qty[j], &c__1);
/*<                   x(j,j) = temp >*/
        i__2 = j + j * x_dim1;
        x[i__2].r = temp.r, x[i__2].i = temp.i;
/*<    80          continue >*/
L80:
/*<    90       continue >*/
/* L90: */
        ;
    }
/*<   100    continue >*/
L100:

/*        set up to compute b, rsd, or xb. */

/*<          if (cb) call zcopy(k,qty,1,b,1) >*/
    if (cb) {
        zcopy_(k, &qty[1], &c__1, &b[1], &c__1);
    }
/*<          kp1 = k + 1 >*/
    kp1 = *k + 1;
/*<          if (cxb) call zcopy(k,qty,1,xb,1) >*/
    if (cxb) {
        zcopy_(k, &qty[1], &c__1, &xb[1], &c__1);
    }
/*<          if (cr .and. k .lt. n) call zcopy(n-k,qty(kp1),1,rsd(kp1),1) >*/
    if (cr && *k < *n) {
        i__1 = *n - *k;
        zcopy_(&i__1, &qty[kp1], &c__1, &rsd[kp1], &c__1);
    }
/*<          if (.not.cxb .or. kp1 .gt. n) go to 120 >*/
    if (! cxb || kp1 > *n) {
        goto L120;
    }
/*<             do 110 i = kp1, n >*/
    i__1 = *n;
    for (i__ = kp1; i__ <= i__1; ++i__) {
/*<                xb(i) = (0.0d0,0.0d0) >*/
        i__2 = i__;
        xb[i__2].r = 0., xb[i__2].i = 0.;
/*<   110       continue >*/
/* L110: */
    }
/*<   120    continue >*/
L120:
/*<          if (.not.cr) go to 140 >*/
    if (! cr) {
        goto L140;
    }
/*<             do 130 i = 1, k >*/
    i__1 = *k;
    for (i__ = 1; i__ <= i__1; ++i__) {
/*<                rsd(i) = (0.0d0,0.0d0) >*/
        i__2 = i__;
        rsd[i__2].r = 0., rsd[i__2].i = 0.;
/*<   130       continue >*/
/* L130: */
    }
/*<   140    continue >*/
L140:
/*<          if (.not.cb) go to 190 >*/
    if (! cb) {
        goto L190;
    }

/*           compute b. */

/*<             do 170 jj = 1, k >*/
    i__1 = *k;
    for (jj = 1; jj <= i__1; ++jj) {
/*<                j = k - jj + 1 >*/
        j = *k - jj + 1;
/*<                if (cabs1(x(j,j)) .ne. 0.0d0) go to 150 >*/
        i__2 = j + j * x_dim1;
        i__3 = j + j * x_dim1;
        z__1.r = x[i__3].r * 0. - x[i__3].i * -1., z__1.i = x[i__3].i * 0. +
                x[i__3].r * -1.;
        if ((d__1 = x[i__2].r, abs(d__1)) + (d__2 = z__1.r, abs(d__2)) != 0.)
                {
            goto L150;
        }
/*<                   info = j >*/
        *info = j;
/*           ......exit */
/*<                   go to 180 >*/
        goto L180;
/*<   150          continue >*/
L150:
/*<                b(j) = b(j)/x(j,j) >*/
        i__2 = j;
        z_div(&z__1, &b[j], &x[j + j * x_dim1]);
        b[i__2].r = z__1.r, b[i__2].i = z__1.i;
/*<                if (j .eq. 1) go to 160 >*/
        if (j == 1) {
            goto L160;
        }
/*<                   t = -b(j) >*/
        i__2 = j;
        z__1.r = -b[i__2].r, z__1.i = -b[i__2].i;
        t.r = z__1.r, t.i = z__1.i;
/*<                   call zaxpy(j-1,t,x(1,j),1,b,1) >*/
        i__2 = j - 1;
        zaxpy_(&i__2, &t, &x[j * x_dim1 + 1], &c__1, &b[1], &c__1);
/*<   160          continue >*/
L160:
/*<   170       continue >*/
/* L170: */
        ;
    }
/*<   180       continue >*/
L180:
/*<   190    continue >*/
L190:
/*<          if (.not.cr .and. .not.cxb) go to 240 >*/
    if (! cr && ! cxb) {
        goto L240;
    }

/*           compute rsd or xb as required. */

/*<             do 230 jj = 1, ju >*/
    i__1 = ju;
    for (jj = 1; jj <= i__1; ++jj) {
/*<                j = ju - jj + 1 >*/
        j = ju - jj + 1;
/*<                if (cabs1(qraux(j)) .eq. 0.0d0) go to 220 >*/
        i__2 = j;
        i__3 = j;
        z__1.r = qraux[i__3].r * 0. - qraux[i__3].i * -1., z__1.i = qraux[
                i__3].i * 0. + qraux[i__3].r * -1.;
        if ((d__1 = qraux[i__2].r, abs(d__1)) + (d__2 = z__1.r, abs(d__2)) ==
                0.) {
            goto L220;
        }
/*<                   temp = x(j,j) >*/
        i__2 = j + j * x_dim1;
        temp.r = x[i__2].r, temp.i = x[i__2].i;
/*<                   x(j,j) = qraux(j) >*/
        i__2 = j + j * x_dim1;
        i__3 = j;
        x[i__2].r = qraux[i__3].r, x[i__2].i = qraux[i__3].i;
/*<                   if (.not.cr) go to 200 >*/
        if (! cr) {
            goto L200;
        }
/*<                      t = -zdotc(n-j+1,x(j,j),1,rsd(j),1)/x(j,j) >*/
        i__2 = *n - j + 1;
        zdotc_(&z__3, &i__2, &x[j + j * x_dim1], &c__1, &rsd[j], &c__1);
        z__2.r = -z__3.r, z__2.i = -z__3.i;
        z_div(&z__1, &z__2, &x[j + j * x_dim1]);
        t.r = z__1.r, t.i = z__1.i;
/*<                      call zaxpy(n-j+1,t,x(j,j),1,rsd(j),1) >*/
        i__2 = *n - j + 1;
        zaxpy_(&i__2, &t, &x[j + j * x_dim1], &c__1, &rsd[j], &c__1);
/*<   200             continue >*/
L200:
/*<                   if (.not.cxb) go to 210 >*/
        if (! cxb) {
            goto L210;
        }
/*<                      t = -zdotc(n-j+1,x(j,j),1,xb(j),1)/x(j,j) >*/
        i__2 = *n - j + 1;
        zdotc_(&z__3, &i__2, &x[j + j * x_dim1], &c__1, &xb[j], &c__1);
        z__2.r = -z__3.r, z__2.i = -z__3.i;
        z_div(&z__1, &z__2, &x[j + j * x_dim1]);
        t.r = z__1.r, t.i = z__1.i;
/*<                      call zaxpy(n-j+1,t,x(j,j),1,xb(j),1) >*/
        i__2 = *n - j + 1;
        zaxpy_(&i__2, &t, &x[j + j * x_dim1], &c__1, &xb[j], &c__1);
/*<   210             continue >*/
L210:
/*<                   x(j,j) = temp >*/
        i__2 = j + j * x_dim1;
        x[i__2].r = temp.r, x[i__2].i = temp.i;
/*<   220          continue >*/
L220:
/*<   230       continue >*/
/* L230: */
        ;
    }
/*<   240    continue >*/
L240:
/*<   250 continue >*/
L250:
/*<       return >*/
    return 0;
/*<       end >*/
} /* zqrsl_ */

#ifdef __cplusplus
        }
#endif