Next: Listing for corrxy2r.c Up: Computer programs for Previous: Listing for statsimp.c

Listing for corrcalc.c

/************************************************************************/
/*                                                                      */
/*      Program: corrcalc.c                                             */
/*      Purpose: To compute the 2-D correlation function in S(x,y) form */
/*              for one or more phases of a cement particle image       */
/*      Programmer: Dale P. Bentz                                       */
/*                  NIST                                                */
/*                  Building 226 Room B-350                             */
/*                  Gaithersburg, MD  20899-0001                        */
/*                  Phone: (301) 975-5865                               */
/*                  E-mail: dale.bentz@nist.gov	                        */
/*                                                                      */
/************************************************************************/
#include <stdio.h>
#include <math.h>

#define CSIZE 60   /* Extent of correlation computation in pixels */

int flag[5];

/* Note that maximum image size is 512*512 */
main()
{
        int mask,i,j,image [512] [512],xsize,ysize;
        int scalef,xoff,yoff,inval,iscale,jscale;
        long int sum,ntot;
        float fsum;
        char filen[80],fileo[80];
        FILE *infile,*outfile;

        do{
                printf("Enter size of image in x and y directions \n");
                scanf("%d %d",&xsize,&ysize);
                printf("%d %d \n",xsize,ysize);
        } while ((xsize>512)||(ysize>512));

/* Assumed file format is a simple list of integers (one per line) */
/* representing phases as indicated below and with y in the inner loop */
        printf("Enter name of image file \n");
        scanf("%s",filen);
        infile=fopen(filen,"r");

        printf("Enter name of correlation file to create\n");
        scanf("%s",fileo);
        outfile=fopen(fileo,"w");
/* Output header consisting of the extent of the correlation matrix */
        fprintf(outfile,"%d %d\n",CSIZE+1,CSIZE+1);

/* If original image is acquired at a resolution much less than */
/* 1 um/pixel, scaling factor can be used to for example skip every other */
/* pixel to adjust to appropriate resolution */
        scalef=1;
        printf("Enter scaling factor to use (default=%d) \n",scalef);
        scanf("%d",&scalef);
        printf("%d \n",scalef);

        printf("Phase assignments are assumed to be as follows: \n\n");
        printf("Phase      Image ID   Mask value \n");
        printf(" \n");
        printf("Pores, etc.    0          1\n");
        printf("C3S            1          2\n");
        printf("C2S            2          4\n");
        printf("C3A            3          8\n");
        printf("C4AF           4         16\n");

/* Composite mask is a simple logical and of the above values */
/* Example: For C3S+C2S, mask would be 2&4 = 6 */
        printf("Enter composite value for mask to employ during this run \n");
	printf("Example- for C3S+C2S, composite value would be 2+4=6 \n");
        scanf("%d",&mask);
        printf("%d\n",mask);

/* Read in the original image and convert values to 2**i format */
/* so that simple and operations may be employed to compute correlation */
        for(i=0;i<xsize;i++){
        for(j=0;j<ysize;j++){
                fscanf(infile,"%d\n",&inval);
                image [i] [j]=(int)(pow(2.,(float)inval));
        }
        }

        fclose(infile);
        printf("Finished reading in file \n");
        fflush(stdout);

/* Now perform calculation */
/* Correlation is computed only over a distance of CSIZE pixels */
        for(i=0;i<=CSIZE;i++){
        for(j=0;j<=CSIZE;j++){
                sum=0;
                ntot=0;

/* iscale and jscale represent limits of correlation calculation since */
/* original image is assumed to be non-periodic */
                iscale=i*scalef;
                jscale=j*scalef;
/* Be sure to skip to every scalef pixel to adjust resolution */
                for(xoff=0;xoff<(xsize-iscale);xoff+=scalef){
                for(yoff=0;yoff<(ysize-jscale);yoff+=scalef){
                        ntot+=1;
/* Increment counter if both pixels contain one of the necessary phases */
                if(((mask&(image[xoff][yoff]))!=0)&&
                       ((mask&(image[xoff+iscale][yoff+jscale]))!=0)){
                       sum+=1;
                }
		
                }
                }
/* Determine and output value for S(x,y) */
        fsum=(float)sum/(float)ntot;
        fprintf(outfile,"%d %d %f \n",i,j,fsum);
        fflush(outfile);

        }
        }

        fclose(outfile);
}