System and method for illumination source identification

Abstract

A light detection module has N optical channels, each with an optical filter, a detector, and an amplifier; and an N1 switch with N input ports each connected to one corresponding output port of each channel to receive an amplified detector output corresponding to a filtered optical intensity incident on that detector. The switch cycles between channels, connecting each amplified detector output in turn to the output port. An ADC samples a time dependent optical intensity signal from the switch, generating a corresponding ADC digital signal output. A microcontroller, connected to the N1 switch and the ADC, controls acquisition by the ADC to provide a digital voltage data stream from each channel; making the average optical intensity value characterizing the voltage data stream available from each channel at a digital output port of the microcontroller, as N data values, characterizing the light incident on the N channels of the module.

Claims

1. A light detection module receptive to incident light; the light detection module comprising: N optical channels, wherein each optical channel comprises an optical filter, a detector, and an amplifier, and each optical channel having an output port; an N1 switch having N input ports and a single output port, wherein each input port is connected to a corresponding one and only one output port of each of said N optical channels to receive a corresponding amplified detector output corresponding to a filtered optical intensity incident on the detector, and wherein the switch is operatively configured to sequentially cycle between channels, to connect each of the amplified detector outputs in turn to the output port of the switch; an analog-to-digital converter (ADC) connected to the output port of the N1 switch, and operatively configured to sample a time dependent optical intensity signal input to the ADC from the switch, and generate a corresponding time dependent digital signal output from the ADC; and a microcontroller, connected to the N1 switch and the ADC, and operatively configured to: control an acquisition sequence by setting a sampling rate and a sample duration used by the ADC to provide a digital voltage data stream from each channel in turn; and make at least an average optical intensity value characterizing the voltage data stream available from each optical channel in turn at a digital output port of the microcontroller; such that the module provides an output at the digital output port of the microcontroller comprising N data values, characterizing the light incident on the N optical channels of the module; wherein the microcontroller is further operatively configured, when cycling through the channel output ports in turn, to, for each channel: compute an estimated flicker spectrum from the digital voltage data stream provided by the ADC for that channel; sample the estimated flicker spectrum for that channel at DC and at H select flicker sampling frequencies where H is an integer greater than 0; and make said H values of the estimated flicker spectrum, available at the digital output port of the microcontroller for each channel; such that the output provided by the module at the digital output port of the microcontroller comprises N*(1+H) data values, characterizing the light incident on the module.

2. The light detection module of claim 1, wherein computing the estimated flicker spectrum comprises executing a Fast Fourier Transform algorithm.

3. The light detection module of claim 1, wherein H=1, and the values of the estimated flicker spectrum provided for each channel are a value of the average light intensity and a value of amplitude at a predetermined flicker frequency.

4. The light detection module of claim 3, wherein the predetermined frequency is one of 100 Hz and 120 Hz.

5. The light detection module of claim 2, wherein the digital output of the microcontroller is connected to a microprocessor module with access to reference data on a plurality of artificial light sources; wherein the microprocessor is operatively configured to execute an unmixing algorithm to produce relative abundance information characterizing the mixture of artificial light sources; and wherein the reference data comprises N*(1+H) reference flicker spectrum intensity values previously measured, by detectors identical to, or calibrated with respect to, the N channel detectors, through filters identical to, or calibrated with respect to, the N channel filters, from received emission from the plurality of artificial light sources, and sampled at flicker spectrum frequency values identical to the H frequencies sampled in the incident light by the light detection module.

Description

BRIEF DESCRIPTIONS OF THE FIGURES

(1) FIG. 1. Illustration of flicker attributes flicker percent and flicker index.

(2) FIG. 2. Block diagram of a commercially available the ambient light sensor (PRIOR ART)

(3) FIG. 3A. An example of a flicker intensity spectrum measured from a high-pressure sodium lamp, with the lamp powered on.

(4) FIG. 3B. An example of a flicker intensity spectrum measured from the lamp of FIG. 3A, with the lamp powered off.

(5) FIG. 4. A light-pollution characterization (LPC) system according to one embodiment of the present invention.

(6) FIG. 5. An LPC module according to one embodiment of the present invention.

(7) FIG. 6. An LPC module according to another embodiment of the present invention.

(8) FIG. 7. Example of high-pressure mercury lamp optical emission spectra and selection of bandpass filters to be used to detect some of the emission peaks.

(9) FIG. 8. Example of high-pressure sodium lamp optical emission spectra and selection of bandpass filters to be used to detect some of the emission peaks.

