Next: Listing for sinter3d.c Up: Computer programs for Previous: Listing for rand3d.f

Listing for stat3d.c

/************************************************************************/
/*                                                                      */
/*      Program: stat3d.c                                               */
/*      Purpose: To read in a 3-D image and output phase volumes        */
/*               and report the volume and pore-exposed surface area    */
/*               fractions                                              */
/*               Assumes that input image contains only the values 0-6  */
/*      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 ISIZE 100


main(){
        static int mic [ISIZE] [ISIZE] [ISIZE];
        int valin,ix,iy,iz;
        int ix1,iy1,iz1,k;
        long int surftot,volume[7],surface [7];
        FILE *infile;
        char filen[80];

        printf("Enter name of file to open \n");
        scanf("%s",filen);
        printf("%s \n",filen);

        for(ix=0;ix<=6;ix++){
                volume[ix]=surface[ix]=0;
        }

        infile=fopen(filen,"r");

/* Read in image and accumulate volume totals */
        for(iz=0;iz<ISIZE;iz++){
        for(iy=0;iy<ISIZE;iy++){
        for(ix=0;ix<ISIZE;ix++){
                fscanf(infile,"%d",&valin);
                mic [ix] [iy] [iz]=valin;
                volume[valin]+=1;
        }
        }
        }
	
        fclose(infile);

        for(iz=0;iz<ISIZE;iz++){
        for(iy=0;iy<ISIZE;iy++){
        for(ix=0;ix<ISIZE;ix++){
                if(mic [ix] [iy] [iz]!=0){
                        valin=mic [ix] [iy] [iz];
                /* Check six neighboring pixels for porosity */
                        for(k=1;k<=6;k++){

                                 switch (k){
                                        case 1:
                                                ix1=ix-1;
                                                if(ix1<0){ix1+=ISIZE;}
                                                iy1=iy;
                                                iz1=iz;
                                                break;
                                        case 2:
                                                ix1=ix+1;
                                                if(ix1>=ISIZE){ix1-=ISIZE;}
                                                iy1=iy;
                                                iz1=iz;
                                                break;
                                        case 3:
                                                iy1=iy-1;
                                                if(iy1<0){iy1+=ISIZE;}
                                                ix1=ix;
                                                iz1=iz;
                                                break;
                                        case 4:
                                                iy1=iy+1;
                                                if(iy1>=ISIZE){iy1-=ISIZE;}
                                                ix1=ix;
                                                iz1=iz;
                                                break;
                                        case 5:
                                                iz1=iz-1;
                                                if(iz1<0){iz1+=ISIZE;}
                                                iy1=iy;
                                                ix1=ix;
                                                break;
                                        case 6:
                                                iz1=iz+1;
                                                if(iz1>=ISIZE){iz1-=ISIZE;}
                                                iy1=iy;
                                                ix1=ix;
                                                break;
                                        default:
                                                break;
                                }
        if((ix1<0)||(iy1<0)||(iz1<0)||(ix1>=ISIZE)||(iy1>=ISIZE)||(iz1>=ISIZE)){
                                        printf("%d %d %d \n",ix1,iy1,iz1);
                                        exit(1);
                                }
                                if(mic[ix1] [iy1] [iz1]==0){
                                        surface[valin]+=1;
                                }
                        }
                }
        }
        }
        }

        printf("Phase    Volume      Surface     Surface Fraction \n");
/* Only include non-gypsum phases in surface area fraction calculation */
        surftot=surface[1]+surface[2]+surface[3]+surface[4];
printf("Total surface area count (not including gypsum and aggregate) is %ld\n",
              surftot);
        for(k=0;k<=6;k++){
                if(k<=4){
        printf("%d      %8ld     %8ld       %.5f\n",k,volume[k],surface[k],
                    (float)surface[k]/(float)surftot);
                }
                else{
                       printf("%d      %8ld     %8ld\n",k,volume[k],surface[k]);
                }
        }
}