1 /* -*- Mode: C; tab-width: 8; indent-tabs-mode: t; c-basic-offset: 8 -*-
2  *
3  * Copyright (C) 2014-2015 Richard Hughes <richard@hughsie.com>
4  *
5  * Licensed under the GNU Lesser General Public License Version 2.1
6  *
7  * This library is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * This library is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with this library; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301 USA
20  */
21 
22 /**
23  * SECTION:cd-spectrum
24  * @short_description: A single set of spectral values
25  *
26  * Functions to manipulate spectral values.
27  */
28 
29 #include "config.h"
30 
31 #include <math.h>
32 #include <glib-object.h>
33 
34 #include "cd-color.h"
35 #include "cd-interp-linear.h"
36 #include "cd-spectrum.h"
37 
38 /* this is private */
39 struct _CdSpectrum {
40 	guint			 reserved_size;
41 	gchar			*id;
42 	gdouble			 start;
43 	gdouble			 end;
44 	gdouble			 norm;
45 	gdouble			 wavelength_cal[3];
46 	GArray			*data;
47 };
48 
49 /**
50  * cd_spectrum_dup:
51  * @spectrum: a #CdSpectrum instance.
52  *
53  * Since: 1.1.6
54  **/
55 CdSpectrum *
cd_spectrum_dup(const CdSpectrum * spectrum)56 cd_spectrum_dup (const CdSpectrum *spectrum)
57 {
58 	CdSpectrum *dest;
59 	gdouble tmp;
60 	guint i;
61 
62 	g_return_val_if_fail (spectrum != NULL, NULL);
63 
64 	dest = cd_spectrum_new ();
65 	dest->id = g_strdup (spectrum->id);
66 	dest->start = spectrum->start;
67 	dest->end = spectrum->end;
68 	dest->norm = spectrum->norm;
69 	for (i = 0; i < spectrum->data->len; i++) {
70 		tmp = cd_spectrum_get_value_raw (spectrum, i);
71 		cd_spectrum_add_value (dest, tmp);
72 	}
73 	for (i = 0; i < 3; i++)
74 		dest->wavelength_cal[i] = spectrum->wavelength_cal[i];
75 	return dest;
76 }
77 
78 /**
79  * cd_spectrum_get_id:
80  * @spectrum: a #CdSpectrum instance.
81  *
82  * Gets the spectral data.
83  *
84  * Return value: the textual ID of the sample
85  *
86  * Since: 1.1.6
87  **/
88 const gchar *
cd_spectrum_get_id(const CdSpectrum * spectrum)89 cd_spectrum_get_id (const CdSpectrum *spectrum)
90 {
91 	g_return_val_if_fail (spectrum != NULL, NULL);
92 	return spectrum->id;
93 }
94 
95 /**
96  * cd_spectrum_get_value:
97  * @spectrum: a #CdSpectrum instance.
98  * @idx: an index into the data
99  *
100  * Gets the spectrum data at a specified index.
101  *
102  * Return value: spectral data value, or -1 for invalid
103  *
104  * Since: 1.1.6
105  **/
106 gdouble
cd_spectrum_get_value(const CdSpectrum * spectrum,guint idx)107 cd_spectrum_get_value (const CdSpectrum *spectrum, guint idx)
108 {
109 	g_return_val_if_fail (spectrum != NULL, -1.0f);
110 	g_return_val_if_fail (idx < spectrum->data->len, -1.0f);
111 	return g_array_index (spectrum->data, gdouble, idx) * spectrum->norm;
112 }
113 
114 /**
115  * cd_spectrum_get_value_max:
116  * @spectrum: a #CdSpectrum instance.
117  *
118  * Gets the largest normalised value in the spectrum.
119  *
120  * Since: 1.3.1
121  **/
122 gdouble
cd_spectrum_get_value_max(const CdSpectrum * spectrum)123 cd_spectrum_get_value_max (const CdSpectrum *spectrum)
124 {
125 	gdouble max = 0.f;
126 	guint i;
127 	for (i = 0; i < cd_spectrum_get_size (spectrum); i++)
128 		max = MAX (max, cd_spectrum_get_value (spectrum, i));
129 	return max;
130 }
131 
132 /**
133  * cd_spectrum_get_value_min:
134  * @spectrum: a #CdSpectrum instance.
135  *
136  * Gets the smallest normalised value in the spectrum.
137  *
138  * Since: 1.3.1
139  **/
140 gdouble
cd_spectrum_get_value_min(const CdSpectrum * spectrum)141 cd_spectrum_get_value_min (const CdSpectrum *spectrum)
142 {
143 	gdouble min = G_MAXDOUBLE;
144 	guint i;
145 	for (i = 0; i < cd_spectrum_get_size (spectrum); i++)
146 		min = MIN (min, cd_spectrum_get_value (spectrum, i));
147 	return min;
148 }
149 
150 /**
151  * cd_spectrum_set_value:
152  * @spectrum: a #CdSpectrum instance.
153  * @idx: an index into the data
154  * @data: a data value
155  *
156  * Overwrites the spectrum data at a specified index.
157  *
158  * Since: 1.2.6
159  **/
160 void
cd_spectrum_set_value(CdSpectrum * spectrum,guint idx,gdouble data)161 cd_spectrum_set_value (CdSpectrum *spectrum, guint idx, gdouble data)
162 {
163 	g_return_if_fail (spectrum != NULL);
164 	g_return_if_fail (idx < spectrum->data->len);
165 	g_array_index (spectrum->data, gdouble, idx) = data / spectrum->norm;
166 }
167 
168 /**
169  * cd_spectrum_get_value_raw:
170  * @spectrum: a #CdSpectrum instance.
171  * @idx: an index into the data
172  *
173  * Gets the spectrum data at a specified index, without any normalization
174  * applied. Most people should use cd_spectrum_get_value() instead.
175  *
176  * Return value: spectral data value, or -1 for invalid
177  *
178  * Since: 1.2.6
179  **/
180 gdouble
cd_spectrum_get_value_raw(const CdSpectrum * spectrum,guint idx)181 cd_spectrum_get_value_raw (const CdSpectrum *spectrum, guint idx)
182 {
183 	g_return_val_if_fail (spectrum != NULL, -1.0f);
184 	g_return_val_if_fail (idx < spectrum->data->len, -1.0f);
185 	return g_array_index (spectrum->data, gdouble, idx);
186 }
187 
188 /**
189  * cd_spectrum_get_wavelength:
190  * @spectrum: a #CdSpectrum instance.
191  * @idx: an index into the data
192  *
193  * Gets the wavelenth that corresponds to the specified index.
194  *
195  * Return value: wavelenth value in nm, or -1 for invalid
196  *
197  * Since: 1.1.6
198  **/
199 gdouble
cd_spectrum_get_wavelength(const CdSpectrum * spectrum,guint idx)200 cd_spectrum_get_wavelength (const CdSpectrum *spectrum, guint idx)
201 {
202 	g_return_val_if_fail (spectrum != NULL, -1.0f);
203 
204 	/* fall back to the old method */
205 	if (spectrum->wavelength_cal[0] < 0) {
206 		gdouble step;
207 		guint number_points;
208 		/* if we used cd_spectrum_size_new() and there is no data we can infer
209 		 * the wavelenth based on the declared initial size */
210 		if (spectrum->reserved_size > 0)
211 			number_points = spectrum->reserved_size;
212 		else
213 			number_points = spectrum->data->len;
214 		step = (spectrum->end - spectrum->start) / (number_points - 1);
215 		return spectrum->start + (step * (gdouble) idx);
216 	}
217 
218 	/* use wavelength_cal to work out wavelength */
219 	return spectrum->start +
220 		spectrum->wavelength_cal[0] * (gdouble) idx +
221 		spectrum->wavelength_cal[1] * pow (idx, 2) +
222 		spectrum->wavelength_cal[2] * pow (idx, 3);
223 }
224 
225 /**
226  * cd_spectrum_get_size:
227  * @spectrum: a #CdSpectrum instance.
228  *
229  * Gets the size of the spectrum data.
230  *
231  * Return value: number of data items in this spectrum
232  *
233  * Since: 1.1.6
234  **/
235 guint
cd_spectrum_get_size(const CdSpectrum * spectrum)236 cd_spectrum_get_size (const CdSpectrum *spectrum)
237 {
238 	g_return_val_if_fail (spectrum != NULL, G_MAXUINT);
239 	return spectrum->data->len;
240 }
241 
242 /**
243  * cd_spectrum_get_data:
244  * @spectrum: a #CdSpectrum instance.
245  *
246  * Gets the spectral data.
247  * NOTE: This is not normalized
248  *
249  * Return value: (transfer none) (element-type gdouble): spectral data
250  *
251  * Since: 1.1.6
252  **/
253 GArray *
cd_spectrum_get_data(const CdSpectrum * spectrum)254 cd_spectrum_get_data (const CdSpectrum *spectrum)
255 {
256 	g_return_val_if_fail (spectrum != NULL, NULL);
257 	return spectrum->data;
258 }
259 
260 /**
261  * cd_spectrum_get_start:
262  * @spectrum: a #CdSpectrum instance.
263  *
264  * Gets the start value of the spectral data.
265  *
266  * Return value: the value in nm
267  *
268  * Since: 1.1.6
269  **/
270 gdouble
cd_spectrum_get_start(const CdSpectrum * spectrum)271 cd_spectrum_get_start (const CdSpectrum *spectrum)
272 {
273 	g_return_val_if_fail (spectrum != NULL, 0.0f);
274 	return spectrum->start;
275 }
276 
277 /**
278  * cd_spectrum_get_end:
279  * @spectrum: a #CdSpectrum instance.
280  *
281  * Gets the end value of the spectral data.
282  *
283  * Return value: the value in nm
284  *
285  * Since: 1.1.6
286  **/
287 gdouble
cd_spectrum_get_end(const CdSpectrum * spectrum)288 cd_spectrum_get_end (const CdSpectrum *spectrum)
289 {
290 	g_return_val_if_fail (spectrum != NULL, 0.0f);
291 	return spectrum->end;
292 }
293 
294 /**
295  * cd_spectrum_get_norm:
296  * @spectrum: a #CdSpectrum instance.
297  *
298  * Gets the normalization value of the spectral data.
299  * NOTE: This affects every value in the spectrum.
300  *
301  * Return value: the value
302  *
303  * Since: 1.1.6
304  **/
305 gdouble
cd_spectrum_get_norm(const CdSpectrum * spectrum)306 cd_spectrum_get_norm (const CdSpectrum *spectrum)
307 {
308 	g_return_val_if_fail (spectrum != NULL, 0.0f);
309 	return spectrum->norm;
310 }
311 
312 /**
313  * cd_spectrum_get_resolution:
314  * @spectrum: a #CdSpectrum instance.
315  *
316  * Gets the divisor of the spectra, for instance a .
317  *
318  * Return value: the value
319  *
320  * Since: 1.2.6
321  **/
322 gdouble
cd_spectrum_get_resolution(const CdSpectrum * spectrum)323 cd_spectrum_get_resolution (const CdSpectrum *spectrum)
324 {
325 	g_return_val_if_fail (spectrum != NULL, 0.0f);
326 	return (spectrum->end - spectrum->start) / (gdouble) spectrum->data->len;
327 }
328 
329 /**
330  * cd_spectrum_get_type:
331  *
332  * Gets a specific type.
333  *
334  * Return value: a #GType
335  *
336  * Since: 1.1.6
337  **/
338 GType
cd_spectrum_get_type(void)339 cd_spectrum_get_type (void)
340 {
341 	static GType type_id = 0;
342 	if (!type_id)
343 		type_id = g_boxed_type_register_static ("CdSpectrum",
344 							(GBoxedCopyFunc) cd_spectrum_dup,
345 							(GBoxedFreeFunc) cd_spectrum_free);
346 	return type_id;
347 }
348 
349 /**
350  * cd_spectrum_new:
351  *
352  * Allocates a spectrum.
353  *
354  * Return value: A newly allocated #CdSpectrum object
355  *
356  * Since: 1.1.6
357  **/
358 CdSpectrum *
cd_spectrum_new(void)359 cd_spectrum_new (void)
360 {
361 	CdSpectrum *spectrum;
362 	spectrum = g_slice_new0 (CdSpectrum);
363 	spectrum->norm = 1.f;
364 	spectrum->data = g_array_new (FALSE, FALSE, sizeof (gdouble));
365 	spectrum->wavelength_cal[0] = -1.f;
366 	return spectrum;
367 }
368 
369 /**
370  * cd_spectrum_sized_new:
371  * @reserved_size: the future size of the spectrum
372  *
373  * Allocates a spectrum with a preallocated size.
374  *
375  * Return value: A newly allocated #CdSpectrum object
376  *
377  * Since: 1.1.6
378  **/
379 CdSpectrum *
cd_spectrum_sized_new(guint reserved_size)380 cd_spectrum_sized_new (guint reserved_size)
381 {
382 	CdSpectrum *spectrum;
383 	spectrum = g_slice_new0 (CdSpectrum);
384 	spectrum->norm = 1.f;
385 	spectrum->reserved_size = reserved_size;
386 	spectrum->data = g_array_sized_new (FALSE, FALSE, sizeof (gdouble), reserved_size);
387 	spectrum->wavelength_cal[0] = -1.f;
388 	return spectrum;
389 }
390 
391 /**
392  * cd_spectrum_planckian_new_full:
393  * @temperature: the temperature in Kelvin
394  * @start: the new spectrum start
395  * @end: the new spectrum end
396  * @resolution: the resolution to use when resampling
397  *
398  * Allocates a Planckian spectrum at a specific temperature.
399  *
400  * Return value: A newly allocated #CdSpectrum object
401  *
402  * Since: 1.3.1
403  **/
404 CdSpectrum *
cd_spectrum_planckian_new_full(gdouble temperature,gdouble start,gdouble end,gdouble resolution)405 cd_spectrum_planckian_new_full (gdouble temperature,
406 				gdouble start,
407 				gdouble end,
408 				gdouble resolution)
409 {
410 	CdSpectrum *s = NULL;
411 	const gdouble c1 = 3.74183e-16;	/* 2pi * h * c^2 */
412 	const gdouble c2 = 1.4388e-2;	/* h * c / k */
413 	gdouble wl;
414 	gdouble norm;
415 	gdouble tmp;
416 	guint i;
417 
418 	/* sanity check */
419 	if (temperature < 1.0 || temperature > 1e6)
420 		return NULL;
421 
422 	/* create spectrum with 1nm resolution */
423 	s = cd_spectrum_sized_new (531);
424 	s->id = g_strdup_printf ("Planckian@%.0fK", temperature);
425 	cd_spectrum_set_start (s, start);
426 	cd_spectrum_set_end (s, end);
427 
428 	/* see http://www.create.uwe.ac.uk/ardtalks/Schanda_paper.pdf, page 42 */
429 	wl = 560 * 1e-9;
430 	norm = 0.01 * (c1 * pow (wl, -5.0)) / (exp (c2 / (wl * temperature)) - 1.0);
431 	for (i = 0; i < s->reserved_size; i++) {
432 		wl = cd_spectrum_get_wavelength (s, i) * 1e-9;
433 		tmp = (c1 * pow (wl, -5.0)) / (exp (c2 / (wl * temperature)) - 1.0);
434 		cd_spectrum_add_value (s, tmp / norm);
435 	}
436 	return s;
437 }
438 
439 /**
440  * cd_spectrum_planckian_new:
441  * @temperature: the temperature in Kelvin
442  *
443  * Allocates a Planckian spectrum at a specific temperature.
444  *
445  * Return value: A newly allocated #CdSpectrum object
446  *
447  * Since: 1.1.6
448  **/
449 CdSpectrum *
cd_spectrum_planckian_new(gdouble temperature)450 cd_spectrum_planckian_new (gdouble temperature)
451 {
452 	return cd_spectrum_planckian_new_full (temperature, 300, 830, 1);
453 }
454 
455 /**
456  * cd_spectrum_add_value:
457  * @spectrum: the spectrum
458  *
459  * Adds a value in nm to the spectrum.
460  *
461  * Since: 1.1.6
462  **/
463 void
cd_spectrum_add_value(CdSpectrum * spectrum,gdouble data)464 cd_spectrum_add_value (CdSpectrum *spectrum, gdouble data)
465 {
466 	g_return_if_fail (spectrum != NULL);
467 	g_array_append_val (spectrum->data, data);
468 }
469 
470 /**
471  * cd_spectrum_free:
472  * @spectrum: the spectrum
473  *
474  * Deallocates a color spectrum.
475  *
476  * Since: 1.1.6
477  **/
478 void
cd_spectrum_free(CdSpectrum * spectrum)479 cd_spectrum_free (CdSpectrum *spectrum)
480 {
481 	if (spectrum == NULL)
482 		return;
483 	g_free (spectrum->id);
484 	g_array_unref (spectrum->data);
485 	g_slice_free (CdSpectrum, spectrum);
486 }
487 
488 /**
489  * cd_spectrum_set_id:
490  * @spectrum: the destination spectrum
491  * @id: component id
492  *
493  * Sets a spectrum id.
494  *
495  * Since: 1.1.6
496  **/
497 void
cd_spectrum_set_id(CdSpectrum * spectrum,const gchar * id)498 cd_spectrum_set_id (CdSpectrum *spectrum, const gchar *id)
499 {
500 	g_return_if_fail (spectrum != NULL);
501 	g_return_if_fail (id != NULL);
502 	g_free (spectrum->id);
503 	spectrum->id = g_strdup (id);
504 }
505 
506 /**
507  * cd_spectrum_set_data:
508  * @spectrum: the destination spectrum
509  * @value: (element-type gdouble): component value
510  *
511  * Sets the spectrum data.
512  *
513  * Since: 1.1.6
514  **/
515 void
cd_spectrum_set_data(CdSpectrum * spectrum,GArray * value)516 cd_spectrum_set_data (CdSpectrum *spectrum, GArray *value)
517 {
518 	g_return_if_fail (spectrum != NULL);
519 	g_return_if_fail (value != NULL);
520 	g_array_unref (spectrum->data);
521 	spectrum->data = g_array_ref (value);
522 }
523 
524 /**
525  * cd_spectrum_set_start:
526  * @spectrum: a #CdSpectrum instance.
527  * @start: the start value of the spectral data
528  *
529  * Set the start value of the spectal data in nm.
530  *
531  * Since: 1.1.6
532  **/
533 void
cd_spectrum_set_start(CdSpectrum * spectrum,gdouble start)534 cd_spectrum_set_start (CdSpectrum *spectrum, gdouble start)
535 {
536 	g_return_if_fail (spectrum != NULL);
537 	spectrum->start = start;
538 }
539 
540 /**
541  * cd_spectrum_set_end:
542  * @spectrum: a #CdSpectrum instance.
543  * @end: the end value of the spectral data
544  *
545  * Set the end value of the spectal data in nm.
546  *
547  * If there is already spectral data, the wavelength calibration will
548  * also be set automatically.
549  *
550  * Since: 1.1.6
551  **/
552 void
cd_spectrum_set_end(CdSpectrum * spectrum,gdouble end)553 cd_spectrum_set_end (CdSpectrum *spectrum, gdouble end)
554 {
555 	g_return_if_fail (spectrum != NULL);
556 
557 	/* calculate the calibration co-efficients */
558 	if (spectrum->data->len > 1) {
559 		spectrum->wavelength_cal[0] = (end - spectrum->start) /
560 						(spectrum->data->len - 1);
561 		spectrum->wavelength_cal[1] = 0.f;
562 		spectrum->wavelength_cal[2] = 0.f;
563 	}
564 
565 	/* set this for later */
566 	spectrum->end = end;
567 }
568 
569 /**
570  * cd_spectrum_set_norm:
571  * @spectrum: a #CdSpectrum instance.
572  * @norm: the end value of the spectral data
573  *
574  * Set the normalization value of the spectrum.
575  * NOTE: This affects every value in the spectrum.
576  *
577  * Since: 1.1.6
578  **/
579 void
cd_spectrum_set_norm(CdSpectrum * spectrum,gdouble norm)580 cd_spectrum_set_norm (CdSpectrum *spectrum, gdouble norm)
581 {
582 	g_return_if_fail (spectrum != NULL);
583 	spectrum->norm = norm;
584 }
585 
586 /**
587  * cd_spectrum_get_value_for_nm:
588  * @spectrum: a #CdSpectrum instance.
589  * @wavelength: the wavelength in nm
590  *
591  * Gets the value from the spectral data for a given wavelength.
592  *
593  * Return value: the value for the wavelength
594  *
595  * Since: 1.1.6
596  **/
597 gdouble
cd_spectrum_get_value_for_nm(const CdSpectrum * spectrum,gdouble wavelength)598 cd_spectrum_get_value_for_nm (const CdSpectrum *spectrum, gdouble wavelength)
599 {
600 	guint i;
601 	guint size;
602 	g_autoptr(CdInterp) interp = NULL;
603 
604 	g_return_val_if_fail (spectrum != NULL, -1.f);
605 
606 	/* out of bounds */
607 	size = cd_spectrum_get_size (spectrum);
608 	if (size == 0)
609 		return 1.f;
610 	if (wavelength < spectrum->start)
611 		return cd_spectrum_get_value (spectrum, 0);
612 	if (wavelength > spectrum->end)
613 		return cd_spectrum_get_value (spectrum, size - 1);
614 
615 	/* add all the data points */
616 	interp = cd_interp_linear_new ();
617 	for (i = 0; i < size; i++) {
618 		cd_interp_insert (interp,
619 				  cd_spectrum_get_wavelength (spectrum, i),
620 				  cd_spectrum_get_value (spectrum, i));
621 	}
622 
623 	/* get the interpolated value */
624 	if (!cd_interp_prepare (interp, NULL))
625 		return -1.f;
626 	return cd_interp_eval (interp, wavelength, NULL);
627 }
628 
629 /**
630  * cd_spectrum_limit_min:
631  * @spectrum: a #CdSpectrum instance
632  * @value: the threshold value to limit the spectrum
633  *
634  * Ensures no values in the spectrum fall below a set limit. If they
635  * are found, set them to @value.
636  *
637  * Since: 1.3.1
638  **/
639 void
cd_spectrum_limit_min(CdSpectrum * spectrum,gdouble value)640 cd_spectrum_limit_min (CdSpectrum *spectrum, gdouble value)
641 {
642 	gdouble tmp;
643 	guint i;
644 	for (i = 0; i < spectrum->data->len; i++) {
645 		tmp = cd_spectrum_get_value (spectrum, i);
646 		if (tmp < value)
647 			cd_spectrum_set_value (spectrum, i, value);
648 	}
649 }
650 
651 /**
652  * cd_spectrum_limit_max:
653  * @spectrum: a #CdSpectrum instance
654  * @value: the threshold value to limit the spectrum
655  *
656  * Ensures no values in the spectrum fall above a set limit. If they
657  * are found, set them to @value.
658  *
659  * Since: 1.3.1
660  **/
661 void
cd_spectrum_limit_max(CdSpectrum * spectrum,gdouble value)662 cd_spectrum_limit_max (CdSpectrum *spectrum, gdouble value)
663 {
664 	gdouble tmp;
665 	guint i;
666 	for (i = 0; i < spectrum->data->len; i++) {
667 		tmp = cd_spectrum_get_value (spectrum, i);
668 		if (tmp > value)
669 			cd_spectrum_set_value (spectrum, i, value);
670 	}
671 }
672 
673 /**
674  * cd_spectrum_normalize:
675  * @spectrum: a #CdSpectrum instance
676  * @wavelength: the wavelength in nm
677  * @value: the value to normalize to
678  *
679  * Normalizes a spectrum to a specific value at a specific wavelength.
680  *
681  * Since: 1.1.6
682  **/
683 void
cd_spectrum_normalize(CdSpectrum * spectrum,gdouble wavelength,gdouble value)684 cd_spectrum_normalize (CdSpectrum *spectrum, gdouble wavelength, gdouble value)
685 {
686 	gdouble tmp;
687 	tmp = cd_spectrum_get_value_for_nm (spectrum, wavelength);
688 	spectrum->norm *= value / tmp;
689 }
690 
691 /**
692  * cd_spectrum_normalize_max:
693  * @spectrum: a #CdSpectrum instance
694  * @value: the value to normalize to
695  *
696  * Normalizes a spectrum to a specific value at its maximum value.
697  *
698  * Since: 1.2.6
699  **/
700 void
cd_spectrum_normalize_max(CdSpectrum * spectrum,gdouble value)701 cd_spectrum_normalize_max (CdSpectrum *spectrum, gdouble value)
702 {
703 	gdouble max = 0.f;
704 	gdouble tmp;
705 	guint i;
706 
707 	for (i = 0; i < spectrum->data->len; i++) {
708 		tmp = cd_spectrum_get_value_raw (spectrum, i);
709 		if (tmp > max)
710 			max = tmp;
711 	}
712 	if (max > 0.f)
713 		spectrum->norm = value / max;
714 }
715 
716 /**
717  * cd_spectrum_multiply:
718  * @s1: a #CdSpectrum instance, possibly an illuminant.
719  * @s2: a #CdSpectrum instance, possibly an absorption spectrum.
720  * @resolution: the step size in nm
721  *
722  * Multiplies two spectra together.
723  *
724  * Return value: a #CdSpectrum instance
725  *
726  * Since: 1.1.6
727  **/
728 CdSpectrum *
cd_spectrum_multiply(CdSpectrum * s1,CdSpectrum * s2,gdouble resolution)729 cd_spectrum_multiply (CdSpectrum *s1, CdSpectrum *s2, gdouble resolution)
730 {
731 	CdSpectrum *s;
732 	gdouble i;
733 
734 	s = cd_spectrum_new ();
735 	s->id = g_strdup_printf ("%s✕%s", s1->id, s2->id);
736 	s->start = MAX (s1->start, s2->start);
737 	s->end = MIN (s1->end, s2->end);
738 	for (i = s->start; i <= s->end; i += resolution) {
739 		cd_spectrum_add_value (s, cd_spectrum_get_value_for_nm (s1, i) *
740 					  cd_spectrum_get_value_for_nm (s2, i));
741 	}
742 	return s;
743 }
744 
745 /**
746  * cd_spectrum_multiply_scalar:
747  * @spectrum: a #CdSpectrum instance
748  * @value: a scalar value
749  *
750  * Multiplies a spectra with a scalar value.
751  *
752  * Return value: a #CdSpectrum instance
753  *
754  * Since: 1.3.5
755  **/
756 CdSpectrum *
cd_spectrum_multiply_scalar(CdSpectrum * spectrum,gdouble value)757 cd_spectrum_multiply_scalar (CdSpectrum *spectrum, gdouble value)
758 {
759 	CdSpectrum *s = cd_spectrum_dup (spectrum);
760 	for (guint i = 0; i < spectrum->data->len; i++)
761 		cd_spectrum_add_value (s, cd_spectrum_get_value (spectrum, i) * value);
762 	return s;
763 }
764 
765 /**
766  * cd_spectrum_subtract:
767  * @s1: a #CdSpectrum instance, e.g. a sample
768  * @s2: a #CdSpectrum instance, e.g. a dark calibration
769  * @resolution: the resolution to use when resampling
770  *
771  * Subtracts one spectral plot from another. If the spectra have the same start,
772  * end and the same number of data points they are not resampled.
773  *
774  * Return value: a #CdSpectrum instance
775  *
776  * Since: 1.3.1
777  **/
778 CdSpectrum *
cd_spectrum_subtract(CdSpectrum * s1,CdSpectrum * s2,gdouble resolution)779 cd_spectrum_subtract (CdSpectrum *s1, CdSpectrum *s2, gdouble resolution)
780 {
781 	CdSpectrum *s;
782 	gdouble max;
783 	gdouble min;
784 	gdouble nm;
785 	guint i;
786 
787 	g_return_val_if_fail (s1 != NULL, NULL);
788 	g_return_val_if_fail (s2 != NULL, NULL);
789 
790 	/* we can do this without resampling */
791 	if (fabs (s1->start - s2->start) < 0.01f &&
792 	    fabs (s1->end - s2->end) < 0.01f &&
793 	    s1->data->len == s2->data->len) {
794 		s = cd_spectrum_sized_new (s1->data->len);
795 		s->id = g_strdup_printf ("%s-%s", s1->id, s2->id);
796 		s->start = s1->start;
797 		s->end = s1->end;
798 		for (i = 0; i < 3; i++)
799 			s->wavelength_cal[i] = s1->wavelength_cal[i];
800 		for (i = 0; i < s1->data->len; i++) {
801 			gdouble tmp;
802 			tmp = cd_spectrum_get_value (s1, i) - cd_spectrum_get_value (s2, i);
803 			cd_spectrum_add_value (s, tmp);
804 		}
805 		return s;
806 	}
807 
808 	/* resample */
809 	min = MIN (cd_spectrum_get_start (s1), cd_spectrum_get_start (s2));
810 	max = MAX (cd_spectrum_get_end (s1), cd_spectrum_get_end (s2));
811 	s = cd_spectrum_new ();
812 	s->id = g_strdup_printf ("%s-%s", s1->id, s2->id);
813 	s->start = min;
814 	s->end = max;
815 	for (nm = min; nm <= max; nm += resolution) {
816 		gdouble tmp;
817 		tmp = cd_spectrum_get_value_for_nm (s1, nm) -
818 			cd_spectrum_get_value_for_nm (s2, nm);
819 		cd_spectrum_add_value (s, tmp);
820 	}
821 	return s;
822 }
823 
824 /**
825  * cd_spectrum_to_string:
826  * @spectrum: a #CdSpectrum instance
827  * @max_width: the terminal width
828  * @max_height: the terminal height
829  *
830  * Returns a graphical representation of the spectrum.
831  *
832  * Return value: a printable ASCII string
833  *
834  * Since: 1.3.1
835  **/
836 gchar *
cd_spectrum_to_string(CdSpectrum * spectrum,guint max_width,guint max_height)837 cd_spectrum_to_string (CdSpectrum *spectrum, guint max_width, guint max_height)
838 {
839 	GString *str = g_string_new ("");
840 	guint i, j;
841 	gdouble val_max;
842 	gdouble nm_scale;
843 
844 	/* make space for the axes */
845 	max_width -= 9;
846 	max_height -= 2;
847 
848 	/* find value maximum */
849 	val_max = cd_spectrum_get_value_max (spectrum);
850 	if (val_max < 0.001)
851 		val_max = 0.001;
852 	nm_scale = (cd_spectrum_get_end (spectrum) -
853 		    cd_spectrum_get_start (spectrum)) / (gdouble) (max_width - 1);
854 
855 	/* draw grid */
856 	for (i = 0; i < max_height; i++) {
857 		gdouble val;
858 		val = val_max / (gdouble) max_height * (max_height - i);
859 		g_string_append_printf (str, "%7.3f |", val);
860 		for (j = 0; j < max_width; j++) {
861 			gdouble nm;
862 			nm = ((gdouble) j * nm_scale) + cd_spectrum_get_start (spectrum);
863 			if (cd_spectrum_get_value_for_nm (spectrum, nm) >= val)
864 				g_string_append (str, "#");
865 			else
866 				g_string_append (str, "_");
867 		}
868 		g_string_append (str, "\n");
869 	}
870 
871 	/* draw x axis */
872 	g_string_append_printf (str, "%7.3f  ", 0.f);
873 	for (j = 0; j < max_width; j++)
874 		g_string_append (str, "-");
875 	g_string_append (str, "\n");
876 
877 	/* draw X labels */
878 	g_string_append_printf (str, "         %.0fnm",
879 				cd_spectrum_get_start (spectrum));
880 	for (j = 0; j < max_width - 10; j++)
881 		g_string_append (str, " ");
882 	g_string_append_printf (str, "%.0fnm",
883 				cd_spectrum_get_end (spectrum));
884 	g_string_append (str, "\n");
885 
886 	/* success */
887 	return g_string_free (str, FALSE);
888 }
889 
890 /**
891  * cd_spectrum_set_wavelength_cal:
892  * @spectrum: a #CdSpectrum instance
893  * @c1: the 1st coefficient
894  * @c2: the 2nd coefficient
895  * @c3: the 3rd coefficient
896  *
897  * Sets the calibration coefficients used to map pixel indexes to
898  * wavelengths.
899  *
900  * This function will set the 'end' wavelength automatically,
901  * potentially overwriting the value set by cd_spectrum_set_end().
902  *
903  * Since: 1.3.1
904  **/
905 void
cd_spectrum_set_wavelength_cal(CdSpectrum * spectrum,gdouble c1,gdouble c2,gdouble c3)906 cd_spectrum_set_wavelength_cal (CdSpectrum *spectrum,
907 				gdouble c1, gdouble c2, gdouble c3)
908 {
909 	spectrum->wavelength_cal[0] = c1;
910 	spectrum->wavelength_cal[1] = c2;
911 	spectrum->wavelength_cal[2] = c3;
912 
913 	/* recalculate the end wavelength */
914 	spectrum->end = cd_spectrum_get_wavelength (spectrum,
915 						    cd_spectrum_get_size (spectrum) - 1);
916 }
917 
918 /**
919  * cd_spectrum_get_wavelength_cal:
920  * @spectrum: a #CdSpectrum instance
921  * @c1: the 1st coefficient
922  * @c2: the 2nd coefficient
923  * @c3: the 3rd coefficient
924  *
925  * Gets the calibration coefficients used to map pixel indexes to
926  * wavelengths.
927  *
928  * Since: 1.3.1
929  **/
930 void
cd_spectrum_get_wavelength_cal(CdSpectrum * spectrum,gdouble * c1,gdouble * c2,gdouble * c3)931 cd_spectrum_get_wavelength_cal (CdSpectrum *spectrum,
932 				gdouble *c1, gdouble *c2, gdouble *c3)
933 {
934 	if (c1 != NULL)
935 		*c1 = spectrum->wavelength_cal[0];
936 	if (c2 != NULL)
937 		*c2 = spectrum->wavelength_cal[1];
938 	if (c3 != NULL)
939 		*c3 = spectrum->wavelength_cal[2];
940 }
941 
942 /**
943  * cd_spectrum_resample:
944  * @spectrum: a #CdSpectrum instance
945  * @start: the new spectrum start
946  * @end: the new spectrum end
947  * @resolution: the resolution to use when resampling
948  *
949  * Resample a new spectrum with linear index to wavelength coefficients.
950  *
951  * Return value: a #CdSpectrum instance
952  *
953  * Since: 1.3.1
954  **/
955 CdSpectrum *
cd_spectrum_resample(CdSpectrum * spectrum,gdouble start,gdouble end,gdouble resolution)956 cd_spectrum_resample (CdSpectrum *spectrum,
957 		      gdouble start,
958 		      gdouble end,
959 		      gdouble resolution)
960 {
961 	gdouble nm;
962 	CdSpectrum *sp;
963 
964 	sp = cd_spectrum_new ();
965 	cd_spectrum_set_start (sp, start);
966 	for (nm = start; nm <= end; nm += resolution) {
967 		gdouble tmp;
968 		tmp = cd_spectrum_get_value_for_nm (spectrum, nm);
969 		cd_spectrum_add_value (sp, tmp);
970 	}
971 	cd_spectrum_set_end (sp, end);
972 	return sp;
973 }
974 
975 /**
976  * cd_spectrum_resample_to_size:
977  * @spectrum: a #CdSpectrum instance
978  * @size: the output spectrum size
979  *
980  * Resample a new spectrum with the desired number of points.
981  *
982  * Return value: a #CdSpectrum instance
983  *
984  * Since: 1.3.4
985  **/
986 CdSpectrum *
cd_spectrum_resample_to_size(CdSpectrum * spectrum,guint size)987 cd_spectrum_resample_to_size (CdSpectrum *spectrum, guint size)
988 {
989 	gdouble inc;
990 	guint i;
991 	CdSpectrum *sp;
992 
993 	sp = cd_spectrum_new ();
994 	cd_spectrum_set_start (sp, spectrum->start);
995 	cd_spectrum_set_end (sp, spectrum->end);
996 
997 	inc = (spectrum->end - spectrum->start) / (gdouble) (size - 1);
998 	for (i = 0; i < size; i++) {
999 		gdouble nm = spectrum->start + ((gdouble) i * inc);
1000 		gdouble tmp = cd_spectrum_get_value_for_nm (spectrum, nm);
1001 		cd_spectrum_add_value (sp, tmp);
1002 	}
1003 	return sp;
1004 }
1005