(10) FIG. 9. Example of light-emitting diode lamp optical emission spectra and selection of bandpass filters to be used to detect some of the emission spectra characteristics.

(11) FIG. 10. An illustration of power balance in an illumination source's emission.

(12) FIG. 11. Calibration mode in a laboratory and normal operation mode of the LPC system in low Earth orbit.

(13) FIG. 12 A non-limiting example of F=20 filters that are a choice for selection.

(14) FIG. 13 Illustration of M=7 illumination sources spectra overlapped with the M=20 filters shown in FIG. 12.

(15) FIG. 14. Method for practicing the invention flow chart.

DETAILED DESCRIPTION OF THE INVENTION

Definitions

(16) The innovative system disclosed in this patent application is referred to as a light-pollution characterization system (LPC system) although it will be clear that this disclosed apparatus and methods employed can be applied to any other illumination unmixing task without departing from the spirit of the invention. The LPC system comprises a light detection module, a processing module, and a memory module with reference information for each illumination source that is to be detected.

(17) Measured light is assumed to comprise a mixture of M known and unknown illumination sources. The former (known) produce mutually uncorrelated light streams and the latter (unknown) are considered to be noise in the measurement. The method assumes that the known sources can be identified in the mixture if their characteristics (illumination source reference data for sources of the same types) are known prior to the unmixing analysis. If there are M known source types in the reference source data bank (matrix S), the result of the unmixing analysis is an array of real numbers x.sub.m (m=0, 1, 2, . . . M), each of which represents the fraction of a known illumination source (for m>0) in the measured mixture. The sum of all coefficients x.sub.m (x.sub.m) is unity.

(18) Single-Pixel Multi-Spectral Detection

(19) FIG. 5. illustrates one embodiment of a light detection module 500 that detects and analyzes incoming light 511 and delivers measurement results at a digital output 512. The illustrated embodiment comprises N optical channels 510, only one outlined in a dashed box for clarity, each optical channel n comprising an optical filter F.sub.n (501), optical lens L.sub.n (514), optical detector D.sub.n (502), and amplifier .sub.n (503), where n=1, 2, 3, . . . N. Only the components in just the first two channels are fully labeled, for clarity. While this embodiment shows the detector and amplifier as discrete components, in practice, some degree of amplification is usually built into photodetectors, and there may be no need for an additional amplifier 503. Incident light 511 is simultaneously and substantially equally incident on all N channel filters 501, lenses 514, and detectors 503. A portion of light incident on filter F.sub.n is transmitted through the filter F.sub.n, focused with lens L.sub.n, and detected on the detector D.sub.n, and converted to a current which is then coupled to the input of amplifier .sub.n resulting in a voltage at the output of amplifier .sub.n that is substantially proportional to the optical intensity detected by detector D.sub.n. The amplifier may be a transimpedance amplifier or an integrator, which may be gated.

(20) The outputs from N amplifiers 503 are coupled to N input ports of an analog switch module 504, the output of the analog switch 504 is then fed into an analog-to-digital converter 505 (ADC) and the digital information from the ADC 505 is fed by a digital interface 506 into a microcontroller 507 that processes the measurement data and delivers the result at the output terminal 512 in digital form. The microcontroller 507 is programmed to control the sampling times of the ADC 505 via electrical connection 513 and the channel selection on the analog switch 504 using electrical connection 508.

(21) FIG. 6 illustrates another embodiment of a light detection module 600 that detects and analyzes incoming light 601 and delivers measurement results at a digital output 612. The illustrated embodiment comprises N optical channels 610, each optical channel n comprising an optical filter F.sub.n (611), an optical lens L.sub.n (614), optical detector D.sub.n (612), and amplifier 613 denoted with .sub.n, where n=1, 2, 3, . . . N. As noted above, in some related embodiments, amplifier 613 may not be required. Incident light 601 is simultaneously and substantially equally incident on all N channel filters 611. A portion of light incident on filter F.sub.n is transmitted through the filter F.sub.n, focused on the detector D.sub.n by the lens L.sub.n. The light captured by each of the detector is converted to a current which is then coupled to the input of amplifier .sub.n resulting in a voltage at the output of amplifier .sub.n that is substantially proportional to the optical intensity detected by detector D.sub.n. The amplifier may be a transimpedance amplifier or an integrator, which may be gated. FIG. 6 illustrates an embodiment of the light measurement module that is similar as the one shown in FIG. 5, while not using an analog selection switch 504. In this embodiment, each of the N channels has its own dedicated analog-to-digital converter 615 so that the selection between the channels is performed by the micro-controller digitally.

