---

pct.c

Go to the documentation of this file.

#include        "sadie.h"

/*-Copyright Information------------------------------------------------------*/
/* Copyright (c) 1988 by the University of Arizona Digital Image Analysis Lab */
/*----------------------------------------------------------------------------*/
/*-General Information--------------------------------------------------------*/
/*                                                                            */
/*   This function performs a principal component transform on a multi-band   */
/*   image.                                                                   */
/*                                                                            */
/*----------------------------------------------------------------------------*/
/*-Interface Information------------------------------------------------------*/
void PCT (
IMAGE *in,      /*  I   Pointer to the input image                            */
short inc,      /*  I   Line and pixel subsample increment to be used for     */
                /*      calculating the covariance matrix.                    */
IMAGE **out     /*  O   Address of a pointer to the output image.             */
                /*      The input image may be overwritten (out = &in).       */
/*----------------------------------------------------------------------------*/
) { register short i, j, k, l;
    char   msg[SLEN];
    double *cov=0, *cor=0, *eig=0, sum, pinc, psum;

    if (TIMES) TIMING(T_PCT);
    if (NAMES) {
        MESSAGE('I',"");
        MESSAGE('I',"PCT");
        MESSAGE('I',"");
        sprintf(msg," Input image:                         %s",in->text);
        MESSAGE('I',msg);
        sprintf(msg," Line and pixel subsample increment:  %d",inc);
        MESSAGE('I',msg);
        MESSAGE('I'," ...............");
    }

    /* check input */
    if (!CHECKIMG(in)) {
        MESSAGE('E'," Can't identify image.");
        goto the_end;
    } else if (in->nbnd == 1) {
        MESSAGE('E'," Image must be multi-band.");
        goto the_end;
    } else if (inc <= 0) {
        MESSAGE('E'," Subsample increment must be greater than zero.");
        goto the_end;
    }

    /* create image of appropriate size */
    if (!CHECKIMG(*out)) GETMEM(in->nbnd,in->nlin,in->npix,out);
    if (!*out) goto the_end;

    /* allocate space for cov, cor, and eig */
    if (!(cov=(double *)malloc(in->nbnd*in->nbnd*sizeof(double)))) {
        MESSAGE('E'," Memory request failed.");
        goto the_end;
    }
    if (!(cor=(double *)malloc(in->nbnd*in->nbnd*sizeof(double)))) {
        MESSAGE('E'," Memory request failed.");
        goto the_end;
    }
    if (!(eig=(double *)malloc(in->nbnd*in->nbnd*sizeof(double)))) {
        MESSAGE('E'," Memory request failed.");
        goto the_end;
    }

    /* initialize progress indicator */
    if (LINES  &&  !PROGRESS(psum=0.0)) goto the_end;
    pinc = 1.0/(double)in->nbnd/(double)in->nlin;

    /* compute covariance matrix of the input image and its eigenvectors */
    COVAR(in,inc,cov,cor);
    if (ABORT) goto the_end;
    EIGEN(cov,in->nbnd,eig);
    if (ABORT) goto the_end;

    /* perform principal component transform */
    for (i=0; i<in->nbnd; i++) {
        for (j=0; j<in->nlin; j++) {
            for (k=0; k<in->npix; k++) {
                for (sum=0.0,l=0; l<in->nbnd; l++) {
                    sum += (double)in->data[l][j][k]*eig[l*in->nbnd+i];
                }
                (*out)->data[i][j][k] = (PIXEL)sum;
            }
            if (LINES  &&  !PROGRESS(psum+=pinc)) goto the_end;
        }
    }

    the_end:
    if (cov)   free(cov);
    if (cor)   free(cor);
    if (eig)   free(eig);
    if (TIMES) TIMING(T_EXIT);
}

---