trisurf-ng.git - Gitblit

Found one bug in the code when calculating mprod. Fixed.

Samo Penic

2021-07-27 3197854e89eb9ffe5ff3137e274d16806c2028a5

commit \| author \| age
7f6076	1	/* vim: set ts=4 sts=4 sw=4 noet : */
SP	2	/* vim: set ts=4 sts=4 sw=4 noet : */
854cb6	3	#include<stdio.h>
SP	4	#include<general.h>
	5	#include<snapshot.h>
22721d	6	#include<stdlib.h>
SP	7	#include<stdarg.h>
	8	#include <zlib.h>
f0bcea	9	#include<inttypes.h>
11fef8	10	#include<config.h>
2e8c7f	11	#include <time.h>
698ae1	12	#include "io.h"
a3dbf0	13	/* a helper function that utilizes ts_string data structure and performs same as sprintf */
22721d	14	ts_uint ts_sprintf(ts_string str, char fmt, ...){
SP	15	va_list ap;
	16	va_start(ap,fmt);
	17	ts_uint n=vsprintf(&(str->string[str->beg]),fmt,ap);
	18	va_end(ap);
	19	str->beg+=n;
	20	return n;
	21	}
	22
a3dbf0	23	/* outputs additional data into paraview xml file */
854cb6	24	ts_bool xml_trisurf_data(FILE fh, ts_vesicle vesicle){
22721d	25
SP	26	ts_string data=(ts_string )malloc(sizeof(ts_sprintf));
02d65c	27	data->string=(char )malloc(5120000sizeof(char)); /TODO: warning, can break if the string is to long /
22721d	28	data->beg=0;
SP	29
854cb6	30	xml_trisurf_header(fh, vesicle);
22721d	31	xml_trisurf_tria(data,vesicle->tlist);
SP	32	xml_trisurf_tria_neigh(data,vesicle->tlist);
	33	xml_trisurf_vtx_neigh(data,vesicle->vlist);
	34	xml_trisurf_vtx_tristar(data,vesicle->vlist);
4891eb	35	xml_trisurf_nucleus(data,vesicle);
49981c	36	xml_trisurf_constvolarea(data,V0,A0);
22721d	37	#ifdef COMPRESSION
f0bcea	38	char *compressed;
588bbb	39	ts_uint nbytes=ts_compress_string64(data->string, data->beg-1, &compressed); //suppress null character at the end with by substracting 1
f0bcea	40	fwrite(compressed, sizeof(unsigned char), nbytes, fh);
SP	41	free (compressed);
22721d	42	#else
SP	43	fprintf(fh,"%s", data->string);
	44	#endif
a3dbf0	45	free(data->string); /* TODO: valgrind is not ok with this! */
22721d	46	free(data);
854cb6	47	xml_trisurf_footer(fh);
SP	48	return TS_SUCCESS;
	49	}
	50
	51	ts_bool xml_trisurf_header(FILE fh, ts_vesicle vesicle){
2e8c7f	52	/* format current time */
SP	53	time_t current_time;
	54	char *c_time_string;
	55	current_time = time(NULL);
	56	c_time_string = ctime(&current_time);
a3dbf0	57	int npoly, nfono;
SP	58
e29756	59	fprintf(fh, "<trisurfversion>Trisurf (commit %s), compiled on %s %s</trisurfversion>\n",TS_VERSION, __DATE__, __TIME__);
2e8c7f	60	fprintf(fh, "<dumpdate>%s</dumpdate>\n", c_time_string);
SP	61
cc95bd	62	fprintf(fh, "<tape>");
698ae1	63	fprintf(fh,"%s",tapetxt);
f0bcea	64	fprintf(fh, "</tape>\n");
a3dbf0	65	if(vesicle->poly_list!=NULL){
SP	66	npoly=vesicle->poly_list->n;
	67	if(npoly!=0){
	68	nfono=vesicle->poly_list->poly[0]->vlist->n;
	69	} else {
	70	nfono=0;
	71	}
	72	} else {
	73	npoly=0;
	74	nfono=0;
	75	}
86b69b	76	fprintf(fh, "<trisurf nvtx=\"%u\" npoly=\"%u\" nmono=\"%u\" compressed=\"false\">\n", vesicle->vlist->n, npoly, nfono);
854cb6	77	return TS_SUCCESS;
SP	78	}
	79
	80	ts_bool xml_trisurf_footer(FILE *fh){
	81	fprintf(fh, "</trisurf>\n");
	82	return TS_SUCCESS;
	83	}
	84