(22) As in the previous embodiment, the system is operatively configured to acquire data by sampling all N channels 610 in time with a sampling rate and over a predetermined exposure time. The microcontroller triggers the measurement using trigger line 619 and acquires the data via the digital interface lines 616.

(23) Both embodiments illustrated in FIG. 5 and FIG. 6 are operatively configured to acquire data by sampling all N channels in time with a sampling rate and over a predetermined time interval which will be referred to as the exposure time. During the exposure time, the microcontrollers 507 and 617 will acquire P values that represent temporal response for each channel. This sequence is referred to as P.sub.k,n, where k is the ordinal number of the sample point in the sequence and n is the channel on which the sequence was taken. The microcontrollers are operatively configured to execute an FFT algorithm on P.sub.k,n sequence for each n of the N optical channels 501 to produce a flicker spectrum in digital form for each channel n and then sample this discrete flicker spectrum at predetermined flicker frequencies or frequency points. The total number of data points selected and stored is K, where it contains K.sub.DC DC values and K.sub.AC AC values which are sampled from the flicker spectrum at any predetermined frequency points. The number of channels corresponds to K.sub.DC.

(24) Optical Spectrum and Creation of the Reference Data Base

(25) The N optical filters are each designed to transmit a segment of the optical spectrum of incident light. The segments of the optical spectrum are selected to provide partially orthogonal measurements of the incident light.

(26) There are many commercially available illumination sources used in the world today which can be loosely categorized into (i) open flame or fire, (ii) incandescent sources in which a filament is heated to produce emission that approximates black-body radiation, (iii) gas discharge and fluorescent lamps, and (iv) solid-state lighting. To create a reference-spectrum base, one needs at least one spectrum of each of the illumination sources. Inasmuch as different types of light sources emit light with distinct spectral features, varying by design and manufacturing technology choices and manufacturing tolerance, present public illumination technologies span a wide range of light qualities and cost structures. This is particularly true with metal-halide lamps where the mixture of gases is often proprietary. This variability also applies to the solid-state sources as the ratio between blue and yellow portion of the spectra define the correlated color temperature and the phosphor choice is still under development and each manufacturer uses a different phosphor and a different blue light-emitting diode. The logic for selecting the wavelengths for the bandpass filters is illustrated with the help of FIGS. 7, 8, and 9.

(27) Measured emission spectra 700 of a high-pressure mercury (HPM) lamp, shown in FIG. 7, exhibits several emission lines 701 which are characteristic for the mercury discharge. Clearly, one approach to detect mercury discharge is to pass the emitted light through a sequence of bandpass filters, each with characteristics illustrated with 751, and detect the presence of energy of that line. In 750 the selected lines are denoted with wavelengths .sub.1, .sub.2, .sub.3, .sub.7, .sub.9.

(28) Measured emission spectra 800 of high-pressure sodium (HPS) lamp, shown in FIG. 8, exhibits several emission lines 801 which are characteristic for the sodium discharge. Similarly, detecting sodium discharge can be performed by passing the light emitted by the HPS lamp through a sequence of bandpass filters, each with characteristics illustrated with 851, and detecting the presence of energy of that line. In 850, the select lines are denoted with wavelengths .sub.6, .sub.8, .sub.10, .sub.11.

(29) Finally, detecting an LED emission is illustrated in FIG. 9 where an example of LED emission 900 is shown. The emission spectrum is continuous and exhibits at least two broad peaks, the blue peak 901 and the yellow peak 902. Clearly, detecting the power spectral density around these two peaks with a filter combination at .sub.4, .sub.5, .sub.10, shown with 950, would provide evidence of a light-emitting diode emission.

(30) With detectors placed behind filters with central wavelengths .sub.1 to .sub.11, one obtains 11 signal channels and, clearly, presence of signals at some wavelengths, as for example, .sub.1 and .sub.11, will uniquely detect the presence of the HPM and HPS lamps. When illumination from any of these light sources is measured using one of the LPC systems illustrated shown in FIG. 4A or 4B, we obtain a set of K measured values which include both DC and AC characteristics. In one embodiment, the filters are selected so that the central wavelengths .sub.n of a bandpass filter with bandwidth .sub.n coincide with the emission features of common illumination sources, as it has been shown in FIGS. 7, 8, and 9. In yet another embodiment, the filters are selected to encompass multiple peaks, while in another embodiment, the passband characteristics of at least one of the filters overlaps with at least one other filter. The number of wavelengths defines the number of the optical channels N.

