fabmodules/src/scripts/math_png_py

657 wiersze
18 KiB
Python
Executable File

#!/usr/bin/env python
#
# math_png
# compile and evaluate .math to PNG
#
# Neil Gershenfeld
# CBA MIT 4/18/11
#
# (c) Massachusetts Institute of Technology 2010
# Permission granted for experimental and personal use;
# license for commercial sale available from MIT.
#
# todo
# rotated view
#
import sys
import string
import Image
import os
from numpy import *
#
# command line args
#
if (not((len(sys.argv) == 3) | (len(sys.argv) == 4) | (len(sys.argv) == 5) | (len(sys.argv) == 6)| (len(sys.argv) == 9))):
print "command line: math_png_py in.math out.png [resolution [number [view [rx ry rx]]]]"
print " in.math = input math string file"
print " out.png = output PNG image"
print " resolution = pixels per mm (optional, default 10)"
print " number = number of z slices to evaluate (optional, default 1)"
print " view = view projection(s) (optional, z|3, default z)"
print " rx ry rz = 3D view angle (optional, degrees, default 70 0 20)"
print "[This command is deprecated; use math_png instead.]"
sys.exit()
resolution = "10"
if (len(sys.argv) > 3):
resolution = sys.argv[3]
number = "1"
if (len(sys.argv) > 4):
number = sys.argv[4]
view = "z"
if (len(sys.argv) > 5):
view = sys.argv[5]
rx = "70"
ry = "10"
rz = "20"
if (len(sys.argv) > 6):
rx = sys.argv[6]
ry = sys.argv[7]
rz = sys.argv[8]
#
# read .math file
#
input_file_name = sys.argv[1]
output_file_name = sys.argv[2]
input_file = open(input_file_name,'r')
file_type = string.split(input_file.readline(), ': ')[1][:-1]
units = float(string.split(input_file.readline(), ': ')[1][:-1])
dx,dy,dz = string.split(input_file.readline())[3:]
xmin,ymin,zmin = string.split(input_file.readline())[3:]
math_string = input_file.readline()[1][:-1]
input_file.close()
print "read "+input_file_name
print " math string type: "+file_type
print " ",units,"mm per unit"
print " dx: "+dx+", dy: "+dy+", dz: "+dz
print " xmin: "+xmin+", ymin: "+ymin+", zmin: "+zmin
#
# write evaluation program
#
program = """//
// math_png_eval.c
// math string to PNG evaluation generated by math_png
//
// Neil Gershenfeld
// CBA MIT 3/6/11
//
// (c) Massachusetts Institute of Technology 2010
// Permission granted for experimental and personal use;
// license for commercial sale available from MIT.
//
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <stdint.h>
#include <png.h>
#include <string.h>
//
// data structures
//
struct fab_vars {
unsigned char empty;
unsigned char interior;
unsigned char edge;
unsigned char north;
unsigned char south;
unsigned char east;
unsigned char west;
unsigned char stop;
unsigned char corner;
unsigned char corner2;
unsigned char direction;
unsigned int nx,ny,nz;
unsigned int bit_depth;
unsigned int word_size;
double dx,dy,dz;
double xmin,ymin,zmin;
uint8_t **red,**green,**blue;
uint16_t **intensity;
short int **xptr, **yptr;
struct fab_path_type *path;
struct fab_mesh_type *mesh;
png_bytep *row_pointers;
png_structp png_ptr;
png_infop info_ptr;
};
void fab_write_png_K16(struct fab_vars *v, char *output_file_name) {
//
// write 16-bit grayscale PNG from fab_vars
//
FILE *output_file;
int x,y;
png_uint_32 res_x,res_y;
int unit_type;
png_byte color_type;
png_byte bit_depth;
png_byte *ptr;
//
// open PNG file
//
output_file = fopen(output_file_name, "wb");
v->png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING,
NULL, NULL, NULL);
v->info_ptr = png_create_info_struct(v->png_ptr);
png_init_io(v->png_ptr, output_file);
//
// set vars
//
bit_depth = 16;
color_type = PNG_COLOR_TYPE_GRAY;
png_set_IHDR(v->png_ptr, v->info_ptr, v->nx, v->ny,
bit_depth, color_type, PNG_INTERLACE_NONE,
PNG_COMPRESSION_TYPE_BASE, PNG_FILTER_TYPE_BASE);
res_x = 1000 * v->nx / v->dx;
res_y = 1000 * v->ny / v->dy;
png_set_pHYs(v->png_ptr, v->info_ptr, res_x, res_y, PNG_RESOLUTION_METER);
png_write_info(v->png_ptr, v->info_ptr);
//
// allocate pixels
//
v->row_pointers = (png_bytep*) malloc(sizeof(png_bytep) * v->ny);
for (y = 0; y < v->ny; ++y)
v->row_pointers[y] = (png_byte*) malloc(v->info_ptr->rowbytes);
//
// set pixels
//
for (y = 0; y < v->ny; ++y)
for (x = 0; x < v->nx; ++x) {
ptr = &(v->row_pointers[y][x*2]);
ptr[0] = (v->intensity[y][x] >> 8) & 255;
ptr[1] = v->intensity[y][x] & 255;
}
//
// write, close, and return
//
png_write_image(v->png_ptr, v->row_pointers);
png_write_end(v->png_ptr, NULL);
fclose(output_file);
printf("write %s\\n",output_file_name);
printf(" x pixels: %d, y pixels: %d\\n",v->nx,v->ny);
printf(" x pixels/m: %d, y pixels/m: %d\\n",(int) res_x,(int) res_y);
printf(" dx: %f mm, dy: %f mm\\n",v->dx,v->dy);
}
void fab_write_png_RGB24(struct fab_vars *v, char *output_file_name) {
//
// write 24-bit RGB PNG from fab_vars
//
FILE *output_file;
int x,y;
png_uint_32 res_x,res_y;
int unit_type;
png_byte color_type;
png_byte bit_depth;
png_byte *ptr;
//
// open PNG file
//
output_file = fopen(output_file_name, "wb");
v->png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING,
NULL, NULL, NULL);
v->info_ptr = png_create_info_struct(v->png_ptr);
png_init_io(v->png_ptr, output_file);
//
// set vars
//
bit_depth = 8;
color_type = PNG_COLOR_TYPE_RGB;
png_set_IHDR(v->png_ptr, v->info_ptr, v->nx, v->ny,
bit_depth, color_type, PNG_INTERLACE_NONE,
PNG_COMPRESSION_TYPE_BASE, PNG_FILTER_TYPE_BASE);
res_x = 1000 * v->nx / v->dx;
res_y = 1000 * v->ny / v->dy;
png_set_pHYs(v->png_ptr, v->info_ptr, res_x, res_y, PNG_RESOLUTION_METER);
png_write_info(v->png_ptr, v->info_ptr);
//
// allocate pixels
//
v->row_pointers = (png_bytep*) malloc(sizeof(png_bytep) * v->ny);
for (y = 0; y < v->ny; ++y)
v->row_pointers[y] = (png_byte*) malloc(v->info_ptr->rowbytes);
//
// set pixels
//
for (y = 0; y < v->ny; ++y)
for (x = 0; x < v->nx; ++x) {
ptr = &(v->row_pointers[y][x*3]);
ptr[0] = v->red[y][x];
ptr[1] = v->green[y][x];
ptr[2] = v->blue[y][x];
}
//
// write, close, and return
//
png_write_image(v->png_ptr, v->row_pointers);
png_write_end(v->png_ptr, NULL);
fclose(output_file);
printf("write %s\\n",output_file_name);
printf(" x pixels: %d, y pixels: %d\\n",v->nx,v->ny);
printf(" x pixels/m: %d, y pixels/m: %d\\n",(int) res_x,(int) res_y);
printf(" dx: %f mm, dy: %f mm\\n",v->dx,v->dy);
}
struct fab_vars init_vars() {
//
// fab_vars constructor
//
struct fab_vars v;
//
v.empty = 0;
v.interior = 1;
v.edge = (1 << 1);
v.north = (1 << 2);
v.west = (2 << 2);
v.east = (3 << 2);
v.south = (4 << 2);
v.stop = (5 << 2);
v.corner = (6 << 2);
v.corner2 = (7 << 2);
v.direction = (7 << 2);
//
return v;
}
"""
if (view == 'z'):
program +="""
main(int argc, char **argv) {
//
// local vars
//
struct fab_vars v = init_vars(&v);
double units = 25.4;
double X,Y,Z;
v.dx = units*"""+dx+""";
v.dy = units*"""+dy+""";
v.dz = units*"""+dz+""";
v.xmin = units*"""+xmin+""";
v.ymin = units*"""+ymin+""";
v.zmin = units*"""+zmin+""";
double resolution="""+resolution+""";
v.nx = resolution*v.dx;
v.ny = resolution*v.dy;
v.nz = """+number+""";
int x,y,z;
uint32_t val;
uint8_t byte;
float scale;
//
// create arrays
//
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
v.intensity = malloc(v.ny*sizeof(uint16_t *));
for (y = 0; y < v.ny; ++y) {
v.intensity[y] = (uint16_t*) malloc(v.nx*sizeof(uint16_t));
}
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
v.red = malloc(v.ny*sizeof(uint16_t *));
v.green = malloc(v.ny*sizeof(uint16_t *));
v.blue = malloc(v.ny*sizeof(uint16_t *));
for (y = 0; y < v.ny; ++y) {
v.red[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
v.green[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
v.blue[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
}
}
else {
printf("math_png: oops -- file type not recognized\\n");
exit(-1);
}
//
// evaluate
//
for (z = 0; z < v.nz; ++z) {
if (v.nz == 1) {
Z = v.zmin/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535;
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = 1;
}
else {
Z = (v.zmin + v.dz*z/(v.nz-1.0))/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535 * z/(v.nz-1.0);
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = z/(v.nz-1.0);
}
printf(" Z = %f\\n",Z);
for (y = 0; y < v.ny; ++y) {
Y = (v.ymin + v.dy*(v.ny-y-1)/(v.ny-1.0))/units;
for (x = 0; x < v.nx; ++x) {
X = (v.xmin + v.dx*x/(v.nx-1.0))/units;
val = """+math_string+""";
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
val = scale*val;
if (val > v.intensity[y][x])
v.intensity[y][x] = val;
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
byte = scale * (val & 255);
if (val > v.red[y][x])
v.red[y][x] = byte;
byte = scale * ((val >> 8) & 255);
if (val > v.green[y][x])
v.green[y][x] = byte;
byte = scale * ((val >> 16) & 255);
if (val > v.blue[y][x])
v.blue[y][x] = byte;
}
}
}
}
//
// write PNG
//
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
fab_write_png_K16(&v,\""""+output_file_name+"""\");
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
fab_write_png_RGB24(&v,\""""+output_file_name+"""\");
}
"""
elif (view == 'r'):
program +="""
main(int argc, char **argv) {
//
// to be completed
//
// local vars
//
struct fab_vars v = init_vars(&v);
double units = 25.4;
double X,Y,Z,X0,Y0,Z0,X1,Y1,Z1;
v.dx = units*"""+dx+""";
v.dy = units*"""+dy+""";
v.dz = units*"""+dz+""";
v.xmin = units*"""+xmin+""";
v.ymin = units*"""+ymin+""";
v.zmin = units*"""+zmin+""";
double resolution="""+resolution+""";
int nx = resolution*v.dx;
int ny = resolution*v.dy;
int nz = resolution*v.dz;
v.nx = resolution*(v.dx);
v.ny = resolution*(v.dy);
v.nz = """+number+""";
int x,y,z;
uint32_t val;
int scale;
double pi = 3.14159;
double rx = pi*"""+rx+"""/180.0;
double ry = pi*"""+ry+"""/180.0;
double rz = pi*"""+rz+"""/180.0;
//
// create RGB array
//
v.rgb = malloc(v.ny*sizeof(uint32_t *));
for (y = 0; y < v.ny; ++y)
v.rgb[y] = (uint32_t*) malloc(2*v.nx*sizeof(uint32_t));
//
// evaluate xyz
//
for (z = 0; z < v.nz; ++z) {
if (v.nz == 1) {
Z0 = v.zmin/units;
scale = 1;
}
else {
Z0 = (v.zmin + v.dz*z/(v.nz-1.0))/units;
scale = 255 * z/(v.nz-1.0);
}
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = scale + (scale << 8) + (scale << 16);
printf(" Z = %f\\n",Z0);
for (y = 0; y < v.ny; ++y) {
Y0 = (v.ymin + v.dy*(v.ny-y-1)/(v.ny-1.0))/units;
for (x = 0; x < v.nx; ++x) {
X0 = (v.xmin + v.dx*x/(v.nx-1.0))/units;
X1 = X0;
Y1 = cos(rx)*Y0 - sin(rx)*Z0;
Z1 = sin(rx)*Y0 + cos(rx)*Z0;
X = cos(rz)*X1 - sin(rz)*Y1;
Y = sin(rz)*X1 + cos(rz)*Y1;
Z = Z1;
val = """+math_string+""";
val = scale*val;
if (val > v.rgb[y][x]) {
v.rgb[y][x] = val;
}
}
}
}
//
// write PNG
//
fab_write_png(&v,\""""+output_file_name+"""\");
}
"""
elif (view == '3'):
program +="""
main(int argc, char **argv) {
//
// local vars
//
struct fab_vars v = init_vars(&v);
double units = 25.4;
double X,Y,Z,X0,Y0,Z0,X1,Y1,Z1;
v.dx = units*"""+dx+""";
v.dy = units*"""+dy+""";
v.dz = units*"""+dz+""";
v.xmin = units*"""+xmin+""";
v.ymin = units*"""+ymin+""";
v.zmin = units*"""+zmin+""";
double resolution="""+resolution+""";
int nx = resolution*v.dx;
int ny = resolution*v.dy;
int nz = resolution*v.dz;
int border = nx/25;
v.nx = border + nx + nz;
v.ny = border + ny + nz;
v.nz = """+number+""";
int x,y,z;
uint32_t val;
uint8_t byte;
float scale;
//
// create arrays
//
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
v.intensity = malloc(v.ny*sizeof(uint16_t *));
for (y = 0; y < v.ny; ++y) {
v.intensity[y] = (uint16_t*) malloc(v.nx*sizeof(uint16_t));
}
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
v.red = malloc(v.ny*sizeof(uint16_t *));
v.green = malloc(v.ny*sizeof(uint16_t *));
v.blue = malloc(v.ny*sizeof(uint16_t *));
for (y = 0; y < v.ny; ++y) {
v.red[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
v.green[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
v.blue[y] = (uint8_t*) malloc(v.nx*sizeof(uint8_t));
}
}
else {
printf("math_png: oops -- file type not recognized\\n");
exit(-1);
}
//
// evaluate xy
//
for (z = 0; z < v.nz; ++z) {
if (v.nz == 1) {
Z = v.zmin/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535;
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = 1;
}
else {
Z = (v.zmin + v.dz*z/(v.nz-1.0))/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535 * z/(v.nz-1.0);
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = z/(v.nz-1.0);
}
printf(" Z = %f\\n",Z);
for (y = 0; y < ny; ++y) {
Y = (v.ymin + v.dy*(ny-y-1)/(ny-1.0))/units;
for (x = 0; x < nx; ++x) {
X = (v.xmin + v.dx*x/(nx-1.0))/units;
val = """+math_string+""";
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
val = scale*val;
if (val > v.intensity[y][x])
v.intensity[y][x] = val;
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
byte = scale * (val & 255);
if (val > v.red[y][x])
v.red[y][x] = byte;
byte = scale * ((val >> 8) & 255);
if (val > v.green[y][x])
v.green[y][x] = byte;
byte = scale * ((val >> 16) & 255);
if (val > v.blue[y][x])
v.blue[y][x] = byte;
}
}
}
}
//
// evaluate xz
//
for (y = 0; y < v.nz; ++y) {
if (v.nz == 1) {
Y = v.ymin/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535;
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = 1;
}
else {
Y = (v.ymin + v.dy*(v.nz-y-1)/(v.nz-1.0))/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535 * y/(v.nz-1.0);
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = y/(v.nz-1.0);
}
printf(" Y = %f\\n",Y);
for (z = 0; z < nz; ++z) {
Z = (v.zmin + v.dz*(nz-z-1)/(nz-1.0))/units;
for (x = 0; x < nx; ++x) {
X = (v.xmin + v.dx*x/(nx-1.0))/units;
val = """+math_string+""";
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
val = scale*val;
if (val > v.intensity[ny+border+z][x])
v.intensity[ny+border+z][x] = val;
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
byte = scale * (val & 255);
if (val > v.red[y][x])
v.red[y][x] = byte;
byte = scale * ((val >> 8) & 255);
if (val > v.green[y][x])
v.green[y][x] = byte;
byte = scale * ((val >> 16) & 255);
if (val > v.blue[y][x])
v.blue[y][x] = byte;
}
}
}
}
//
// evaluate zy
//
for (x = 0; x < v.nz; ++x) {
if (v.nz == 1) {
X = v.xmin/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535;
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = 1;
}
else {
X = (v.xmin + v.dx*x/(v.nz-1.0))/units;
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
scale = 65535 * x/(v.nz-1.0);
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
scale = x/(v.nz-1.0);
}
printf(" X = %f\\n",X);
for (y = 0; y < ny; ++y) {
Y = (v.ymin + v.dy*(ny-y-1)/(ny-1.0))/units;
for (z = 0; z < nz; ++z) {
Z = (v.zmin + v.dz*(nz-z-1)/(nz-1.0))/units;
val = """+math_string+""";
if (0 == strcmp(\""""+file_type+"""\","Boolean")) {
val = scale*val;
if (val > v.intensity[y][nx+border+z])
v.intensity[y][nx+border+z] = val;
}
else if (0 == strcmp(\""""+file_type+"""\","RGB")) {
byte = scale * (val & 255);
if (val > v.red[y][x])
v.red[y][x] = byte;
byte = scale * ((val >> 8) & 255);
if (val > v.green[y][x])
v.green[y][x] = byte;
byte = scale * ((val >> 16) & 255);
if (val > v.blue[y][x])
v.blue[y][x] = byte;
}
}
}
}
//
// write PNG
//
if (0 == strcmp(\""""+file_type+"""\","Boolean"))
fab_write_png_K16(&v,\""""+output_file_name+"""\");
else if (0 == strcmp(\""""+file_type+"""\","RGB"))
fab_write_png_RGB24(&v,\""""+output_file_name+"""\");
}
"""
else:
print "view "+view+" not recognized"
sys.exit()
eval_file = open("math_png_eval.c",'w')
eval_file.write(program)
eval_file.close()
#
# compile
#
print "compile"
if os.uname()[0] == 'Darwin':
os.system("gcc -O math_png_eval.c -o math_png_eval -lm -I/usr/X11/include -L/usr/X11/lib -lpng")
else:
os.system("gcc -O math_png_eval.c -o math_png_eval -lm -lpng")
#
# execute
#
print "execute"
os.system("./math_png_eval")
#os.system("rm ./math_png_eval.c")
os.system("rm ./math_png_eval")