Estimation of the System Transfer Function for Certain Linear Systems

1. A method implemented on a computer system comprising a programmable processor and memory, the method for producing an estimate  of a system transfer matrix A for a linear system, wherein the linear system is characterized by y=Ax+e, where x is an input vector, e is a noise vector and y is the output vector resulting from input vector x, the method comprising: selecting a set of calibration inputs X cal =[x 1 , . . . x m ] where each of the x m is an input vector that is a calibration pattern; applying the calibration inputs X cal to the linear system; observing the calibration outputs Y cal resulting from the calibration inputs X cal according to Y=AX+E, where Y cal −=[y 1 , . . . , y m ] where each of the y m is an output vector resulting from the corresponding input calibration pattern x m , and E−=[e 1 , . . . , e m ] where each of the e m is a noise vector; identifying, by the computer system, a subspace of A, assuming that a column span of A equals a row span of A; calculating, by the computer system, an estimate  within the subspace of A, the estimate  calculated from the calibration outputs Y cal and calibration inputs X cal and subject to a constraint that a row span of the estimate  equals a column span of the estimate Â; observing outputs Y 1 produced by the linear system; and based on solving an inverse problem for the calculated estimate Â, estimating, by the computer system, corresponding inputs X 1 that produced the outputs Y 1 .

2. The computer-implemented method of claim 1 wherein the linear system is a plenoptic imaging system, the system transfer matrix A is a pupil image function (PIF) matrix for the plenoptic imaging system, the input x is an object to be imaged by the plenoptic imaging system, the output y is a plenoptic image resulting from object x, the calibration inputs X are calibration patterns, and the outputs Y are plenoptic images resulting from the calibration patterns X.

3. The computer-implemented method of claim 1 wherein identifying the subspace of A comprises performing a singular value decomposition (SVD) of Y.

4. The computer-implemented method of claim 3 wherein SVD(Y)=(WΛZ) and estimating  within the subspace is given by Â=Y({circumflex over (V)} T X) † {circumflex over (V)} T , where {circumflex over (V)}=W(:, 1:{circumflex over (k)}) and {circumflex over (k)} is greater than or equal to the rank of A.

5. The computer-implemented method of claim 1 wherein the calibration inputs X comprise random calibration inputs.

6. The computer-implemented method of claim 5 wherein the calibration inputs X is populated with i.i.d. Gaussian entries ˜N(0,1/m) where N( ) is the normal distribution and m is the number of calibration inputs.

7. The computer-implemented method of claim 5 wherein the calibration inputs X is populated with i.i.d symmetric Bernoulli entries taking ±sqrt(1/m) with probability ½.

8. The computer-implemented method of claim 5 wherein the calibration inputs X is populated with i.i.d. random variables taking ±sqrt(3/m) with probability ⅙ and zero with probability ⅔.

9. The computer-implemented method of claim 5 wherein the calibration inputs X is populated with entries that satisfy the Restricted Isometry Property (RIP) for matrices and/or vectors.

10. The computer-implemented method of claim 1 wherein  is an N×N matrix, and the number of calibration inputs in X is less than N/4.

11. The computer-implemented method of claim 1 wherein  is an N×N matrix, and the number of calibration inputs in X grows sub-linearly with N.

12. The computer-implemented method of claim 1 wherein calculating the estimate  assumes that A has four-quadrant symmetry.

13. The computer-implemented method of claim 1 wherein calculating the estimate  assumes that A has low rank.

14. The computer-implemented method of claim 1 wherein: the calibration inputs X is populated with entries that satisfy the Restricted Isometry Property (RIP) for matrices; and calculating the estimate  comprises: performing a singular value decomposition (SVD) of Y where SVD(Y)=(WΛZ); and calculating Â=Y({circumflex over (V)} T X) † {circumflex over (V)} T , where {circumflex over (V)}=W(:, 1:{circumflex over (k)}) and {circumflex over (k)} is greater than or equal to the rank of A.

15. A non-transitory computer-readable medium containing instructions that, when loaded into memory and executed by a processor, cause execution of a method for producing an estimate  of a system transfer matrix A for a linear system, wherein the linear system is characterized by y=Ax+e, where x is an input vector, e is a noise vector and y is the output vector resulting from input vector x, the method comprising: selecting a set of calibration inputs X cal =[x 1 , . . . , x m ] where each of the x m is an input vector that is a calibration pattern; applying the calibration inputs X cal to the linear system; observing the calibration outputs Y cal resulting from the calibration inputs X cal according to Y=AX+E, where Y cal −=[y 1 , . . . , y m ] where each of the y m is an output vector resulting from the corresponding input calibration pattern x m , and E−=[e 1 , . . . , e m ] where each of the e m is a noise vector; identifying, by the processor executing the instructions, a subspace of A, assuming that a column span of A equals a row span of A; calculating, by the processor executing the instructions, an estimate  within the subspace of A, the estimate  calculated from the calibration outputs Y cal and calibration inputs X cal and subject to a constraint that a row span of the estimate  equals a column span of the estimate Â; observing outputs Y 1 produced by the linear system; and based on solving an inverse problem for the calculated estimate Â, estimating, by the processor executing the instructions, corresponding inputs X 1 that produced the outputs Y 1 .

16. The non-transitory computer-readable medium of claim 15 wherein: identifying the subspace of A comprises performing a singular value decomposition (SVD) of Y according to SVD(Y)=(WΛZ); and estimating  within the subspace is given by Â=Y({circumflex over (V)} T X) † {circumflex over (V)} T , where {circumflex over (V)}=W(:, 1:{circumflex over (k)}) and {circumflex over (k)} is greater than or equal to the rank of A.

17. The non-transitory computer-readable medium of claim 15 wherein the calibration inputs X is populated with random entries that satisfy the Restricted Isometry Property (RIP) for matrices.

18. A system for producing an estimate  of a system transfer matrix A for a linear system, wherein the linear system is characterized by y=Ax+e, where x is an input vector, e is a noise vector and y is the output vector resulting from input vector x, the system comprising: a computer system comprising a memory and a processor; a selection module that selects a set of calibration inputs X cal [x 1 , . . . , x m ] where each of the x m is an input vector that is a calibration pattern; an observation module that captures the calibration outputs Y cal resulting from the application of the calibration inputs X cal to the linear system according to Y=AX+E, where Y cal −=[y 1 , . . . , y m ] where each of the is an output vector resulting from the corresponding input calibration pattern x m , and E−=[e 1 , . . . , e m ] where each of the e m is a noise vector wherein the processor executes instruction stored in the memory, that identifies a subspace of A, assuming that a column span of A equals a row span of A; and that further calculates an estimate A within the subspace of A, the estimate  calculated from the calibration outputs Y cal and calibration inputs X cal and subject to a constraint that a row span of the estimate  equals a column span of the estimate Â; the observation module further capturing outputs Y 1 produced by the linear system; and the processor executing further instruction stored in the memory, for solving an inverse problem for the calculated estimate  to estimate corresponding inputs X 1 that produced the outputs Y 1 .

19. The system of claim 18 wherein calculating the estimate A comprises the processor executing further instruction stored in the memory for: identifying the subspace of A comprises performing a singular value decomposition (SVD) of Y according to SVD(Y)=(WΛZ), and estimating  within the subspace is given by Â=Y({circumflex over (V)} T X) † {circumflex over (V)} T , where {circumflex over (V)}=W(:, 1:{circumflex over (k)}) and {circumflex over (k)} is greater than or equal to the rank of A.