(31) Absolute Emission Power and Power Consumption

(32) The power balance in a light source's emission is illustrated in FIG. 10. All of the electrical power P.sub.ELEC delivered to an illumination source is converted to power emitted in the electromagnetic spectrum P.sub.EM and to power dissipated as heat P.sub.HEAT. For example, a tungsten wire in an incandescent lamp reaches temperature above 2000 K and emits light in the visible part of the electromagnetic spectrum, while the lamp body heats to temperatures substantially below 1000 K. (For simplicity, we assume that the sections of the lamp that have temperature between these two extrema contribute negligibly to the electromagnetic emission spectrum.) Under this assumption, we divide the emission spectrum from any lamp into the light emitting part P.sub.EM [W] with temperature greater than 1000 K and the heat emitting fraction P.sub.HEAT [W] for which the body temperature is well below 1000 K and hence is referred to as heat dissipation. Clearly, P.sub.ELEC=P.sub.EM+P.sub.HEAT.

(33) A fraction of the P.sub.EM [W] power is emitted into the visible part of the spectrum, and we can refer to this fraction as the radiance [W] in the visible part of the spectrum or, more efficiently, we quantify the visible part of the spectrum as luminous flux [lm]. Luminous efficacy K.sub.r of an illumination source is measured in lumens per watt [lm/W] and it is the ratio of luminous flux L [lm] per every watt of electromagnetic power P.sub.EM generated for the purpose of light generation.

(34) Every illumination source (i.e., lamp) also dissipates heat P.sub.HEAT [W] while providing the power delivered to the electromagnetic spectrum P.sub.EM [W] for illumination. Since P.sub.HEAT is actually heat and by Planck's law any material body at elevated temperature emits radiation, we shall maintain that the heat can also be drawn as an electromagnetic spectrum emitted from the lamp.

(35) These metrics are illustrated in FIG. 10 where P.sub.EM is the area under the electromagnetic part of the spectrum that is emitted by the light source, and the area under the second curve at longer wavelengths equals the dissipated power P.sub.HEAT. The graphs are not to scale.

(36) The two metrics used to quantify these luminous-flux and power conversions for illumination sources are luminous efficacy of the radiation given by K.sub.r:=L/P.sub.EM [lm/W] (IEC 845-21-090), while luminous efficiency of a lamp is the ratio P.sub.EM/P.sub.ELEC. As an example, a halogen incandescent lamp has luminous efficacy of the radiation K.sub.r5.5 lm/W, while its luminous efficiency is less than 1%. An LED lamp can have K.sub.r130 lm/W and 13%.

(37) A goal of the LPC system is relate the measured quantity to at least one of the following illumination sources quantities: (a) the luminous flux coming from the Earth, (b) overall radiance coming from the Earth and (c) electrical power used to deliver the measured luminous flux. This would enable the production of a global map of the estimated quantity.

(38) Suppose we have an illumination source m that emits P.sub.EM(m,t) into space and we capture N different spectral segments of its spectrum each centered at .sub.n with bandwidth .sub.n. Each of the filters F.sub.n has its own insertion loss, each lens L.sub.n its own scattering losses and inefficiency, each detector D.sub.n behind its own responsivity, all of which are different for a different center wavelength n, and finally the transimpedance of the amplifiers A.sub.n. We shall merge the insertion loss, filter shape, lens losses, detector responsivities and amplifier gains into one effective responsivity R.sub.n that converts the incident power [W] into a dimensionless signal at the output of the amplifier. We use dimensionless quantity because once the hardware is built, it can always be related to a voltage or a current. Each of the elements Y.sub.n(t) of the measured vector Y is related to the power spectral density via:

(39) $\begin{array}{cc}{Y}_n\left(t\right)=\underset{0}{\overset{}{}}{R}_n\left(\right).Math.P_{EM}^{}\left(,t\right)d& \left(2\right)\end{array}$
where R.sub.n() [1/W] is the effective responsivity on channel n versus wavelength [nm] and P.sub.EM(,t) [W/nm] is the power spectral density of the electromagnetic radiation produced by the illumination source versus wavelength. The factor (beta) depends on how far away from the source we make the measurement of the emitted light. It will be high in the laboratory, but very small when the light from the Earth is measured on a satellite. Note that Y.sub.n(t) is a radiometric quantity and is not directly relatable to the luminous flux because it includes emission at the wavelengths that are outside of the visible range.

