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