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