(40) With the simultaneous measurement of the light incident on all n channels, we obtain a vector Y(t). Every illumination source hence has a unique vector Y(t) that captures the shape of the spectral emission P.sub.EM(,t) incident on the sensor. We now estimate the flicker spectrum on each of the components of the vector Y.sub.n(t) and obtain Fourier components for each of n components. For each of the N channels of flicker spectra and for the K-component measurement vector, we select a certain number of desired components. The measurement vector now contains N elements that correspond to the DC values of detector signal, while we can select any number of the harmonic amplitudes of the signal from each channel N to further describe the captured signal. If for each channel we use the DC value and H harmonics, the total number of data points taken per illumination source will be K=N.Math.(1+H). The intensity of measured light is proportional to the sum of N DC components of the vector Y because the AC components do not contribute to the overall power of the light captured. We refer to this intensity measure as the N-element norm of the vector Y.

(41) The intensity of the light captured is proportional to the N-element norm of measurement vector Y, while the relative values of the components of the normalized vector y are used to determine the composition, namely, break the measured light into known illumination source components using the unmixing algorithm disclosed in this application. The N-element norm of measurement vector Y can then be related to the power dissipation on the surface of the Earth needed to achieve the measured intensity. The N-element norm and the normalized vectors are given as:

(42) $\begin{array}{cc}\begin{array}{cc}{.Math.Y.Math.}_{\left(N\right)}=\underset{n=1}{\overset{N}{.Math.}}{Y}_n& y_k=\frac{Y_k}{\underset{n=1}{\overset{N}{.Math.}}{Y}_n}\end{array}& \left(3\right)\end{array}$

(43) Note that N is the number of channels, but that k is the number of elements in the measurement vector Y and that the number of elements K in this vector is larger than the number of channels by the number of harmonics of interest in the captured light, hence 1kK and that K>N. Normalized measurement vector y has components denoted with y.sub.k.

(44) Note furthermore, that the DC values (elements of vector Y) are non-negative by definition, but that the Fourier amplitudes of the flicker spectra will generally be complex numbers and the argument of the first harmonic will depend on the timing when the spectrum was acquired relative to the phase of the power grid at the surface of the Earth. For each light-source all of the complex harmonic amplitudes should be rotated by the phase of the first harmonic amplitude to ensure that the first harmonic amplitude is a real positive number, while the amplitudes of harmonics with harmonic numbers greater than unity should be converted from complex numbers to real numbers with values equal to the absolute value of the phase-corrected harmonic amplitude.

(45) Variables

(46) N=number of channels (number of parallel filters, lenses and detectors). M=number of illumination sources considered and characterized to make the reference matrix; number of columns in the reference matrix. H=number of frequencies in the measured flicker spectrum where harmonic amplitude values are taken (H does not include the DC value). K=number of elements in the measurement vector Y: K=N(1+H), where 1 stands for the DC value; number of rows in the reference matrix
Unmixing Method Description

(47) The measurement of light emitted from the Earth results in a normalized vector y. The assumption is that this vector y is a linear superposition of known (reference) illumination sources y.sub.m (m=0, 1, 2, . . . M). This is generally referred to as linear mixing.

(48) $\begin{array}{cc}y=\underset{m=1}{\overset{M}{.Math.}}{s}_m{x}_m& \left(4\right)\end{array}$
where the elements of the m-row vector x=[x.sub.1 x.sub.2 . . . x.sub.M].sup.T are target abundances (the superscript T means transposed). The reference matrix S contains m columns, one column for each illumination source measured during the calibration. The reference matrix is shown below:

(49) $\begin{array}{cc}{s}_m=\left[\begin{array}{c}{s}_{1,m}\\ s_{2,m}\\ .Math.\\ s_{K,m}\end{array}\right].fwdarw.s=\left[\begin{array}{cccc}{s}_{1,1}& s_{1,2}& .Math.& s_{1,M}\\ .Math.& .Math.& & .Math.\\ s_{N,1}& s_{N,2}& .Math.& S_{N,M}\\ s_{N+1,1}& s_{N+1,2}& .Math.& s_{N+1,M}\\ .Math.& .Math.& & .Math.\\ s_{K,1}& s_{K,2}& .Math.& s_{K,M}\end{array}\right]\begin{array}{c}\}\mathrm{DC}\mathrm{components}\\ \}\mathrm{AC}\mathrm{components}\end{array}& \left(5\right)\end{array}$

(50) The method is a spectral unmixing algorithm to be used on the processing module 403, which relies on a reference matrix data bank (within memory module 404). The data are vectors of numbers, each of a length K, stored for each individual light source to which the device has been calibrated. As such, the reference data bank contains KM data points, where M is the number of light sources (number of targets).

