1## Copyright (C) 2013 Martin Vogel <octave@martin-vogel.info>
2##
3## This program is free software; you can redistribute it and/or modify it
4## under the terms of the GNU General Public License as published by the
5## Free Software Foundation; either version 3 of the License, or (at your
6## option) any later version.
7##
8## This program is distributed in the hope that it will be useful, but WITHOUT
9## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
10## FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
11## for more details.
12##
13## You should have received a copy of the GNU General Public License
14## along with this program; see the file COPYING.  If not, see
15## <http://www.gnu.org/licenses/>.
16
17## -*- texinfo -*-
18## @deftypefn  {Function File} {@var{M} =} jones_lindiattenuator()
19## @deftypefnx {Function File} {@var{M} =} jones_lindiattenuator(@var{d})
20## @deftypefnx {Function File} {@var{M} =} jones_lindiattenuator(@var{px},@var{py})
21## @deftypefnx {Function File} {@var{M} =} jones_lindiattenuator(..., @var{mode})
22## Return the Jones matrix for a linear diattenuator at zero
23## rotation.
24##
25## @itemize @minus
26## @item @var{d} is the diattenuation of the element, i.e.
27## @code{d=(px-py)/(px+py)}. Reversibly, transmission in y direction
28## is @code{(1-d)/(1+d)}, if transmission in x direction is 1.
29## @item @var{px} is the transmittance in x direction.
30## @item @var{py} is the transmittance in y direction.
31## @item @var{mode} is a string defining the interpretation of
32## transmittance values: 'intensity' (default) or 'amplitude'.
33## @end itemize
34##
35## Arguments @var{d}, @var{px} or @var{py} can be passed as a scalar
36## or as a matrix or as a cell array. In the two latter cases, a cell
37## array @var{M} of Jones matrices is returned. The size of @var{M}
38## is set to the maximum of the parameters' size.
39##
40## References:
41##
42## @enumerate
43## @item E. Collett, Field Guide to Polarization,
44##       SPIE Field Guides vol. FG05, SPIE (2005). ISBN 0-8194-5868-6.
45## @item R. A. Chipman, "Polarimetry," chapter 22 in Handbook of Optics II,
46##       2nd Ed, M. Bass, editor in chief (McGraw-Hill, New York, 1995)
47## @item @url{http://en.wikipedia.org/wiki/Jones_calculus, "Jones calculus"},
48##       last retrieved on Jan 13, 2014.
49## @end enumerate
50##
51## @seealso{jones_linpolarizer}
52## @end deftypefn
53
54function JM = jones_lindiattenuator(varargin)
55
56  if nargin<1
57
58    px = 1;
59    py = 1;
60    was_cell = false;
61
62  elseif nargin==1
63
64    px = varargin{1};
65    [px, was_cell] = __c2n__(px, 0);
66    py = (1-px)./(1+px);
67    px(:) = 1;
68
69  else
70
71    px = varargin{1};
72    [px, was_cellx] = __c2n__(px, 1);
73
74    py = varargin{2};
75    [py, was_celly] = __c2n__(py, 1);
76
77    was_cell = was_cellx || was_celly;
78
79  end
80
81  % check mode
82  s_function = @s_lindiattenuator_int;
83  if nargin>=2 && ischar(varargin{end})
84    if strncmpi(varargin{end},'amp',3)
85      s_function = @s_lindiattenuator_amp;
86    end
87  end
88
89  % any matrix in parameters?
90  if (any([numel(px),numel(py)] > 1)) || was_cell
91
92    % adjust dimensions, i.e. fill missing dimensions with 1
93    spx = size(px);
94    spy = size(py);
95
96    maxdim = max([length(spx),length(spy)]);
97    if length(spx) < maxdim
98      spx = [spx, ones(1,maxdim-length(spx))];
99    end
100    if length(spy) < maxdim
101      spy = [spy, ones(1,maxdim-length(spy))];
102    end
103
104    % generate Jones matrices
105    maxsize = max([spx;spy]);
106    JM = cell(maxsize);
107    JM_subs = cell(1,ndims(JM));
108    numelJM = numel(JM);
109
110    % flatten parameter arrays
111    px = px(:);
112    py = py(:);
113    numelpx = numel(px);
114    numelpy = numel(py);
115
116    for jmi=1:numelJM
117      [JM_subs{:}] = ind2sub(size(JM),jmi);
118      JM{JM_subs{:}} = s_function(px(mod(jmi-1,numelpx)+1), py(mod(jmi-1,numelpy)+1));
119    end
120
121  else
122
123    JM = s_function(px, py);
124
125  end
126
127end
128
129% helper function
130function JM = s_lindiattenuator_amp(px_in_amplitude_domain,py_in_amplitude_domain)
131
132  JM = zeros(2,2);
133  JM(1,1) = px_in_amplitude_domain;
134  JM(2,2) = py_in_amplitude_domain;
135
136end
137
138% helper function 2
139function JM = s_lindiattenuator_int(kx_in_intensity_domain,ky_in_intensity_domain)
140
141  JM = zeros(2,2);
142  JM(1,1) = sqrt(kx_in_intensity_domain);
143  JM(2,2) = sqrt(ky_in_intensity_domain);
144
145end
146
147
148%!test
149%! % test default return value
150%! A = jones_lindiattenuator();
151%! R = A*A-A;
152%! assert(norm(R,inf), 0, 1e-9);
153%!
154%!test
155%! % test equality of providing diattenuation or kx and ky
156%! d = rand(1, 1);
157%! kx = 1;
158%! ky = (1-d)./(1+d);
159%! A1 = jones_lindiattenuator(d);
160%! A2 = jones_lindiattenuator(kx, ky);
161%! R = A1-A2;
162%! assert(norm(R,inf), 0, 1e-9);
163%!
164%!test
165%! % test serial application of linear diattenuators
166%! kx = rand(1, 1);
167%! ky = rand(1, 1);
168%! A1 = jones_lindiattenuator(kx,ky);
169%! A2 = jones_lindiattenuator(kx*kx,ky*ky);
170%! R = A1*A1-A2;
171%! assert(norm(R,inf), 0, 1e-9);
172%!
173%!test
174%! % another test of serial application of retarder elements
175%! kx1 = rand(1, 1);
176%! kx2 = rand(1, 1);
177%! ky1 = rand(1, 1);
178%! ky2 = rand(1, 1);
179%! A1 = jones_lindiattenuator(kx1,ky1);
180%! A2 = jones_lindiattenuator(kx2,ky2);
181%! A12 = jones_lindiattenuator(kx1*kx2,ky1*ky2);
182%! R = A1*A2-A12;
183%! assert(norm(R,inf), 0, 1e-9);
184%!
185%!test
186%! % test mode
187%! kx = rand(1, 1);
188%! ky = rand(1, 1);
189%! A1 = jones_lindiattenuator(kx,ky);
190%! A2 = jones_lindiattenuator(kx,ky,'int');
191%! A3 = jones_lindiattenuator(sqrt(kx),sqrt(ky),'amp');
192%! R1 = A1-A2;
193%! R2 = A1-A3;
194%! assert(norm(R1,inf)+norm(R2,inf), 0, 1e-9);
195%!
196%!test
197%! % test correct size of return values
198%! for dim = 1:5
199%!   asize = randi([1 4], 1, dim);
200%!   R = rand(asize);
201%!   if numel(R) == 1
202%!     R = {R};
203%!   end
204%!   rsize = size(R);
205%!   C = jones_lindiattenuator(R);
206%!   csize = size(C);
207%!   assert(rsize == csize);
208%! end
209%!
210%!test
211%! % another test correct size of return values
212%! kx = rand(3,4,5);
213%! ky = rand(5,4,3);
214%! C = jones_lindiattenuator(kx,ky);
215%! csize = size(C);
216%! assert(csize == [5,4,5]);
217
218
219