22721d	85	ts_bool xml_trisurf_tria(ts_string data, ts_triangle_list tlist){
854cb6	86	ts_uint i;
588bbb	87	ts_sprintf(data,"<tria>");
854cb6	88	for(i=0; i<tlist->n;i++){
e87e4e	89	ts_sprintf(data,"%u %u %u ",tlist->tria[i]->vertex[0]->idx, tlist->tria[i]->vertex[1]->idx, tlist->tria[i]->vertex[2]->idx);
854cb6	90	}
588bbb	91	ts_sprintf(data,"</tria>");
854cb6	92	return TS_SUCCESS;
SP	93	}
	94
22721d	95	ts_bool xml_trisurf_tria_neigh(ts_string data, ts_triangle_list tlist){
854cb6	96	ts_uint i;
22721d	97	ts_sprintf(data,"<trianeigh>\n");
854cb6	98	for(i=0; i<tlist->n;i++){
e87e4e	99	ts_sprintf(data,"%u %u %u ",tlist->tria[i]->neigh[0]->idx, tlist->tria[i]->neigh[1]->idx, tlist->tria[i]->neigh[2]->idx);
854cb6	100	}
22721d	101	ts_sprintf(data,"</trianeigh>\n");
854cb6	102	return TS_SUCCESS;
SP	103	}
	104
22721d	105	ts_bool xml_trisurf_vtx_neigh(ts_string data, ts_vertex_list vlist){
854cb6	106	ts_uint i,j;
SP	107	for(i=0;i<vlist->n;i++){
22721d	108	ts_sprintf(data,"<vtxn idx=\"%u\">",vlist->vtx[i]->idx);
854cb6	109	for(j=0;j<vlist->vtx[i]->neigh_no;j++){
22721d	110	ts_sprintf(data,"%u ",vlist->vtx[i]->neigh[j]->idx);
854cb6	111	}
588bbb	112	ts_sprintf(data, "</vtxn>");
854cb6	113	}
SP	114	return TS_SUCCESS;
	115	}
	116
22721d	117	ts_bool xml_trisurf_vtx_tristar(ts_string data, ts_vertex_list vlist){
854cb6	118	ts_uint i,j;
SP	119	for(i=0;i<vlist->n;i++){
22721d	120	ts_sprintf(data,"<tristar idx=\"%u\">",vlist->vtx[i]->idx);
854cb6	121	for(j=0;j<vlist->vtx[i]->tristar_no;j++){
22721d	122	ts_sprintf(data,"%u ",vlist->vtx[i]->tristar[j]->idx);
854cb6	123	}
588bbb	124	ts_sprintf(data, "</tristar>");
854cb6	125	}
SP	126	return TS_SUCCESS;
	127	}
	128
4891eb	129	ts_bool xml_trisurf_nucleus(ts_string data, ts_vesicle vesicle){
SP	130	if(vesicle->R_nucleus>0.0 \|\| (vesicle->R_nucleusX>0.0 && vesicle->R_nucleusY>0.0 && vesicle->R_nucleusZ>0.0)){
49981c	131	ts_sprintf(data,"<nucleus>%.17e %.17e %.17e</nucleus>",vesicle->nucleus_center[0], vesicle->nucleus_center[1], vesicle->nucleus_center[2]);
4891eb	132	}
SP	133	return TS_SUCCESS;
	134	}
49981c	135	ts_bool xml_trisurf_constvolarea(ts_string *data, ts_double volume, ts_double area){
SP	136	ts_sprintf(data,"<constant_volume>%.17e</constant_volume>",volume);
	137	ts_sprintf(data,"<constant_area>%.17e</constant_area>",area);
	138
	139	return TS_SUCCESS;
	140	}
f0bcea	141
SP	142
	143	/* UTILITIES */
	144
	145	/* zlib compression base64 encoded */
	146	/* compressed must not be pre-malloced */
8c1bb1	147	/* taken from https://gist.github.com/arq5x/5315739 */
f0bcea	148	ts_uint ts_compress_string64(char data, ts_uint data_len, char *compressed){
SP	149	z_stream defstream;
	150	defstream.zalloc = Z_NULL;
	151	defstream.zfree = Z_NULL;
	152	defstream.opaque = Z_NULL;
	153	defstream.avail_in = data_len+1;
	154	defstream.next_in = (unsigned char *)data;
	155	char compr=(char )malloc(data_lensizeof(char ));
	156	defstream.avail_out = data_len+1;
	157	defstream.next_out = (unsigned char *)compr;
	158	deflateInit(&defstream, Z_BEST_COMPRESSION);
	159	deflate(&defstream, Z_FINISH);
	160	deflateEnd(&defstream);
	161	/base64 encode/
	162	size_t nbase;
	163	compressed=base64_encode((unsigned char )compr,(size_t)defstream.total_out,&nbase);
	164	//fwrite(base64, sizeof(unsigned char), nbase, fh);
	165	free(compr);
	166	return nbase;
	167	}
	168
8c1bb1	169	ts_uint ts_decompress_string64(char b64, ts_uint data_len, char *decompressed){
SP	170	return TS_SUCCESS;
	171
	172	}