(51) The matrix S is formed with M column vectors s.sub.m, where each column vector contains K elements (rows), each of the K elements is an output from the N channels measured on illumination source m of the M different illumination sources and H selected higher flicker frequency components of the temporal response, hence K=N(1+H). The components (rows) of the column vector s.sub.m are denoted s.sub.k,m, where 1kK. The elements of these vectors are non-negative and the N-element norm of the vector components of s.sub.m is equal to unity:

(52) $\begin{array}{cc}m\underset{n=1}{\overset{N}{.Math.}}{s}_{n,m}=1.& \left(5\right)\end{array}$

(53) Note that the s.sub.n,m are positive real numbers as was described above. The reference measurements are done for each m of the M illumination sources and a matrix S with K rows and M columns is formed. The reference illumination sources are also referred to as targets (to be consistent with the hyperspectral sub-pixel imaging terminology).

(54) The assumption is that the mixtures of light coming from the many sources are linear superpositions in the sense that the power from each of the available sources is weighted by some target abundance coefficient and that the individual contributions are added. For each measurement, the target abundance is given by the elements in an M-row vector x=[x.sub.1 x.sub.2 . . . x.sub.M].sup.T. The total intensity of the measurement is denoted |Y|.sub.N, while the target abundances x.sub.m sum to unity (all have a taxicab norm equal to 1):

(55) $\begin{array}{cc}\underset{m=1}{\overset{M}{.Math.}}{x}_m=1.& \left(6\right)\end{array}$

(56) Any measurement of an illumination mixture is expressed in terms of the total intensity P and the normalized distribution of channels given by the vector y. This is expressed ideally with the expression y=S.Math.x:

(57) $\begin{array}{cc}\left[\begin{array}{c}{y}_1\\ y_2\\ .Math.\\ y_k\end{array}\right]=\left[\begin{array}{cccc}{s}_{1,1}& s_{1,2}& .Math.& s_{1,M}\\ s_{2,1}& s_{2,2}& .Math.& s_{2,M}\\ .Math.& .Math.& & .Math.\\ s_{K,1}& s_{K,2}& .Math.& s_{K,M}\end{array}\right].Math.\left[\begin{array}{c}{x}_1\\ x_2\\ .Math.\\ x_M\end{array}\right]& \left(7\right)\end{array}$

(58) With the normalizations, the measured vector y must have its N-element norm equal to unity:

(59) $\begin{array}{cc}{y}_k=\underset{m=1}{\overset{M}{.Math.}}{s}_{k,m}{x}_m.fwdarw.\underset{k=1}{\overset{N}{.Math.}}{y}_k=\underset{m=1}{\overset{M}{.Math.}}{x}_m\left(\underset{n=1}{\overset{N}{.Math.}}{s}_{n,m}\right)=1.fwdarw.\underset{n=1}{\overset{N}{.Math.}}{s}_{n,m}=1& \left(8\right)\end{array}$

(60) Therefore, a system with K channels produces a (K1) vector y,

(61) 0$\begin{array}{cc}y=Sx& \left(9\right)\end{array}$

(62) Since y (K1) and S (KM) are measured quantities, we are seeking an estimate of the unknown vector of target abundances x (M1). Vector z is the best fit for x in the least square sense |ySx|>|ySz|, where | . . . | is the Euclidian norm of the vector between the | signs. The optimal target abundance vector z is S.sup.+y, where S.sup.+ (MK) is the Moore-Penrose inverse of matrix S. The error in the method is | S.sup.+yx|. This approach to linear square optimization is well-known in the art and can be found in publicly available literature.sup.4. .sup.4 Such as, for example, J. M. Ortega, Matrix Theory published by Plenum Press in New York (1987)

(63) The above-described approach is a least squares unconstrained method, but other constrained approaches can be found in publicly available literature and are shortly summarized below. All of the approaches are publicly available. Adding constraints to the unmixing approach is a natural addition since the target abundance quantities as well as reference spectra used in the unmixing are normalized and constrained. Namely, the target abundances x.sub.m sum to unity (all have a taxicab norm equal to 1), see Equation (6). Therefore, the optimal target abundance vector z coordinates also sum to unity and this vector can be found from the unconstrained optimization problem by seeking the vector of target abundances subject to the linear restriction: a.sup.Tx=1, where a is a unity column vector (M1) where each coordinate is equal to one.

