added cbox creation tool

[physik/nlsop.git] / nlsop_create_cbox.c

diff --git a/nlsop_create_cbox.c b/nlsop_create_cbox.c
new file mode 100644 (file)
index 0000000..962abc3
--- /dev/null
+++ b/nlsop_create_cbox.c
@@ -0,0 +1,102 @@
+#define _GNU_SOURCE
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <fcntl.h>
+#include <unistd.h>
+
+#include <math.h>
+
+#include <signal.h>
+
+#include "nlsop.h"
+#include "dfbapi.h"
+
+int main(int argc,char **argv) {
+
+       int fd;
+       double atp;
+       int area,size,x,y,z,d;
+       int i,j,k;
+       d3_lattice d3l;
+       info info;
+       unsigned char *status;
+       int *conc;
+       double sigma;
+       int C,C0;
+
+       if(argc!=5) {
+               printf("usage: %s <x> <y> <z> <conc>\n",argv[0]);
+               return -23;
+       }
+
+       if((fd=open("cbox.save",O_WRONLY))<0) {
+               perror("open");
+               printf("unable to open save file!\n");
+               return -23;
+       }
+
+       x=atoi(argv[1]);
+       y=atoi(argv[2]);
+       z=atoi(argv[3]);
+       d=500/3; /* maximum 180keV, so ... */
+       atp=atof(argv[4]);
+       area=x*y;
+       size=area*z;
+       C=(1.0*atp*SI_PER_VOLUME)/(1-atp);
+       sigma=(1.0*((700-500)/3)*((700-500)/3))/(log(C));
+
+       printf("sanity checks ...\n");
+
+       if(z<d) {
+               printf("stupid!\n");
+               return -23;
+       }
+
+       if(atp>=1) {
+               printf("even more stupid!\n");
+               return -23;
+       }
+
+       printf("malloc ...\n");
+
+       status=malloc(size*sizeof(unsigned char));
+       memset(status,0,size*sizeof(unsigned char));
+       conc=malloc(size*sizeof(int));
+
+       printf("computing c profile ...\n");
+
+       i=area*(d+1);
+       C0=C;
+       printf("  first constant part -> %d\n",C0);
+       for(k=0;k<i;k++) *(conc+k)=C0;
+       printf("  last gauss part:\n");
+       for(k=d+1;k<z;k++) {
+               C=C0*exp(-1.0*(k-167)*(k-167)/sigma);
+               printf("    c = %d\n",C);
+               for(j=0;j<area;j++) *(conc+k*area+j)=C;
+       }
+
+       printf("nlsp header stuff ...\n");
+
+       d3l.max_x=x;
+       d3l.max_y=y;
+       d3l.max_z=z;
+       memset(&info,0,sizeof(info));
+
+       printf("writing file ...\n");
+
+       write(fd,&d3l,sizeof(d3_lattice));
+       write(fd,&info,sizeof(info));
+       write(fd,status,size*sizeof(unsigned char));
+       write(fd,conc,size*sizeof(int));
+
+       close(fd);
+
+       printf("done!\n");
+
+       return 1;
+}
+