f0bcea	173
SP	174	/* base64 encoding, taken from http://stackoverflow.com/questions/342409/how-do-i-base64-encode-decode-in-c */
	175	static char encoding_table[] = {'A', 'B', 'C', 'D', 'E', 'F', 'G', 'H',
	176	'I', 'J', 'K', 'L', 'M', 'N', 'O', 'P',
	177	'Q', 'R', 'S', 'T', 'U', 'V', 'W', 'X',
	178	'Y', 'Z', 'a', 'b', 'c', 'd', 'e', 'f',
	179	'g', 'h', 'i', 'j', 'k', 'l', 'm', 'n',
	180	'o', 'p', 'q', 'r', 's', 't', 'u', 'v',
	181	'w', 'x', 'y', 'z', '0', '1', '2', '3',
	182	'4', '5', '6', '7', '8', '9', '+', '/'};
	183	static char *decoding_table = NULL;
	184	static int mod_table[] = {0, 2, 1};
	185
	186
7f6076	187	char base64_encode(const unsigned char data, size_t input_length, size_t *output_length) {
f0bcea	188	output_length = 4 ((input_length + 2) / 3);
SP	189	int i,j;
	190	char encoded_data = malloc(output_length);
	191	if (encoded_data == NULL) return NULL;
	192
	193	for (i = 0, j = 0; i < input_length;) {
	194
7f6076	195	uint32_t octet_a = i < input_length ? (unsigned char)data[i++] : 0;
SP	196	uint32_t octet_b = i < input_length ? (unsigned char)data[i++] : 0;
	197	uint32_t octet_c = i < input_length ? (unsigned char)data[i++] : 0;
	198	uint32_t triple = (octet_a << 0x10) + (octet_b << 0x08) + octet_c;
f0bcea	199
7f6076	200	encoded_data[j++] = encoding_table[(triple >> 3 * 6) & 0x3F];
SP	201	encoded_data[j++] = encoding_table[(triple >> 2 * 6) & 0x3F];
	202	encoded_data[j++] = encoding_table[(triple >> 1 * 6) & 0x3F];
	203	encoded_data[j++] = encoding_table[(triple >> 0 * 6) & 0x3F];
f0bcea	204	}
SP	205
	206	for (i = 0; i < mod_table[input_length % 3]; i++)
7f6076	207	encoded_data[*output_length - 1 - i] = '=';
f0bcea	208
SP	209	return encoded_data;
	210	}
	211
	212
7f6076	213	unsigned char base64_decode(const char data, size_t input_length, size_t *output_length) {
f0bcea	214	int i,j;
SP	215	if (decoding_table == NULL) build_decoding_table();
	216
	217	if (input_length % 4 != 0) return NULL;
	218
	219	output_length = input_length / 4 3;
	220	if (data[input_length - 1] == '=') (*output_length)--;
	221	if (data[input_length - 2] == '=') (*output_length)--;
	222
	223	unsigned char decoded_data = malloc(output_length);
	224	if (decoded_data == NULL) return NULL;
	225
	226	for (i = 0, j = 0; i < input_length;) {
7f6076	227	uint32_t sextet_a = data[i] == '=' ? 0 & i++ : decoding_table[(int)data[i++]];
SP	228	uint32_t sextet_b = data[i] == '=' ? 0 & i++ : decoding_table[(int)data[i++]];
	229	uint32_t sextet_c = data[i] == '=' ? 0 & i++ : decoding_table[(int)data[i++]];
	230	uint32_t sextet_d = data[i] == '=' ? 0 & i++ : decoding_table[(int)data[i++]];
f0bcea	231
7f6076	232	uint32_t triple = (sextet_a << 3 * 6)
f0bcea	233	+ (sextet_b << 2 * 6)
SP	234	+ (sextet_c << 1 * 6)
	235	+ (sextet_d << 0 * 6);
	236
7f6076	237	if (j < output_length) decoded_data[j++] = (triple >> 2 8) & 0xFF;
SP	238	if (j < output_length) decoded_data[j++] = (triple >> 1 8) & 0xFF;
	239	if (j < output_length) decoded_data[j++] = (triple >> 0 8) & 0xFF;
f0bcea	240	}
8c1bb1	241	if(decoding_table !=NULL) free(decoding_table);
f0bcea	242	return decoded_data;
SP	243	}
	244
	245
	246	void build_decoding_table() {
	247
	248	decoding_table = malloc(256);
	249	int i;
	250	for (i = 0; i < 64; i++)
7f6076	251	decoding_table[(unsigned char) encoding_table[i]] = i;
f0bcea	252	}
SP	253
	254
	255	void base64_cleanup() {
	256	free(decoding_table);
	257	}
	258
	259
	260
	261
	262
	263
	264