(64) The second constraint is that the coordinates in the target abundance vector x must be non-negative numbers. Adding this constraint complicates finding the abundance vector since no closed-form solution has been found and an iterative algorithm must be employed to find the optimal target abundance vector z that is the best fit for x in the least square sense subject to the constraint that each coordinate of z is a non-negative number. Note that both the non-negativity constraint and the sum to one constraint can be included. The constrained method is a special case of a more general convex optimization problem.

(65) Interior-point methods are a certain class of iterative algorithms commonly employed to solve general convex optimization problems. It is known that the iterative approach converges towards the optimal solution as the methods satisfy Karush-Kuhn-Tucker (KKT) conditions, which are necessary derivative tests for the optimal solution.

(66) The most important advantage of using flicker harmonics in the inversion disclosed above is that it provides more information on each illumination source than was available from just DC data. This is evident from the following consideration: in the case when only DC light intensity is being captured by N filters, the reference matrix size is (NM) and there are three options: When N<M, the problem is underdefined and the estimate computation using the pseudoinverse will produce a large error largely independent of the choice of filters. We shall not use this case. When N=M, the problem is well defined provided that the determinant of the reference matrix is not zero. This is generally true because it is unlikely that captured the spectra are scaled versions of any other spectrum in the reference matrix. If there is no noise added to the measured data, zero error is achieved. When N>M, the system is over-defined and zero error is achieved if there is no noise.

(67) The relative size of N and M is bound by two opposing conditions. In order to achieve high unmixing accuracy, we need that NM. On the other hand, we need M to be as large as possible to capture largest possible number of illumination sources. Unfortunately, the number of filters N is limited by the size and weight of the module that needs to fit onto a nanosatellite. The ideal approach to resolving this problem is to effectively increase N without increasing the number of channels. This is precisely what adding the harmonic information does: the condition for high accuracy of the unmixing algorithm is in fact ensuring that KM. Consequently, adding H harmonics, effectively increases the number of channels and thereby provides higher accuracy relatively to DC only approaches. For example, if only the first flicker harmonic is captured in addition to the DC value on each of the N channels, the measurement vector is doubled and hence the number of reference sources in the reference matrix can be doubled. In this way, the reference matrix has double the number of rows and columns. This is a preferable approach in the case where the number of filters is limited, which is the case with the disclosed LPCM. As an example, using just one harmonic (H=1) will effectively double the number of channels, which is a great advantage for a module with a limited number of channels.

(68) Calibration

(69) The first step in practicing the invention is to create a reference matrix with reference illumination source data. This may be in the form of a matrix that has K rows and M columns and is populated with numbers that have been obtained by characterizing M different illumination sources using an N-channel LPC system that captures L harmonics of the incident light in each channel. The laboratory experimental environment is illustrated in FIG. 11(a) where the light-pollution characterization module (LPCM), which is a part of the LPC system is located in a laboratory and some distance .sub.1 from the illumination source also positioned in the laboratory and a measurement is being made of the filtered light in N-channels and flicker spectra of each channel analyzed to provide N DC values plus N.Math.H harmonics of the emitted light. In addition, the emitted light may also be measured using a color-meter to provide calibrated photometric characterization and the electrical power dissipation noted.

(70) Now, in FIG. 11(b) we analyze the situation when the LPCM is in orbit at a distance .sub.2 facing towards the Earth (nadir pointing) and the same light source is position at the surface of the Earth. Then, for any position of the satellite, we need an of the free-space loss . This can be done analytically by assuming that the light dispersion obeys ray optics and that we know the approximate emission pattern, for example, whether the emitter can be treated as a Lambertian emitter. In one embodiment, the filter selection is done in such a way that the filter passbands avoid wavelengths in which the atmosphere has highest absorption.

(71) Filter Selection

(72) The accuracy of the unmixing algorithm further depends on the filter selection, namely, how are the central wavelengths and full-width half-maxima selected. This is preferably done in the following manner: for each source m, we manually select a set .sub.m of wavelengths that are characteristic for that light sources .sub.m={.sub.1, .sub.2, . . . .sub.m} along with filter bandpass widths .sub.m={.sub.1, .sub.2, . . . .sub.m}. For each source m, hence, there is a number of filters associated with the source. We add filters at wavelengths where no source appears and obtain a complete set of possible filters with the total number of filters denoted with F. Each filter in this union has a defined center wavelength and full-width half-maximum (FWHM). FIG. 12 shows an exemplary list of 20 filters (F=20) from which one can select N filters. FIG. 13 shows the optical emission spectra of 7 different illumination sources (M=7) with the locations of the filter center wavelengths.

