ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/gen/genglaze.c
Revision: 2.2
Committed: Tue Jun 10 23:13:13 2025 UTC (3 days, 20 hours ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 2.1: +37 -9 lines
Log Message:
fix(genglaze): TW fixed a memory leak and some minor improvements

File Contents

# User Rev Content
1 greg 2.1 /*
2     Glazing system multi-layer optics calculation.
3     Source equations: LBNL WINDOW technical documention.
4     System calcuation: Section 7.6.1
5     Anuglar calculation Section 7.7
6     T. Wang
7     */
8    
9     #include "rtmath.h"
10     #include "color.h"
11     #include "data.h"
12    
13     enum {
14     NTHETA = 59,
15     MAX_NAME = 256,
16     };
17    
18     const double DBL_EPSILON = 1E-9;
19 greg 2.2 const double LARGE_DOUBE = 3.40282347E+38;
20 greg 2.1 const double RAD_PER_DEG = 0.0174532925199;
21     const double COEFFS_TAU_CLEAR[] = {-0.0015, 3.355, -3.84, 1.46, 0.0288};
22     const double COEFFS_TAU_BRONZE[] = {0.002, 2.813, -2.341, -0.05725, 0.599};
23     const double COEFFS_RHO_CLEAR[] = {0.999, -0.563, 2.043, -2.532, 1.054};
24     const double COEFFS_RHO_BRONZE[] = {0.997, -1.868, 6.513, -7.862, 3.225};
25    
26     const double THETAS[NTHETA] = {
27     0., 5., 10., 15., 20., 25., 30., 35., 40.,
28     41., 42., 43., 44., 45., 46., 47., 48., 49, 50.,
29     51., 52., 53., 54., 55., 56., 57., 58., 59, 60.,
30     61., 62., 63., 64., 65., 66., 67., 68., 69, 70.,
31     71., 72., 73., 74., 75., 76., 77., 78., 79, 80.,
32     81., 82., 83., 84., 85., 86., 87., 88., 89, 90.,
33     };
34    
35    
36     char *progname;
37    
38     typedef struct {
39     char *filename;
40     int is_mono;
41     double thickness_m;
42    
43     /* Interpolated data at standard wavelengths */
44     double *t_0;
45     double *rf_0;
46     double *rb_0;
47    
48     /* Calculated angular data */
49     double *t_lambda_theta;
50     double *rf_lambda_theta;
51     double *rb_lambda_theta;
52    
53     } GlazingLayer;
54    
55    
56 greg 2.2 static inline double polynomial_5(double x, const double coeffs[5]) {
57 greg 2.1 return x * ( x * ( x * (coeffs[4] * x + coeffs[3]) + coeffs[2]) + coeffs[1]) + coeffs[0];
58     }
59    
60    
61     /* Calculate rho0 (unpolarized reflectance of a single surface) */
62     static void rho0_calc(const double *t_0, const double *r_0, const int nwvl, double *rho0) {
63     for (int i = 0; i < nwvl; ++i) {
64     const double beta = t_0[i] * t_0[i] - r_0[i] * r_0[i] + 2.0 * r_0[i] + 1.0;
65     double denominator = 2.0 * (2.0 - r_0[i]);
66     if (fabs(denominator) < DBL_EPSILON)
67     denominator = DBL_EPSILON;
68     double discriminant = beta * beta - 2.0 * denominator * r_0[i];
69     if (discriminant < 0.0)
70     discriminant = 0.0;
71     rho0[i] = (beta - sqrt(discriminant)) / denominator;
72     if (rho0[i] < 0.0)
73     rho0[i] = 0.0;
74     if (rho0[i] > 1.0)
75     rho0[i] = 1.0;
76     }
77     }
78    
79     /* refraction index calc */
80     static void n_calc(const double *rho_0, const int nwvl, double *n)
81     {
82     for (int i=0; i < nwvl; ++i) {
83     const double sqrt_rho = sqrt(rho_0[i]);
84     n[i] = ((1.0 - sqrt_rho) < DBL_EPSILON) ? 1.0 : (1.0 + sqrt_rho) / (1.0 - sqrt_rho);
85     }
86     }
87    
88     /* Calculate absoprtion coefficient */
89     static void alpha_calc(const double *t0, const double *r0, const double *rho0,
90     const double thickness_m, const int nwvl, double *abs_coeff)
91     {
92     for (int i = 0; i < nwvl; ++i) {
93     const double numerator = r0[i] - rho0[i];
94     const double denominator = rho0[i] * t0[i];
95     if (denominator > DBL_EPSILON && numerator > DBL_EPSILON) {
96     const double log_val = log(numerator / denominator);
97     abs_coeff[i] = -log_val / thickness_m;
98     } else {
99     abs_coeff[i] = 1.0;
100     }
101     if (abs_coeff[i] < 0.0 || isnan(abs_coeff[i]) || isinf(abs_coeff[i])) {
102     abs_coeff[i] = 0.0;
103     }
104     }
105     }
106    
107    
108     void angular_monolithic(GlazingLayer *layer, const int nwvl) {
109     double *rho_0 = malloc(nwvl * sizeof(double));
110     double *refrac = malloc(nwvl * sizeof(double));
111     double *abs_coeff = malloc(nwvl * sizeof(double));
112    
113     rho0_calc(layer->t_0, layer->rf_0, nwvl, rho_0);
114     n_calc(rho_0, nwvl, refrac);
115     alpha_calc(layer->t_0, layer->rf_0, rho_0, layer->thickness_m, nwvl, abs_coeff);
116    
117     for (int itheta = 0; itheta < NTHETA; ++itheta) {
118     double phi = THETAS[itheta] * RAD_PER_DEG;
119     double cos_phi = cos(phi);
120     double sin_phi = sin(phi);
121    
122     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
123     double sin_phi_prime_sq = sin_phi * sin_phi / (refrac[iwvl] * refrac[iwvl]);
124    
125     double cos_phi_prime = (sin_phi_prime_sq >= 1.0) ? 0.0 : sqrt(1.0 - sin_phi_prime_sq);
126    
127     double rho_lambda;
128     if (cos_phi_prime <= DBL_EPSILON) {
129     rho_lambda = 1.0;
130     } else {
131     double nf_cos_phi = refrac[iwvl] * cos_phi;
132     double nf_cos_phi_p = refrac[iwvl] * cos_phi_prime;
133     double par_num = nf_cos_phi - cos_phi_prime;
134     double par_den = nf_cos_phi + cos_phi_prime;
135     double per_num = nf_cos_phi_p - cos_phi;
136     double per_den = nf_cos_phi_p + cos_phi;
137     if (fabs(par_den) < DBL_EPSILON || fabs(per_den) < DBL_EPSILON) {
138     rho_lambda = 1.0;
139     } else {
140     double par_comp = pow(par_num / par_den, 2.0);
141     double per_comp = pow(per_num / per_den, 2.0);
142     rho_lambda = 0.5 * (par_comp + per_comp);
143     }
144     }
145     if (rho_lambda < 0.0)
146     rho_lambda = 0.0;
147     if (rho_lambda > 1.0)
148     rho_lambda = 1.0;
149    
150 greg 2.2 double exp_arg = (cos_phi_prime < DBL_EPSILON) ? -LARGE_DOUBE : (-abs_coeff[iwvl] * layer->thickness_m / cos_phi_prime);
151 greg 2.1 /* Clamp exponent argument to prevent overflow in exp() */
152     if (exp_arg < -700.0)
153     exp_arg = -700.0;
154    
155     double T_internal = exp(exp_arg);
156    
157     double denominator = 1.0 - rho_lambda * rho_lambda * T_internal * T_internal;
158    
159     double t_lambda, r_lambda;
160    
161     if (fabs(denominator) < DBL_EPSILON || rho_lambda >= 1.0) {
162     t_lambda = 0.0;
163     r_lambda = 1.0;
164     } else {
165     double tau_lambda = 1.0 - rho_lambda;
166     t_lambda = tau_lambda * tau_lambda * T_internal / denominator;
167     r_lambda = rho_lambda * (1.0 + t_lambda * T_internal);
168    
169     }
170     if (t_lambda < 0.0)
171     t_lambda = 0.0;
172     if (t_lambda > 1.0)
173     t_lambda = 1.0;
174    
175     const int flat_idx = itheta * nwvl + iwvl;
176     layer->t_lambda_theta[flat_idx] = t_lambda;
177     layer->rf_lambda_theta[flat_idx] = r_lambda;
178     layer->rb_lambda_theta[flat_idx] = r_lambda;
179     }
180     }
181 greg 2.2 free(rho_0);
182     free(refrac);
183     free(abs_coeff);
184 greg 2.1 }
185    
186    
187     void angular_coated(GlazingLayer *layer, const int nwvl) {
188     for (int itheta = 0; itheta < NTHETA; ++itheta) {
189     const double phi = THETAS[itheta] * RAD_PER_DEG;
190     const double cos_phi = cos(phi);
191    
192     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
193     const double t_0_lambda = layer->t_0[iwvl];
194     const double *coeffs_tau = (t_0_lambda > 0.645) ? COEFFS_TAU_CLEAR : COEFFS_TAU_BRONZE;
195     const double *coeffs_rho = (t_0_lambda > 0.645) ? COEFFS_RHO_CLEAR : COEFFS_RHO_BRONZE;
196    
197     const double tau_bar = polynomial_5(cos_phi, coeffs_tau);
198     const double rho_bar_term = polynomial_5(cos_phi, coeffs_rho);
199     const double rho_bar = rho_bar_term - tau_bar;
200    
201     const int flat_idx = itheta * nwvl + iwvl;
202     layer->t_lambda_theta[flat_idx] = t_0_lambda * tau_bar;
203     layer->rf_lambda_theta[flat_idx] = layer->rf_0[iwvl] * (1.0 - rho_bar) + rho_bar;
204     layer->rb_lambda_theta[flat_idx] = layer->rb_0[iwvl] * (1.0 - rho_bar) + rho_bar;
205     }
206     }
207     }
208    
209    
210     void multi_layer_calc(
211     GlazingLayer *layers,
212     int nlayers,
213     int nwvl,
214     double *total_t,
215     double *total_rf,
216     double *total_rb)
217     {
218     if (nlayers <= 0)
219     return;
220    
221     size_t total_size = (size_t)NTHETA * nwvl * sizeof(double);
222    
223     memcpy(total_t, layers[0].t_lambda_theta, total_size);
224     memcpy(total_rf, layers[0].rf_lambda_theta, total_size);
225     memcpy(total_rb, layers[0].rb_lambda_theta, total_size);
226    
227     if (nlayers == 1)
228     return;
229    
230     double *prev_t = malloc(total_size);
231     double *prev_rf = malloc(total_size);
232     double *prev_rb = malloc(total_size);
233    
234     if (!prev_t || !prev_rf || !prev_rb)
235     perror("Failed to allocate temporary storage in multi_layer");
236    
237     for (int j = 1; j < nlayers; ++j) {
238     memcpy(prev_t, total_t, total_size);
239     memcpy(prev_rf, total_rf, total_size);
240     memcpy(prev_rb, total_rb, total_size);
241    
242     double *t_j = layers[j].t_lambda_theta;
243     double *rf_j = layers[j].rf_lambda_theta;
244     double *rb_j = layers[j].rb_lambda_theta;
245    
246     for (int itheta = 0; itheta < NTHETA; ++itheta) {
247     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
248     const int idx = itheta * nwvl + iwvl;
249     double denominator = 1.0 - rf_j[idx] * prev_rb[idx];
250     if (fabs(denominator) < DBL_EPSILON) {
251     denominator = DBL_EPSILON;
252     }
253     total_t[idx] = prev_t[idx] * t_j[idx] / denominator;
254     total_rf[idx] = prev_rf[idx] + prev_t[idx] * prev_t[idx] * rf_j[idx] / denominator;
255     total_rb[idx] = rb_j[idx] + t_j[idx] * t_j[idx] * prev_rb[idx] / denominator;
256     }
257     }
258     }
259     free(prev_t);
260     free(prev_rf);
261     free(prev_rb);
262     }
263    
264    
265     int add_layer(GlazingLayer **layers_array, int *count, int *capacity,
266     const char *filename, double thickness_m, int is_monolithic) {
267     if (*count >= *capacity) {
268     *capacity = (*capacity == 0) ? 4 : *capacity * 2;
269     GlazingLayer *new_layers = realloc(*layers_array, *capacity * sizeof(GlazingLayer));
270     if (!new_layers) {
271     perror("Failed to reallocate memory for layers");
272     return 0;
273     }
274     *layers_array = new_layers;
275     }
276    
277     memset(&((*layers_array)[*count]), 0, sizeof(GlazingLayer));
278    
279     #ifdef _WIN32
280     (*layers_array)[*count].filename = _strdup(filename);
281     #else
282     (*layers_array)[*count].filename = strdup(filename);
283     #endif
284    
285     if (!(*layers_array)[*count].filename) {
286     perror("Failed to duplicate filename");
287     return 0;
288     }
289     (*layers_array)[*count].is_mono = is_monolithic;
290     (*layers_array)[*count].thickness_m = thickness_m;
291     (*layers_array)[*count].t_0 = NULL;
292     (*layers_array)[*count].rf_0 = NULL;
293     (*layers_array)[*count].rb_0 = NULL;
294     (*layers_array)[*count].t_lambda_theta = NULL;
295     (*layers_array)[*count].rf_lambda_theta = NULL;
296     (*layers_array)[*count].rb_lambda_theta = NULL;
297    
298     (*count)++;
299     return 1;
300     }
301    
302    
303     int interpolate(GlazingLayer *layer, const double wvl_start, const double wvl_end, const double wvl_interval, const int nwvl) {
304     DATARRAY *dp = getdata(layer->filename);
305     if (!dp) {
306     fprintf(stderr, "Error: Cannot open file '%s'\n", layer->filename);
307     return 0;
308     }
309    
310     layer->t_0 = malloc(nwvl * sizeof(double));
311     layer->rf_0 = malloc(nwvl * sizeof(double));
312     layer->rb_0 = malloc(nwvl * sizeof(double));
313    
314     double wvl = wvl_start;
315     int i = 0;
316     while (wvl <= wvl_end) {
317     double t_pt[2] = {2., wvl};
318     layer->t_0[i] = datavalue(dp, t_pt);
319     double rf_pt[2] = {0., wvl};
320     layer->rf_0[i] = datavalue(dp, rf_pt);
321     double rb_pt[2] = {1., wvl};
322     layer->rb_0[i] = datavalue(dp, rb_pt);
323     wvl = wvl + wvl_interval;
324     i = i + 1;
325    
326     }
327     freedata(dp);
328     return 1;
329     }
330    
331    
332     int write_output_file(
333     const char *tfilename,
334     const char *rfilename,
335     const double *tdata,
336     const double *rfdata,
337     const double *rbdata,
338     const double wvl_start,
339     const double wvl_end,
340     const int nwvl,
341     int argc, char *argv[])
342     {
343     FILE *tfp = fopen(tfilename, "w");
344     FILE *rfp = fopen(rfilename, "w");
345     if (!tfp || !rfp) {
346     fprintf(stderr, "Error: Cannot open output files\n");
347     return(0);
348     }
349    
350     fprintf(tfp, "# ");
351     for (int i = 0; i< argc; ++i) {
352     fprintf(tfp, "%s ", argv[i]);
353     }
354     fprintf(tfp, "\n2\n0 0 %d", NTHETA);
355     for (int i = 0; i < NTHETA; i++) {
356     fprintf(tfp, " %d", (int)THETAS[i]);
357     }
358     fprintf(tfp, "\n%.0f %.0f %d\n", wvl_start, wvl_end, nwvl);
359    
360     for (int itheta = 0; itheta < NTHETA; ++itheta) {
361     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
362     fprintf(tfp, "%.6f\n", tdata[itheta*nwvl+iwvl]);
363     }
364     }
365     fclose(tfp);
366    
367     fprintf(rfp, "# ");
368     for (int i = 0; i< argc; ++i) {
369     fprintf(rfp, "%s ", argv[i]);
370     }
371     fprintf(rfp, "\n2\n0 0 %d", NTHETA * 2 - 1);
372     for (int i = 0; i < NTHETA; i++) {
373     fprintf(rfp, " %d", (int)THETAS[i]);
374     }
375     for (int i = NTHETA-2; i >= 0; --i) {
376     fprintf(rfp, " %d", (int)(180 - THETAS[i]));
377     }
378     fprintf(rfp, "\n%.0f %.0f %d\n", wvl_start, wvl_end, nwvl);
379    
380     for (int itheta = 0; itheta < NTHETA; ++itheta) {
381     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
382     fprintf(rfp, "%.6f\n", rfdata[itheta * nwvl + iwvl]);
383     }
384     }
385     for (int itheta = NTHETA - 2; itheta >= 0 ; --itheta) {
386     for (int iwvl = 0; iwvl < nwvl; ++iwvl) {
387     fprintf(rfp, "%.6f\n", rbdata[itheta * nwvl + iwvl]);
388     }
389     }
390    
391     fclose(rfp);
392     return(1);
393     }
394    
395    
396     void print_usage() {
397     fprintf(stderr, "Usage: %s [-m monolithic_layer.dat thickness | -c coated_layer.dat] ...\n", progname);
398     fprintf(stderr, "Calculate multi-layer glazing optics from spectral data files.\n");
399     fprintf(stderr, "Options:\n");
400     fprintf(stderr, " -m <filename> thickness Specify an uncoated (monolithic) glazing layer .dat file and its thickness (meter).\n");
401     fprintf(stderr, " -c <filename> Specify a coated or laminate glazing layer .dat file.\n");
402     fprintf(stderr, " -s start_wvl end_wvl interval Specify wavelength range and interval.\n");
403     fprintf(stderr, " -p prefix Specify prefix name to the output files.\n");
404     fprintf(stderr, " -h, --help Show this help message.\n");
405     fprintf(stderr, "Layer order is determined by the sequence of options on the command line.\n");
406     fprintf(stderr, "Output Files: prefix_t.dat, prefix_r.dat\n");
407     }
408    
409 greg 2.2
410     static void free_layer_data(GlazingLayer *layer) {
411     if (!layer) return;
412    
413     free(layer->filename);
414     free(layer->t_0);
415     free(layer->rf_0);
416     free(layer->rb_0);
417     free(layer->t_lambda_theta);
418     free(layer->rf_lambda_theta);
419     free(layer->rb_lambda_theta);
420    
421     // Set pointers to NULL to prevent double-free errors
422     layer->filename = NULL;
423     layer->t_0 = NULL;
424     layer->rf_0 = NULL;
425     layer->rb_0 = NULL;
426     layer->t_lambda_theta = NULL;
427     layer->rf_lambda_theta = NULL;
428     layer->rb_lambda_theta = NULL;
429     }
430    
431    
432 greg 2.1 void cleanup_layers(GlazingLayer *layers, int num_layers) {
433     if (!layers) return;
434     for (int i = 0; i < num_layers; ++i) {
435     if (layers[i].filename)
436 greg 2.2 free_layer_data(&layers[i]);
437 greg 2.1 }
438     free(layers);
439     }
440    
441    
442     int main(int argc, char *argv[])
443     {
444     GlazingLayer *layers = NULL;
445     int num_layers = 0;
446     int layer_capacity = 0;
447     int success = 1;
448     double wvl_start_nm = 380.;
449     double wvl_end_nm = 780.;
450     double wvl_interval_nm = 5.;
451     int nwvl = 81;
452     double thickness_m = 0.003;
453     char *filename;
454     char *prefix = "unnamed";
455     char file_t[MAX_NAME];
456     char file_r[MAX_NAME];
457     progname = argv[0];
458     if (argc <= 1) {
459     print_usage();
460     return EXIT_FAILURE;
461     }
462    
463     for (int i=1; i < argc; ++i) {
464     switch (argv[i][1]) {
465     case 'm':
466     filename = argv[++i];
467     thickness_m = atof(argv[++i]);
468     if (!add_layer(&layers, &num_layers, &layer_capacity, filename, thickness_m, 1)) {
469 greg 2.2 cleanup_layers(layers, num_layers);
470     return EXIT_FAILURE;
471 greg 2.1 }
472     break;
473     case 'c':
474     filename = argv[++i];
475     if (!add_layer(&layers, &num_layers, &layer_capacity, filename, thickness_m, 0)) {
476     cleanup_layers(layers, num_layers);
477     return EXIT_FAILURE;
478     }
479     break;
480     case 'p':
481     prefix = argv[++i];
482     break;
483     case 's':
484     wvl_start_nm = atof(argv[++i]);
485     wvl_end_nm = atof(argv[++i]);
486     wvl_interval_nm = atof(argv[++i]);
487     if (wvl_start_nm > wvl_end_nm) {
488     fprintf(stderr, "Starting wavelength > End wavelength\n");
489 greg 2.2 cleanup_layers(layers, num_layers);
490 greg 2.1 return EXIT_FAILURE;
491     }
492     if (((int)(wvl_end_nm - wvl_start_nm) % (int)wvl_interval_nm) > 0) {
493     fprintf(stderr,
494     "Error: Wavelength range (%f to %f nm) must be evenly divisible by the interval (%f nm).\n",
495     wvl_start_nm, wvl_end_nm, wvl_interval_nm);
496 greg 2.2 cleanup_layers(layers, num_layers);
497 greg 2.1 return EXIT_FAILURE;
498     }
499     break;
500     case 'h':
501     print_usage();
502     return 1;
503     default:
504     break;
505     }
506     }
507     snprintf(file_t, MAX_NAME, "%s_t.dat", prefix);
508     snprintf(file_r, MAX_NAME, "%s_r.dat", prefix);
509    
510     if (fabs(wvl_start_nm - wvl_end_nm) < DBL_EPSILON && wvl_interval_nm > DBL_EPSILON) {
511     nwvl = 1;
512     } else {
513     nwvl = (int)((wvl_end_nm - wvl_start_nm) / wvl_interval_nm + 1.5);
514     }
515    
516     if (num_layers == 0) {
517     fprintf(stderr, "Error: No layers specified.\n");
518     print_usage();
519     return EXIT_FAILURE;
520     }
521    
522     for (int i = 0; i < num_layers; ++i) {
523     if (!interpolate(&layers[i], wvl_start_nm, wvl_end_nm, wvl_interval_nm, nwvl)) {
524     fprintf(stderr, "Error: Failed to parse or interpolate data for layer %d.\n", i + 1);
525     success = 0;
526     break;
527     }
528     layers[i].t_lambda_theta = malloc(NTHETA * nwvl * sizeof(double));
529     layers[i].rf_lambda_theta = malloc(NTHETA * nwvl * sizeof(double));
530     layers[i].rb_lambda_theta = malloc(NTHETA * nwvl * sizeof(double));
531    
532     if (layers[i].is_mono) {
533     angular_monolithic(&layers[i], nwvl);
534     } else {
535     angular_coated(&layers[i], nwvl);
536     }
537     }
538    
539     if (!success) {
540     cleanup_layers(layers, num_layers);
541     return EXIT_FAILURE;
542     }
543    
544    
545     size_t total_flat_size = (size_t)NTHETA * nwvl * sizeof(double);
546     double *total_t = malloc(total_flat_size);
547     double *total_rf = malloc(total_flat_size);
548     double *total_rb = malloc(total_flat_size);
549    
550     if (!total_t || !total_rf || !total_rb) {
551     perror("Failed to allocate memory for final results");
552     cleanup_layers(layers, num_layers);
553     free(total_t);
554     free(total_rf);
555     free(total_rb);
556     return EXIT_FAILURE;
557     }
558    
559     multi_layer_calc(layers, num_layers, nwvl, total_t, total_rf, total_rb);
560    
561     success &= write_output_file(file_t, file_r, total_t, total_rf, total_rb, wvl_start_nm, wvl_end_nm, nwvl, argc, argv);
562    
563     if (success) {
564     printf("# ");
565     for (int i = 0;i < argc; ++i) {
566     printf("%s ", argv[i]);
567     }
568 greg 2.2 printf("\nvoid specdata refl_spec_%s\n", prefix);
569 greg 2.1 printf("4 noop %s . 'Acos(Rdot)/DEGREE'\n0\n0\n\n", file_r);
570 greg 2.2 printf("void specdata trans_spec_%s\n", prefix);
571 greg 2.1 printf("4 noop %s . 'Acos(abs(Rdot))/DEGREE'\n0\n0\n\n", file_t);
572 greg 2.2 printf("void WGMDfunc glaze_mat_%s\n13\n\trefl_spec_%s 1 0 0\n\ttrans_spec_%s 1 0 0\n\tvoid\n\t0 0 1 .\n0\n9 0 0 0 0 0 0 0 0 0\n", prefix, prefix, prefix);
573 greg 2.1 } else {
574     fprintf(stderr, "Error: Failed to write one or more output files.\n");
575     }
576    
577     cleanup_layers(layers, num_layers);
578     free(total_t);
579     free(total_rf);
580     free(total_rb);
581    
582     return success ? EXIT_SUCCESS : EXIT_FAILURE;
583     }