/* Calculate N */
/* Calculate G */
/* Calculate H */
/* Quadratic Arguments */
#include <stdio.h>
#include <math.h>
float quadratic(float a, float b, float c);
float calculations(int i, int j);
int getval(float highval, float lowval, float val);
int halftoning(float lowval, float highval);
int bitset(int i, int j, int gval);
float imagemap(void);
/* Global Variables */
float xeye=0, yeye=0, zeye=3;
float xmax=250, ymax=250;
float xwindow=0, ywindow=0, zwindow=2;
float xlight=1, ylight=0, zlight=0;
float x=0, y=0, z=0, r=1, t, k;
float dirx=0, diry=1, dirz=0;
float pointx, pointy, pointz;
float f, f1, f2, f3;
float greys[250][250];
int bits[500][500];
int inbits[10][10];
float h1, h2, h3, g1, g2, g3, m, n, am1, am2, am3, am;
float bm1, bm2, bm3, bm, an1, an2, an3, an, md1, md2, md3;
int flag;
int main(int argc, char** argv) /* main part of program */
{
int i, j, val, biti, bitj;
char P, P_;
FILE *input, *output;
char file1[20], file2[20];
float norm1, norm2, norm3;
float lightn1, lightn2, lightn3;
float dotnorm1, dotnorm2, dotnorm3, dotprod;
float intlev, lowval=0, highval=0;
double normval, normval2;
/* Print an error message if not incorrect command line is given */
if((argc!=3)||(sscanf(argv[1], "%s", &file1)!=1)
||(sscanf(argv[2], "%s", &file2)!=1))
{
fprintf(stderr, "Incorrect command line given in %s\n", argv[0]);
exit(1);
}
/* Check that the input file is of the correct type and that the Eye */
/* Position is greater than the Window Position */
input=fopen(file1, "r"); /*Open file "file1" to be read */
output=fopen(file2, "w");
/* Open the output file "file2" to be written into */
fprintf(output, "P1\n500 500\n");
/* Read the input file */
fscanf(input, "%c", &P);
fscanf(input, "%c", &P_);
/* If the incorrect input file is being used then print an error message */
if((P!='P')||(P_!='1')){
printf("*** Error: Incorrect Input File ***\n");
exit(1);
}
/* Scan the Bits In from the Input File */
fscanf(input, "%d %d", &i, &j);
printf("%d %d\n", i, j);
for(i=0; i<10; i++);
{
for(j=0; j<10; j++);
{
fscanf(input, "%d", &inbits[i][j]);
}
}
biti=0;
bitj=0;
flag=1;
highval=50;
/* Main loop here */
for(i=0; i<=ymax; i++)
{
for(j=0; j<=xmax; j++)
{
t=calculations(i, j); /* Set t equal to the value Calculated*/
/* from using i and j for the values in calculations part of program */
if (t==0) val=0; /* If t equals zero, set val equal to 0 */
else
{
/* Calculating point on Sphere */
pointz=xeye+(t*(xwindow-xeye));
pointy=yeye+(t*(ywindow-yeye));
pointz=zeye+(t*(zwindow-zeye));
/* Calculating cylinder normal vector */
f1=(pointx-x-(k*dirx))/r;
f2=(pointy-y-(k*diry))/r;
f3=(pointz-z-(k*dirz))/r;
/* Calculating the Normals */
norm1=(pointx-x)/r;
norm2=(pointy-y)/r;
norm3=(pointz-z)/r;
/* Calculate the Unit Vector */
normval=sqrt((norm1*norm1)+(norm2*norm2)+(norm3*norm3));
/* Calculate the Vector of Light Source */
normval2=sqrt((xlight*xlight)+(ylight*ylight)+(zlight*zlight));
lightn1=xlight/normval2;
lightn2=ylight/normval2;
lightn3=zlight/normval2;
/* Check unit Vector */
normval=sqrt((lightn1*lightn1)+(lightn2*lightn2)+(lightn3*lightn3));
/* Calculate Dot Product */
dotnorm1=f1*lightn1;
dotnorm2=f2*lightn2;
dotnorm3=f3*lightn3;
dotprod=dotnorm1+dotnorm2+dotnorm3;
/* Calculate light Levels */
if (dotprod<0) dotprod=0;
intlev=dotprod+0.7;
/* Initialise high and low Levels */
if (flag==1){
highval=0;
lowval=t;
flag=0;
}
val=t*intlev;
if (val>highval) highval=val;
if (val
greys[i][j]=val;
}
}
halftoning(lowval, highval);
for(i=0; i<(xmax*2); i++)
{
for(j=0; j<(ymax*2); j++)
{
fprintf(output, "%d", bits[i][j]);
}
fprintf(output, "\n");
}
}
int halftoning(float lowval, float highval)
{
int i, j;
float greylev[5];
float valrange=highval-lowval;
for(i=0; i<5; i++);
{
greylev[i]=((i+1)*(valrange/5))+lowval;
}
for(i=0; i
for(j=0; j
if(greys[i][j]>greylev[3]) bitset(i,j,0);
else
if(greys[i][j]>greylev[2]) bitset(i,j,1);
else
if(greys[i][j]>greylev[1]) bitset(i,j,2);
else
if(greys[i][j]>greylev[4]) bitset(i,j,4);
else bitset(i,j,4);
}
}
}
int bitset(int i, int j, int glev)
{
switch(glev){
case 4:
bits[(i*2)-1][(j*2)-1]=1;
bits[(i*2)][(j*2)-1]=1;
bits[(i*2)-1][(j*2)]=1;
bits[(i*2)][(j*2)]=1;
break;
case 3:
bits[(i*2)-1][(j*2)-1]=0;
bits[(i*2)][(j*2)-1]=1;
bits[(i*2)-1][(j*2)]=1;
bits[(i*2)][(j*2)]=1;
break;
case 2:
bits[(i*2)-1][(j*2)-1]=0;
bits[(i*2)][(j*2)-1]=1;
bits[(i*2)-1][(j*2)]=0;
bits[(i*2)][(j*2)]=1;
break;
case 1:
bits[(i*2)-1][(j*2)-1]=0;
bits[(i*2)][(j*2)-1]=0;
bits[(i*2)-1][(j*2)]=1;
bits[(i*2)][(j*2)]=0;
break;
case 0:
bits[(i*2)-1][(j*2)-1]=0;
bits[(i*2)][(j*2)-1]=0;
bits[(i*2)-1][(j*2)]=0;
bits[(i*2)][(j*2)]=0;
break;
}
}
int getval(float highval, float lowval, float val)
{
int i,j;
float greylev[5];
float valrange=highval-lowval;
for(i=0; i<5; i++)
{
greylev[i]=((i+1)*(valrange/5))+lowval;
}
for(i=0; i
for(j=0; j
if(val>greylev[3]) {printf("Returning 5\n");return 5;}
else
if(val>greylev[2]) {printf("Returning 4\n");return 4;}
else
if(val>greylev[1]) return 3;
else
if(val>greylev[0]) return 2;
else return 1;
}
}
}
float calculations(int i, int j)
{
float a,b,c;
float dx,dy,dz;
ywindow=((i-(ymax/2))/(ymax/2));
xwindow=((j-(xmax/2))/(xmax/2));
/* Normalise Direction Vector */
dx=xwindow-xeye;
dy=ywindow-yeye;
dz=zwindow-zeye;
/* Calculate M */
am1=(xeye-x)*dirx;
am2=(yeye-y)*diry;
am3=(zeye-z)*dirz;
am=am1+am2+am3;
bm1=dirx*dirx;
bm2=diry*diry;
bm3=dirz*dirz;
bm=bm1+bm2+bm3;
m=am/bm;
an1=dx*dirx;
an2=dy*diry;
an3=dz*dirz;
an=an1+an2+an3;
n=an/bm;
md1=m*dirx;
md2=m*diry;
md3=m*dirz;
g1=(xeye-x-(md1));
g2=(yeye-y-(md2));
g3=(zeye-z-(md3));
h1=dx-(n*dirx);
h2=dy-(n*diry);
h3=dz-(n*dirz);
a=((g1*g2)+(g2*g2)+(g3*g3))-(r*r);
b=2*((g1*h1)+(g2*h2)+(g3*h3));
c=((h1*h1)+(h2*h2)+(h3*h3));
t=quadratic(a,b,c);
k=m+(n*t);
return t;
}
float quadratic(float a, float b, float c)
{
float p,q,square;
square=(b*b)-(4*a*c);
if (square<0) return 0;
p=((0-b)+(sqrt((b*b)-(4*a*c))))/(2*a);
q=((0-b)-(sqrt((b*b)-(4*a*c))))/(2*a);
if (p>q) return q;
else return p;
}
e-mail: nolan.harley@ntlworld.com