(73) Now select N filters from the defined F (the are .sub.NCF combination for this selection .sub.NC.sub.F is a binomial coefficient n choose k) for each combination perform a simple Monte Carlo method in selecting a random target abundance vector and computing the error |S.sup.+yx| between the estimated S.sup.+y and the initial target abundance x vectors. The algorithm then finds the filter combination that gives the lowest error |S.sup.+yx|.

(74) The far-right column of the table shown in FIG. 12 shows which combination of N=7 filters provides the lowest error when applies to the M=7 illumination sources shown in FIG. 12.

(75) Method of Practicing the Invention

(76) FIG. 14 illustrates the flow chart of the method of identifying illumination source in captured light. This method uses both the average and characteristic of the periodic fluctuation in the illumination to perform the unmixing. The steps 351 to 357 are referred to as calibration steps, while the steps 358 to 363 are referred to as normal operation steps. As indicated with the step 362, normal operation steps may be repeated any number of times, but require the calibration steps to be executed at least once for M>1 illumination sources. In calibration step 351, M different illumination sources are provided. The set of provided illumination sources contains representative categories of illumination sources that are expected to be detected in the normal operation mode. At the output from step 351, the illumination source counter m is set to the first illumination source: m=1. In step 352, illumination source m is powered on, and the electrical power delivered to the illumination source P.sub.m is noted. In step 353, the temporal variation of N optical spectrum segments is captured and converted into N electrical signals. In step 354, the N electrical signals are converted to digital form. Both the DC average and flicker spectrum of the captured data is estimated. The flicker spectrum is sampled at H predetermined frequencies plus the amplitude value at zero frequency to obtain N.Math.(1+H) amplitudes to be used in the data analysis. In step 355, the N.Math.(1+H) digital data points of illumination source m are stored in the reference data bank. In step 356, the method decides whether all of the M illumination sources have been characterized. If the answer is no, then the illumination source counter m is incremented in step 357, and the calibration method is continued by starting from the step 352, with now incremented illumination source counter m. If the calibration has completed with m=M, then the LPC system can be used in normal operation mode. In step 358, the LPC system is in an operating environment, which may be on a satellite measuring light-pollution from an orbit, an airplane, a vehicle or may be stationary. The LPC is exposed to light that is to be analyzed. In step 359, the LPC captures temporal fluctuations of the optical emission captured in at least in one of N spectral segments and converts this into electrical signal. In step 360, the electrical signal is converted into digital form. Each of the DC average and flicker spectrum of the captured data is estimated. The flicker spectrum is sampled at H predetermined frequencies plus one amplitude at zero frequency to obtain 1+H amplitudes to be used in the data analysis. In step 361, the N.Math.(1+H) digital measurement data and the measurement number or location where the measurement was done is saved in digital form in LPC memory. In step 362, the LPCM is ready to make another measurement and if this is required (YES), the method will start over from step 358 at a different location or measurement ordinal number. If LPC memory is full or there is an opportunity to download the memory contents (NO), the data is downloaded to a remote processor or the same processor in the LPC for further processing. In step 363, the downloaded data is processed for every location or measurement in that the unmixing algorithm is executed and target abundance coefficients x.sub.m for 1mM are determined and, optionally, electrical power P.sub.m consumed by each of the sources is used to compute electrical power dissipation per illumination source in the mixture. The method uses both DC and AC characteristics of conventional illumination sources to improve the accuracy of unmixing.

(77) There are variations on the presented method and the LPC system that can be implemented to optimize the unmixing without departing from the invention. Some of these variations are discussed below.

(78) Any number of harmonics of 100 Hz and 120 Hz may be selected for sampling during measurement stage in step 354 and 360, and the set of flicker-spectrum samples does not have to be a harmonic, but rather suitable selected frequency which is off interest. Natural light sources, such as, fire, volcanic eruption or Aurora Borealis will fluctuate and will not have peaks in the flicker spectrum, but may be detected by selecting an appropriate set of frequencies that capture the power spectral density of the light intensity fluctuation and from which the form of the power spectral density can be inferred.

(79) The filter passbands generally do not overlap, but it is also possible that for the purpose of improving signal to noise ratio, some of the filters have overlapping bands.

(80) For the purposes of this application, flicker frequency spectrum refers to oscillations or fluctuations in the light intensity with frequency components between 0 Hz and at least 20 kHz.

(81) It should be understood that there are many variations of the disclosed methods and systems described in the embodiments above without limiting the present invention, and that the scope of the present invention is to be determined by the following claims.