ARTS  2.4.0(git:4fb77825)
m_sensor.cc
Go to the documentation of this file.
1 /* Copyright (C) 2003-2012
2  Patrick Eriksson <Patrick.Eriksson@chalmers.se>
3  Stefan Buehler <sbuehler(at)ltu.se>
4 
5  This program is free software; you can redistribute it and/or modify it
6  under the terms of the GNU General Public License as published by the
7  Free Software Foundation; either version 2, or (at your option) any
8  later version.
9 
10  This program is distributed in the hope that it will be useful,
11  but WITHOUT ANY WARRANTY; without even the implied warranty of
12  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13  GNU General Public License for more details.
14 
15  You should have received a copy of the GNU General Public License
16  along with this program; if not, write to the Free Software
17  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
18  USA. */
19 
20 /*===========================================================================
21  === File description
22  ===========================================================================*/
23 
35 /*===========================================================================
36  === External declarations
37  ===========================================================================*/
38 
39 #include <algorithm>
40 #include <cmath>
41 #include <stdexcept>
42 #include <string>
43 #include "arts.h"
44 #include "auto_md.h"
45 #include "check_input.h"
46 #include "interpolation_poly.h"
47 #include "m_select.h"
48 #include "math_funcs.h"
49 #include "messages.h"
50 #include "ppath.h"
51 #include "sensor.h"
52 #include "sorting.h"
53 #include "special_interp.h"
54 #include "xml_io.h"
55 
56 extern const Numeric PI;
57 extern const Numeric NAT_LOG_2;
58 extern const Numeric DEG2RAD;
59 extern const Numeric RAD2DEG;
60 extern const Index GFIELD1_F_GRID;
61 extern const Index GFIELD4_FIELD_NAMES;
62 extern const Index GFIELD4_F_GRID;
63 extern const Index GFIELD4_ZA_GRID;
64 extern const Index GFIELD4_AA_GRID;
65 extern const Numeric SPEED_OF_LIGHT;
66 
67 /*===========================================================================
68  === The functions (in alphabetical order)
69  ===========================================================================*/
70 
71 /* Workspace method: Doxygen documentation will be auto-generated */
74  GriddedField4& r,
76  const Index& n_za_grid,
77  const Numeric& fwhm,
78  const Numeric& xwidth_si,
79  const Numeric& dx_si,
80  const Verbosity& verbosity) {
81  if (dx_si > xwidth_si)
82  throw runtime_error("It is demanded that dx_si <= xwidth_si.");
83 
84  antenna_dim = 1;
85 
86  antenna_dlos.resize(1, 1);
87  antenna_dlos(0, 0) = 0.0;
88 
89  antenna_responseGaussian(r, fwhm, xwidth_si, dx_si, 0, verbosity);
90 
91  // za grid for response
93  const Index nr = r_za_grid.nelem();
94 
95  // Cumulative integral of response (factor /2 skipped, but does not matter)
96  Vector cumtrapz(nr);
97  cumtrapz[0] = 0;
98  for (Index i = 1; i < nr; i++) {
99  cumtrapz[i] = cumtrapz[i - 1] + r.data(0, 0, i - 1, 0) + r.data(0, 0, i, 0);
100  }
101 
102  // Equally spaced vector between end points of cumulative sum
103  Vector csp;
104  nlinspace(csp, cumtrapz[0], cumtrapz[nr - 1], n_za_grid);
105 
106  // Get mblock_za_grid by interpolation
107  mblock_dlos_grid.resize(n_za_grid, 1);
108  ArrayOfGridPos gp(n_za_grid);
109  gridpos(gp, cumtrapz, csp);
110  Matrix itw(n_za_grid, 2);
111  interpweights(itw, gp);
112  interp(mblock_dlos_grid(joker, 0), itw, r_za_grid, gp);
113 }
114 
115 
116 
117 /* Workspace method: Doxygen documentation will be auto-generated */
123  const Index& atmosphere_dim,
124  const Verbosity& verbosity) {
125  // Sizes
126  const Index nmblock = sensor_pos.nrows();
127  const Index nant = antenna_dlos.nrows();
128 
129  // Input checks
130  chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
131  chk_if_in_range("antenna_dim", antenna_dim, 1, 2);
132  //
133  if (sensor_pos.ncols() != atmosphere_dim)
134  throw runtime_error(
135  "The number of columns of sensor_pos must be "
136  "equal to the atmospheric dimensionality.");
137  if (atmosphere_dim <= 2 && sensor_los.ncols() != 1)
138  throw runtime_error("For 1D and 2D, sensor_los shall have one column.");
139  if (atmosphere_dim == 3 && sensor_los.ncols() != 2)
140  throw runtime_error("For 3D, sensor_los shall have two columns.");
141  if (sensor_los.nrows() != nmblock) {
142  ostringstream os;
143  os << "The number of rows of sensor_pos and sensor_los must be "
144  << "identical, but sensor_pos has " << nmblock << " rows,\n"
145  << "while sensor_los has " << sensor_los.nrows() << " rows.";
146  throw runtime_error(os.str());
147  }
148  if (antenna_dim == 2 && atmosphere_dim < 3) {
149  throw runtime_error("If *antenna_dim* is 2, *atmosphere_dim* must be 3.");
150  }
151  if (antenna_dlos.empty()) throw runtime_error("*antenna_dlos* is empty.");
152  if (antenna_dlos.ncols() < 1 || antenna_dlos.ncols() > 2)
153  throw runtime_error("*antenna_dlos* must have one or 2 columns.");
154  if (atmosphere_dim < 3 && antenna_dlos.ncols() == 2)
155  throw runtime_error(
156  "*antenna_dlos* can only have two columns for 3D atmosphers.");
157 
158  // Claculate new sensor_pos and sensor_los
159  const Matrix pos_copy = sensor_pos;
160  const Matrix los_copy = sensor_los;
161  //
162  sensor_pos.resize(nmblock * nant, pos_copy.ncols());
163  sensor_los.resize(nmblock * nant, los_copy.ncols());
164  //
165  for (Index ib = 0; ib < nmblock; ib++) {
166  for (Index ia = 0; ia < nant; ia++) {
167  const Index i = ib * nant + ia;
168 
169  sensor_pos(i, joker) = pos_copy(ib, joker);
170  sensor_los(i, joker) = los_copy(ib, joker);
171 
172  sensor_los(i, 0) += antenna_dlos(ia, 0);
173  if (antenna_dlos.ncols() == 2) sensor_los(i, 1) += antenna_dlos(ia, 1);
174  }
175  }
176 
177  // Set other variables
179  //
180  antenna_dlos.resize(1, 1);
181  antenna_dlos = 0;
182 }
183 
184 
185 
186 /* Workspace method: Doxygen documentation will be auto-generated */
189  const Verbosity& verbosity) {
190  CREATE_OUT2;
191 
192  out2 << " Sets the antenna dimensionality to 1.\n";
193  antenna_dim = 1;
194 
195  out2 << " Sets *mblock_dlos_grid* to have one row with value 0.\n";
196  mblock_dlos_grid.resize(1, 1);
197  mblock_dlos_grid = 0;
198 }
199 
200 
201 
202 /* Workspace method: Doxygen documentation will be auto-generated */
204  const Numeric& fwhm,
205  const Numeric& xwidth_si,
206  const Numeric& dx_si,
207  const Index& do_2d,
208  const Verbosity&) {
209  if (dx_si > xwidth_si)
210  throw runtime_error("It is demanded that dx_si <= xwidth_si.");
211 
212  Vector x, y;
213  gaussian_response_autogrid(x, y, 0, fwhm, xwidth_si, dx_si);
214 
215  r.set_name("Antenna response");
216 
217  r.set_grid_name(0, "Polarisation");
218  r.set_grid(0, {"NaN"});
219 
220  r.set_grid_name(1, "Frequency");
221  r.set_grid(1, Vector(1, -999));
222 
223  r.set_grid_name(2, "Zenith angle");
224  r.set_grid(2, x);
225 
226  // common for 1D and 2D
227  r.set_grid_name(3, "Azimuth angle");
228  const Index n = y.nelem();
229  //
230  if (do_2d) {
231  r.set_grid(3, x);
232  r.data.resize(1, 1, n, n);
233 
234  // The code below follows the function *gaussian_response*
235  const Numeric si = fwhm / (2 * sqrt(2 * NAT_LOG_2));
236  const Numeric a = 1 / (si * sqrt(2 * PI));
237 
238  for (Index z = 0; z < n; z++) {
239  for (Index b = 0; b < n; b++) {
240  r.data(0, 0, z, b) =
241  a * exp(-0.5 * pow(sqrt(x[z] * x[z] + x[b] * x[b]) / si, 2.0));
242  }
243  }
244  } else {
245  r.set_grid(3, Vector(1, 0));
246  r.data.resize(1, 1, n, 1);
247  r.data(0, 0, joker, 0) = y;
248  }
249 }
250 
251 
252 
253 /* Workspace method: Doxygen documentation will be auto-generated */
255  const Numeric& leff,
256  const Numeric& xwidth_si,
257  const Numeric& dx_si,
258  const Index& nf,
259  const Numeric& fstart,
260  const Numeric& fstop,
261  const Index& do_2d,
262  const Verbosity& verbosity) {
263  if (dx_si > xwidth_si)
264  throw runtime_error("It is demanded that dx_si <= xwidth_si.");
265 
266  // Calculate response for highest frequency, with xwidth_si scaled from
267  // fstart
268  Vector x, y;
269  Numeric fwhm = RAD2DEG * SPEED_OF_LIGHT / (leff * fstop);
271  x, y, 0, fwhm, (fstop / fstart) * xwidth_si, dx_si);
272 
273  r.set_name("Antenna response");
274 
275  r.set_grid_name(0, "Polarisation");
276  r.set_grid(0, {"NaN"});
277 
278  Vector f_grid;
279  VectorNLogSpace(f_grid, nf, fstart, fstop, verbosity);
280  r.set_grid_name(1, "Frequency");
281  r.set_grid(1, f_grid);
282 
283  r.set_grid_name(2, "Zenith angle");
284  r.set_grid(2, x);
285 
286  // common for 1D and 2D
287  r.set_grid_name(3, "Azimuth angle");
288  const Index n = y.nelem();
289  //
290  if (do_2d) {
291  r.set_grid(3, x);
292  r.data.resize(1, nf, n, n);
293 
294  for (Index i = 0; i < nf; i++) {
295  fwhm = RAD2DEG * SPEED_OF_LIGHT / (leff * f_grid[i]);
296  const Numeric si = fwhm / (2 * sqrt(2 * NAT_LOG_2));
297  const Numeric a = 1 / (si * sqrt(2 * PI));
298 
299  for (Index z = 0; z < n; z++) {
300  for (Index b = 0; b < n; b++) {
301  r.data(0, i, z, b) =
302  a * exp(-0.5 * pow(sqrt(x[z] * x[z] + x[b] * x[b]) / si, 2.0));
303  }
304  }
305  }
306  } else {
307  r.set_grid(3, Vector(1, 0));
308  r.data.resize(1, nf, n, 1);
309  //
310  r.data(0, nf - 1, joker, 0) = y;
311  //
312  for (Index i = 0; i < nf - 1; i++) {
313  fwhm = RAD2DEG * SPEED_OF_LIGHT / (leff * f_grid[i]);
314  gaussian_response(y, x, 0, fwhm);
315  r.data(0, i, joker, 0) = y;
316  }
317  }
318 }
319 
320 
321 
322 /* Workspace method: Doxygen documentation will be auto-generated */
324  const Numeric& resolution,
325  const Verbosity&) {
326  r.resize(1);
327  r[0].set_name("Backend channel response function");
328 
329  Vector x(2);
330 
331  r[0].set_grid_name(0, "Frequency");
332  x[1] = resolution / 2.0;
333  x[0] = -x[1];
334  r[0].set_grid(0, x);
335 
336  r[0].data.resize(2);
337  r[0].data[0] = 1 / resolution;
338  r[0].data[1] = r[0].data[0];
339 }
340 
341 
342 
343 /* Workspace method: Doxygen documentation will be auto-generated */
345  const Vector& fwhm,
346  const Vector& xwidth_si,
347  const Vector& dx_si,
348  const Verbosity&) {
349  if ((fwhm.nelem() != xwidth_si.nelem() && xwidth_si.nelem() != 1) ||
350  (fwhm.nelem() != dx_si.nelem() && dx_si.nelem() != 1)) {
351  std::ostringstream os;
352  os << "*xwidth_si* and *dx_si* must have one element or the same number of\n";
353  os << "elements as *fwhm*.";
354  throw std::runtime_error(os.str());
355  }
356 
357  Index nchannels = fwhm.nelem();
358  r.resize(nchannels);
359 
360  Vector x, y;
361  Numeric this_xwidth_si = xwidth_si[0];
362  Numeric this_dx_si = dx_si[0];
363  for (Index i = 0; i < nchannels; i++) {
364  if (xwidth_si.nelem() > 1) this_xwidth_si = xwidth_si[i];
365 
366  if (dx_si.nelem() > 1) this_dx_si = dx_si[i];
367 
368  gaussian_response_autogrid(x, y, 0, fwhm[i], this_xwidth_si, this_dx_si);
369 
370  r[i].set_name("Backend channel response function");
371 
372  r[i].set_grid_name(0, "Frequency");
373  r[i].set_grid(0, x);
374 
375  const Index n = y.nelem();
376  r[i].data.resize(n);
377  for (Index j = 0; j < n; j++) r[i].data[j] = y[j];
378  }
379 }
380 
381 
382 
383 /* Workspace method: Doxygen documentation will be auto-generated */
384 void f_gridFromSensorAMSU( // WS Output:
385  Vector& f_grid,
386  // WS Input:
387  const Vector& lo,
388  const ArrayOfVector& f_backend,
390  // Control Parameters:
391  const Numeric& spacing,
392  const Verbosity& verbosity) {
393  CREATE_OUT2;
394  CREATE_OUT3;
395 
396  // Find out how many channels we have in total:
397  // f_backend is an array of vectors, containing the band frequencies for each Mixer.
398  Index n_chan = 0;
399  for (Index i = 0; i < f_backend.nelem(); ++i)
400  for (Index s = 0; s < f_backend[i].nelem(); ++s) ++n_chan;
401 
402  // Checks on input quantities:
403 
404  // There must be at least one channel:
405  if (n_chan < 1) {
406  ostringstream os;
407  os << "There must be at least one channel.\n"
408  << "(The vector *lo* must have at least one element.)";
409  throw runtime_error(os.str());
410  }
411 
412  // Is number of LOs consistent in all input variables?
413  if ((f_backend.nelem() != lo.nelem()) ||
414  (backend_channel_response.nelem() != lo.nelem())) {
415  ostringstream os;
416  os << "Variables *lo_multi*, *f_backend_multi* and *backend_channel_response_multi*\n"
417  << "must have same number of elements (number of LOs).";
418  throw runtime_error(os.str());
419  }
420 
421  // Is number of bands consistent for each LO?
422  for (Index i = 0; i < f_backend.nelem(); ++i)
423  if (f_backend[i].nelem() != backend_channel_response[i].nelem()) {
424  ostringstream os;
425  os << "Variables *f_backend_multi* and *backend_channel_response_multi*\n"
426  << "must have same number of bands for each LO.";
427  throw runtime_error(os.str());
428  }
429 
430  // Start the actual work.
431 
432  // We construct the necessary input for function find_effective_channel_boundaries,
433  // which will identify channel boundaries, taking care of overlaping channels.
434 
435  // A "flat" vector of nominal band frequencies (two for each AMSU channel):
436  Vector f_backend_flat(2 * n_chan);
437 
438  // A "flat" list of channel response functions (two for each AMSU channel)
439  ArrayOfGriddedField1 backend_channel_response_flat(2 * n_chan);
440 
441  // Counts position inside the flat arrays during construction:
442  Index j = 0;
443 
444  for (Index i = 0; i < f_backend.nelem(); ++i)
445  for (Index s = 0; s < f_backend[i].nelem(); ++s) {
446  const GriddedField1& this_grid = backend_channel_response[i][s];
447  const Numeric this_f_backend = f_backend[i][s];
448 
449  // Signal sideband:
450  f_backend_flat[j] = this_f_backend;
451  backend_channel_response_flat[j] = this_grid;
452  j++;
453 
454  // Image sideband:
455  Numeric offset = this_f_backend - lo[i];
456  Numeric f_image = lo[i] - offset;
457  f_backend_flat[j] = f_image;
458  backend_channel_response_flat[j] = this_grid;
459  j++;
460  }
461 
462  // We build up a total list of absolute frequency ranges for
463  // both signal and image sidebands:
464  Vector fmin(2 * n_chan), fmax(2 * n_chan);
465 
466  // We have to add some additional margin at the band edges,
467  // otherwise the instrument functions are not happy. Define
468  // this in terms of the grid spacing:
469  Numeric delta = 1 * spacing;
470 
471  // Call subfunction to do the actual work of merging overlapping
472  // channels and identifying channel boundaries:
474  fmax,
475  f_backend_flat,
476  backend_channel_response_flat,
477  delta,
478  verbosity);
479 
480  // Create f_grid_array. This is an array of Numeric, so that we
481  // can use the STL push_back function.
482  ArrayOfNumeric f_grid_array;
483 
484  for (Index i = 0; i < fmin.nelem(); ++i) {
485  // Band width:
486  const Numeric bw = fmax[i] - fmin[i];
487 
488  // How many grid intervals do I need?
489  const Numeric npf = ceil(bw / spacing);
490 
491  // How many grid points to store? - Number of grid intervals
492  // plus 1.
493  const Index npi = (Index)npf + 1;
494 
495  // Create the grid for this band:
496  Vector grid;
497  nlinspace(grid, fmin[i], fmax[i], npi);
498 
499  out3 << " Band range " << i << ": " << grid << "\n";
500 
501  // Append to f_grid_array:
502  f_grid_array.reserve(f_grid_array.nelem() + npi);
503  for (Index s = 0; s < npi; ++s) f_grid_array.push_back(grid[s]);
504  }
505 
506  // Copy result to output vector:
507  f_grid = f_grid_array;
508 
509  out2 << " Total number of frequencies in f_grid: " << f_grid.nelem() << "\n";
510  // cout << "f_grid: " << f_grid << "\n";
511 }
512 
513 
514 
515 /* Workspace method: Doxygen documentation will be auto-generated */
516 void f_gridFromSensorAMSUgeneric( // WS Output:
517  Vector& f_grid,
518  // WS Input:
519  const ArrayOfVector& f_backend_multi, // Center frequency for each channel
521  // Control Parameters:
522  const Numeric& spacing,
523  const Vector& verbosityVect,
524  const Verbosity& verbosity) {
525  CREATE_OUT2;
526  CREATE_OUT3;
527  Index numFpoints = 0;
528  // Calculate the total number of frequency points
529  for (Index idx = 0; idx < backend_channel_response_multi.nelem(); ++idx) {
530  for (Index idy = 0; idy < backend_channel_response_multi[idx].nelem();
531  ++idy) {
532  numFpoints += backend_channel_response_multi[idx][idy].get_grid_size(0);
533  }
534  }
535 
536  Index numChan = backend_channel_response_multi.nelem(); // number of channels
537 
538  // Checks on input quantities:
539 
540  // There must be at least one channel:
541  if (numChan < 1) {
542  ostringstream os;
543  os << "There must be at least one channel.\n"
544  << "(The vector *lo* must have at least one element.)";
545  throw runtime_error(os.str());
546  }
547 
548  // Start the actual work.
549  // We construct the necessary input for function find_effective_channel_boundaries,
550  // which will identify channel boundaries, taking care of overlaping channels.
551 
552  // A "flat" vector of nominal band frequencies (one for each AMSU channel):
553  Vector f_backend_flat(numChan);
554  // A "flat" list of channel response functions (one for each AMSU channel)
555  ArrayOfGriddedField1 backend_channel_response_nonflat(numChan);
556 
557  // Save the correct verbosity value to each passband
558  Vector FminVerbosityVect(numFpoints);
559  Vector FmaxVerbosityVect(numFpoints);
560  Vector VerbosityValVect(numFpoints);
561  Index VerbVectIdx = 0;
562 
563  for (Index i = 0; i < f_backend_multi.nelem(); ++i) {
564  for (Index ii = 0; ii < f_backend_multi[i].nelem(); ++ii) {
565  const GriddedField1& this_grid = backend_channel_response_multi[i][ii];
566  const Numeric this_f_backend = f_backend_multi[i][ii];
567  // Signal passband :
568  f_backend_flat[i] = this_f_backend;
569  backend_channel_response_nonflat[i] = this_grid;
570  for (Index idy = 1;
571  idy < backend_channel_response_multi[i][ii].get_grid_size(0);
572  ++idy) {
573  if ((backend_channel_response_multi[i][ii].data[idy - 1] == 0) &&
574  (backend_channel_response_multi[i][ii].data[idy] > 0)) {
575  FminVerbosityVect[VerbVectIdx] =
576  f_backend_multi[i][ii] +
577  backend_channel_response_multi[i][ii].get_numeric_grid(0)[idy];
578  VerbosityValVect[VerbVectIdx] = verbosityVect[i];
579  }
580  if ((backend_channel_response_multi[i][ii].data[idy - 1] > 0) &&
581  (backend_channel_response_multi[i][ii].data[idy] == 0)) {
582  FmaxVerbosityVect[VerbVectIdx] =
583  f_backend_multi[i][ii] +
584  backend_channel_response_multi[i][ii].get_numeric_grid(0)[idy];
585  VerbVectIdx++;
586  }
587  }
588  }
589  }
590  // We build up a total list of absolute frequency ranges for all passbands
591  // Code reused from the function "f_gridFromSensorAMSU"
592  Vector fmin,
593  fmax; // - these variables will be resized, therefore len(1) is enough for now.,
594 
595  // We have to add some additional margin at the band edges,
596  // otherwise the instrument functions are not happy. Define
597  // this in terms of the grid spacing:
598  const Numeric delta = 10;
599 
600  // Call subfunction to do the actual work of merging overlapping
601  // channels and identifying channel boundaries:
603  fmax,
604  f_backend_flat,
605  backend_channel_response_nonflat,
606  delta,
607  verbosity);
608 
609  // Create f_grid_array. This is an array of Numeric, so that we
610  // can use the STL push_back function.
611  ArrayOfNumeric f_grid_array;
612 
613  for (Index i = 0; i < fmin.nelem(); ++i) {
614  // Bandwidth:
615  const Numeric bw = fmax[i] - fmin[i];
616  Numeric npf = ceil(bw / spacing); // Set a default value
617 
618  // How many grid intervals do I need?
619  Index verbIdx = 0;
620  if (verbosityVect.nelem() > 0) {
621  // find the grid needed for the particular part of passband
622  for (Index ii = 0; ii < VerbVectIdx; ++ii) {
623  if ((FminVerbosityVect[ii] >= fmin[i]) &&
624  (FmaxVerbosityVect[ii] <= fmax[i])) {
625  if (verbIdx == 0) {
626  verbIdx = ii;
627  } else {
628  if (VerbosityValVect[ii] < VerbosityValVect[verbIdx]) {
629  verbIdx = ii;
630  }
631  }
632  }
633  }
634  if (spacing > VerbosityValVect[verbIdx]) {
635  npf = ceil(
636  bw / VerbosityValVect[verbIdx]); // is the default value to coarse?
637  } else {
638  npf = ceil(bw / spacing); // Default value
639  }
640  }
641 
642  // How many grid points to store? - Number of grid intervals
643  // plus 1.
644  const Index npi = (Index)npf + 1;
645 
646  // What is the actual grid spacing inside the band?
647  const Numeric gs = bw / npf;
648 
649  // Create the grid for this band:
650  Vector grid(fmin[i], npi, gs);
651 
652  out3 << " Band range " << i << ": " << grid << "\n";
653 
654  // Append to f_grid_array:
655  f_grid_array.reserve(f_grid_array.nelem() + npi);
656  for (Index s = 0; s < npi; ++s) f_grid_array.push_back(grid[s]);
657  }
658 
659  // Copy result to output vector:
660  f_grid = f_grid_array;
661 
662  out2 << " Total number of frequencies in f_grid: " << f_grid.nelem() << "\n";
663 }
664 
665 
666 
667 /* Workspace method: Doxygen documentation will be auto-generated */
668 void f_gridFromSensorHIRS( // WS Output:
669  Vector& f_grid,
670  // WS Input:
671  const Vector& f_backend,
673  // Control Parameters:
674  const Numeric& spacing,
675  const Verbosity& verbosity) {
676  CREATE_OUT2;
677  CREATE_OUT3;
678 
679  // Check input
680  if (spacing <= 0) {
681  ostringstream os;
682  os << "Expected positive spacing. Found spacing to be: " << spacing << "\n";
683  throw runtime_error(os.str());
684  }
685  // Call subfunction to get channel boundaries. Also does input
686  // consistency checking for us.
687  Vector fmin, fmax;
688 
689  // We have to add some additional margin at the band edges,
690  // otherwise the instrument functions are not happy. Define
691  // this in terms of the grid spacing:
692  Numeric delta = 1 * spacing;
693 
695  fmin, fmax, f_backend, backend_channel_response, delta, verbosity);
696 
697  // Ok, now we just have to create a frequency grid for each of the
698  // fmin/fmax ranges.
699 
700  // Create f_grid_array. This is an array of Numeric, so that we
701  // can use the STL push_back function.
702  ArrayOfNumeric f_grid_array;
703 
704  for (Index i = 0; i < fmin.nelem(); ++i) {
705  // Band width:
706  const Numeric bw = fmax[i] - fmin[i];
707 
708  // How many grid intervals do I need?
709  const Numeric npf = ceil(bw / spacing);
710 
711  // How many grid points to store? - Number of grid intervals
712  // plus 1.
713  const Index npi = (Index)npf + 1;
714 
715  // Create the grid for this band:
716  Vector grid;
717  nlinspace(grid, fmin[i], fmax[i], npi);
718 
719  out3 << " Band range " << i << ": " << grid << "\n";
720 
721  // Append to f_grid_array:
722  f_grid_array.reserve(f_grid_array.nelem() + npi);
723  for (Index s = 0; s < npi; ++s) f_grid_array.push_back(grid[s]);
724  }
725 
726  // Copy result to output vector:
727  f_grid = f_grid_array;
728 
729  out2 << " Total number of frequencies in f_grid: " << f_grid.nelem() << "\n";
730 }
731 
732 
733 
734 /* Workspace method: Doxygen documentation will be auto-generated */
736  // WS Output:
737  Vector& f_grid,
738  Vector& f_backend,
741  // WS Input:
742  const Matrix& mm_back, /* met_mm_backend */
743  // Control Parameters:
744  const Vector& freq_spacing,
745  const ArrayOfIndex& freq_number,
746  const Numeric& freq_merge_threshold,
747  const Verbosity&) {
748  // Some sizes
749  const Index nchannels = mm_back.nrows();
750 
751  // Checks of input
752  //
753  chk_met_mm_backend(mm_back);
754  //
755  if (freq_spacing.nelem() != 1 && freq_spacing.nelem() != nchannels)
756  throw std::runtime_error(
757  "Length of *freq_spacing* vector must be either 1 or correspond\n"
758  "to the number of rows in *met_mm_backend*.");
759  //
760  if (freq_number.nelem() != 1 && freq_number.nelem() != nchannels)
761  throw std::runtime_error(
762  "Length of *freq_number* array must be either 1 or correspond\n"
763  "to the number of rows in *met_mm_backend*.");
764  //
765  if (freq_merge_threshold <= 0)
766  throw std::runtime_error("The *freq_merge_threshold* must be > 0.\n");
767  //
768  if (freq_merge_threshold > 100.)
769  throw std::runtime_error(
770  "The *freq_merge_threshold* is only meant to merge frequencies\n"
771  "that are basically identical and only differ slightly due to\n"
772  "numerical inaccuracies. Setting it to >100Hz is not reasonable.");
773 
774  ArrayOfNumeric f_grid_unsorted;
775  ArrayOfIndex nf_per_channel(nchannels);
776  ArrayOfIndex index_in_unsorted;
777 
778  f_backend.resize(nchannels);
779 
780  for (Index ch = 0; ch < nchannels; ch++) {
781  const Numeric lo = mm_back(ch, 0);
782  const Numeric offset1 = mm_back(ch, 1);
783  const Numeric offset2 = mm_back(ch, 2);
784  const Numeric bandwidth = mm_back(ch, 3);
785 
786  const Index this_fnumber =
787  (freq_number.nelem() == 1) ? freq_number[0] : freq_number[ch];
788  const Numeric this_spacing =
789  (freq_spacing.nelem() == 1) ? freq_spacing[0] : freq_spacing[ch];
790 
791  if (this_spacing <= 0)
792  throw std::runtime_error("*freq_spacing must be > 0.");
793 
794  if (this_fnumber == 0) {
795  ostringstream os;
796  os << "*freq_number* must be -1 or greater zero:"
797  << "freq_number[" << ch << "] = " << this_fnumber;
798  std::runtime_error(os.str());
799  }
800 
801  // Number of passbands
802  const Index npassb =
803  1 + ((Index)(offset1 > 0)) + (2 * (Index)(offset2 > 0));
804 
805  // Number of frequencies per passband
806  Index nfperband = this_fnumber;
807  //
808  if (this_fnumber == -1 ||
809  bandwidth / (Numeric)this_fnumber > this_spacing) {
810  nfperband = (Index)ceil(bandwidth / this_spacing);
811  }
812 
813  // Fill the data known so far
814  nf_per_channel[ch] = npassb * nfperband;
815  f_backend[ch] = lo;
816 
817  // Distance between each new f_grid value
818  const Numeric df = bandwidth / (Numeric)nfperband;
819 
820  // Loop passbands and determine f_grid values
821  for (Index b = 0; b < npassb; b++) {
822  // Centre frequency of passband
823  Numeric fc = lo;
824  if (npassb == 2) {
825  fc += (-1 + 2 * (Numeric)b) * offset1;
826  } else if (npassb == 4) {
827  if (b <= 1) {
828  fc -= offset1;
829  } else {
830  fc += offset1;
831  }
832  if (b == 0 || b == 2) {
833  fc -= offset2;
834  } else {
835  fc += offset2;
836  }
837  }
838 
839  // Loop frequencies to add
840  for (Index f_index = 0; f_index < nfperband; f_index++) {
841  const Numeric fnew = fc - bandwidth / 2 + (0.5 + (Numeric)f_index) * df;
842 
843  // Does this frequency exist or not?
844  bool found = false;
845  for (size_t f_try = 0; f_try < f_grid_unsorted.size(); f_try++) {
846  if (abs(fnew - f_grid_unsorted[f_try]) < freq_merge_threshold) {
847  found = true;
848  index_in_unsorted.push_back(f_try);
849  break;
850  }
851  }
852  if (!found) {
853  f_grid_unsorted.push_back(fnew);
854  index_in_unsorted.push_back(f_grid_unsorted.size() - 1);
855  }
856  }
857  }
858  }
859 
860  // Determine sort order for f_grid
861  const size_t nf = f_grid_unsorted.size();
862  ArrayOfIndex move2index(nf);
863 
864  // Create frequency position vector (1...nf)
865  ArrayOfIndex sorted_indices;
866  sorted_indices.resize(nf);
867  for (size_t i = 0; i < nf; i++) {
868  sorted_indices[i] = i;
869  }
870 
871  // Sort frequency position vector by frequency
872  std::sort(sorted_indices.begin(),
873  sorted_indices.end(),
874  CmpArrayOfNumeric(f_grid_unsorted));
875 
876  // Create vector with indices in target vector
877  for (size_t i = 0; i < nf; i++) {
878  move2index[sorted_indices[i]] = i;
879  }
880 
881  // Create f_grid
882  f_grid.resize(nf);
883  //
884  for (size_t f_index = 0; f_index < nf; f_index++) {
885  f_grid[move2index[f_index]] = f_grid_unsorted[f_index];
886  }
887 
888  // Create channel2 fgrid variables
889  channel2fgrid_indexes.resize(nchannels);
890  channel2fgrid_weights.resize(nchannels);
891  Index i = 0;
892  for (Index ch = 0; ch < nchannels; ch++) {
893  channel2fgrid_indexes[ch].resize(nf_per_channel[ch]);
894  channel2fgrid_weights[ch].resize(nf_per_channel[ch]);
895  //
896  for (Index j = 0; j < nf_per_channel[ch]; j++) {
897  channel2fgrid_indexes[ch][j] = move2index[index_in_unsorted[i]];
898  channel2fgrid_weights[ch][j] = 1 / (Numeric)nf_per_channel[ch];
899  i++;
900  }
901  }
902 }
903 
904 
905 
906 /* Workspace method: Doxygen documentation will be auto-generated */
908  const Numeric& spacing,
909  const Numeric& width,
910  const Index& centre,
911  const Verbosity&) {
912  // Create linear, equidistant grid (same in zenith and azimuth)
913  Vector x;
914  Numeric w;
915  if (centre) {
916  w = spacing * ceil(width / spacing);
917  } else {
918  w = spacing * (0.5 + floor(width / spacing));
919  }
920  linspace(x, -w, w, spacing);
921  //
922  const Index l = x.nelem();
923 
924  Matrix dlos_try(l * l, 2, 0);
925  Index n_in = 0;
926  const Numeric c = width * width;
927 
928  for (Index i = 0; i < l; i++) {
929  const Numeric a = x[i] * x[i];
930 
931  for (Index j = 0; j < l; j++) {
932  if (a + x[j] * x[j] <= c) {
933  dlos_try(n_in, 0) = x[i];
934  dlos_try(n_in, 1) = x[j];
935  n_in++;
936  }
937  }
938  }
939 
940  mblock_dlos_grid = dlos_try(Range(0, n_in), joker);
941 }
942 
943 
944 
945 /* Workspace method: Doxygen documentation will be auto-generated */
947  const Numeric& spacing,
948  const Numeric& za_width,
949  const Numeric& aa_width,
950  const Index& centre,
951  const Verbosity&) {
952  // Create za-grid
953  Vector za;
954  Numeric w;
955  if (centre) {
956  w = spacing * ceil(za_width / spacing);
957  } else {
958  w = spacing * (0.5 + floor(za_width / spacing));
959  }
960  linspace(za, -w, w, spacing);
961 
962  // Create aa-grid
963  Vector aa;
964  if (centre) {
965  w = spacing * ceil(aa_width / spacing);
966  } else {
967  w = spacing * (0.5 + floor(aa_width / spacing));
968  }
969  linspace(aa, -w, w, spacing);
970 
971  const Index nza = za.nelem();
972  const Index naa = aa.nelem();
973 
974  mblock_dlos_grid.resize(nza * naa, 2);
975 
976  Index n = 0;
977 
978  for (Index a = 0; a < naa; a++) {
979  for (Index z = 0; z < nza; z++) {
980  mblock_dlos_grid(n, 0) = za[z];
981  mblock_dlos_grid(n, 1) = aa[a];
982  n++;
983  }
984  }
985 }
986 
987 
988 
989 /* Workspace method: Doxygen documentation will be auto-generated */
997  const Index& atmosphere_dim,
998  const Index& antenna_dim,
999  const Matrix& antenna_dlos,
1001  const Index& sensor_norm,
1002  const String& option_2d,
1003  const Verbosity& verbosity) {
1004  CREATE_OUT3;
1005 
1006  // Basic checks
1007  chk_if_in_range("atmosphere_dim", atmosphere_dim, 1, 3);
1008  chk_if_in_range("antenna_dim", antenna_dim, 1, 2);
1009  chk_if_bool("sensor_norm", sensor_norm);
1010 
1011  // Some sizes
1012  const Index nf = sensor_response_f_grid.nelem();
1013  const Index npol = sensor_response_pol_grid.nelem();
1014  const Index nlos = sensor_response_dlos_grid.nrows();
1015  const Index nin = nf * npol * nlos;
1016 
1017  // Initialise a output stream for runtime errors and a flag for errors
1018  ostringstream os;
1019  bool error_found = false;
1020 
1021  // Check that sensor_response variables are consistent in size
1022  if (sensor_response_f.nelem() != nin) {
1023  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
1024  << "grid variables (sensor_response_f_grid etc.).\n";
1025  error_found = true;
1026  }
1027  if (sensor_response.nrows() != nin) {
1028  os << "The sensor block response matrix *sensor_response* does not have\n"
1029  << "right size compared to the sensor grid variables\n"
1030  << "(sensor_response_f_grid etc.).\n";
1031  error_found = true;
1032  }
1033 
1034  // Checks related to antenna dimension
1035  if (antenna_dim == 2 && atmosphere_dim < 3) {
1036  os << "If *antenna_dim* is 2, *atmosphere_dim* must be 3.\n";
1037  error_found = true;
1038  }
1039 
1040  // Basic checks of antenna_dlos
1041  if (antenna_dlos.empty()) throw runtime_error("*antenna_dlos* is empty.");
1042  if (antenna_dlos.ncols() < 1 || antenna_dlos.ncols() > 2)
1043  throw runtime_error("*antenna_dlos* must have one or 2 columns.");
1044  if (atmosphere_dim < 3 && antenna_dlos.ncols() == 2)
1045  throw runtime_error(
1046  "*antenna_dlos* can only have two columns for 3D atmosphers.");
1047 
1048  // We allow angles in antenna_los to be unsorted
1049 
1050  // Checks of antenna_response polarisation dimension
1051  //
1052  const Index lpolgrid =
1053  antenna_response.get_string_grid(GFIELD4_FIELD_NAMES).nelem();
1054  //
1055  if (lpolgrid != 1 && lpolgrid != npol) {
1056  os << "The number of polarisation in *antenna_response* must be 1 or be\n"
1057  << "equal to the number of polarisations used (determined by\n"
1058  << "*stokes_dim* or *instrument_pol*).\n";
1059  error_found = true;
1060  }
1061 
1062  // Checks of antenna_response frequency dimension
1063  //
1064  ConstVectorView aresponse_f_grid =
1065  antenna_response.get_numeric_grid(GFIELD4_F_GRID);
1066  //
1067  chk_if_increasing("f_grid of antenna_response", aresponse_f_grid);
1068  //
1069  Numeric f_dlow = 0.0;
1070  Numeric f_dhigh = 0.0;
1071  //
1072  f_dlow = min(sensor_response_f_grid) - aresponse_f_grid[0];
1073  f_dhigh = last(aresponse_f_grid) - max(sensor_response_f_grid);
1074  //
1075  if (aresponse_f_grid.nelem() > 1) {
1076  if (f_dlow < 0) {
1077  os << "The frequency grid of *antenna_response is too narrow. It must\n"
1078  << "cover all considered frequencies (*f_grid*), if the length\n"
1079  << "is > 1. The grid needs to be expanded with " << -f_dlow
1080  << " Hz in\n"
1081  << "the lower end.\n";
1082  error_found = true;
1083  }
1084  if (f_dhigh < 0) {
1085  os << "The frequency grid of *antenna_response is too narrow. It must\n"
1086  << "cover all considered frequencies (*f_grid*), if the length\n"
1087  << "is > 1. The grid needs to be expanded with " << -f_dhigh
1088  << " Hz in\n"
1089  << "the upper end.\n";
1090  error_found = true;
1091  }
1092  }
1093 
1094  // Checks of antenna_response za dimension
1095  //
1096  ConstVectorView aresponse_za_grid =
1097  antenna_response.get_numeric_grid(GFIELD4_ZA_GRID);
1098  //
1099  chk_if_increasing("za_grid of *antenna_response*", aresponse_za_grid);
1100  //
1101  if (aresponse_za_grid.nelem() < 2) {
1102  os << "The zenith angle grid of *antenna_response* must have >= 2 values.\n";
1103  error_found = true;
1104  }
1105 
1106  // Checks of antenna_response aa dimension
1107  //
1108  ConstVectorView aresponse_aa_grid =
1109  antenna_response.get_numeric_grid(GFIELD4_AA_GRID);
1110  //
1111  if (antenna_dim == 1) {
1112  if (aresponse_aa_grid.nelem() != 1) {
1113  os << "The azimuthal dimension of *antenna_response* must be 1 if\n"
1114  << "*antenna_dim* equals 1.\n";
1115  error_found = true;
1116  }
1117  } else {
1118  chk_if_increasing("aa_grid of antenna_response", aresponse_aa_grid);
1119  //
1120  if (aresponse_za_grid.nelem() < 2) {
1121  os << "The zenith angle grid of *antenna_response* must have >= 2\n"
1122  << "values.\n";
1123  error_found = true;
1124  }
1125  }
1126 
1127  // Check of angular grids. These checks differ with antenna_dim
1128  if (antenna_dim == 1) {
1131  os << "For 1D antennas, the zenith angles in *sensor_response_dlos_grid*\n"
1132  << "must be sorted, either in increasing or decreasing order.\n"
1133  << "The original problem is probably found in *mblock_dlos_grid*.\n";
1134  error_found = true;
1135  }
1136 
1137  else {
1138  // Check if the za relative grid is outside sensor_response_dlos_grid.
1139  Numeric za_dlow = 0.0;
1140  Numeric za_dhigh = 0.0;
1141  //
1142  za_dlow = antenna_dlos(0, 0) + aresponse_za_grid[0] -
1144  za_dhigh = max(sensor_response_dlos_grid(joker, 0)) -
1145  (last(antenna_dlos(joker, 0)) + last(aresponse_za_grid));
1146  //
1147  if (za_dlow < 0) {
1148  os << "The WSV zenith angle part of *sensor_response_dlos_grid* is too narrow.\n"
1149  << "It should be expanded with " << -za_dlow
1150  << " deg in the lower end.\n"
1151  << "This change should be probably applied to *mblock_dlos_grid*.\n";
1152  error_found = true;
1153  }
1154  if (za_dhigh < 0) {
1155  os << "The WSV zenith angle part of *sensor_response_dlos_grid* is too narrow.\n"
1156  << "It should be expanded with " << -za_dhigh
1157  << " deg in the upper end.\n"
1158  << "This change should be probably applied to *mblock_dlos_grid*.\n";
1159  error_found = true;
1160  }
1161  }
1162  } else {
1163  // Demands differs between the options and checks are done inside
1164  // sub-functions
1165  }
1166 
1167  // If errors where found throw runtime_error with the collected error
1168  // message.
1169  if (error_found) throw runtime_error(os.str());
1170 
1171  // And finally check if grids and data size match
1172  antenna_response.checksize_strict();
1173 
1174  // Call the core function
1175  //
1176  Sparse hantenna;
1177  //
1178  if (antenna_dim == 1)
1179  antenna1d_matrix(hantenna,
1180  antenna_dim,
1181  antenna_dlos(joker, 0),
1185  npol,
1186  sensor_norm);
1187  else {
1188 
1189  if (option_2d == "interp_response" ) {
1190  antenna2d_interp_response(hantenna,
1191  antenna_dim,
1192  antenna_dlos,
1196  npol);
1197  }
1198  else if (option_2d == "gridded_dlos" ) {
1199  antenna2d_gridded_dlos(hantenna,
1200  antenna_dim,
1201  antenna_dlos,
1205  npol);
1206  }
1207 
1208  else {
1209  throw runtime_error( "Unrecognised choice for *option_2d*." );
1210  }
1211  }
1212 
1213  // Here we need a temporary sparse that is copy of the sensor_response
1214  // sparse matrix. We need it since the multiplication function can not
1215  // take the same object as both input and output.
1216  Sparse htmp = sensor_response;
1217  sensor_response.resize(hantenna.nrows(), htmp.ncols());
1218  mult(sensor_response, hantenna, htmp);
1219 
1220  // Some extra output.
1221  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1222  << sensor_response.ncols() << "\n";
1223 
1224  // Update sensor_response_dlos_grid
1226 
1227  // Set aux variables
1234 }
1235 
1236 
1237 
1238 /* Workspace method: Doxygen documentation will be auto-generated */
1247  const Vector& f_backend,
1249  const Index& sensor_norm,
1250  const Verbosity& verbosity) {
1251  CREATE_OUT3;
1252 
1253  // Some sizes
1254  const Index nf = sensor_response_f_grid.nelem();
1255  const Index npol = sensor_response_pol_grid.nelem();
1256  const Index nlos = sensor_response_dlos_grid.nrows();
1257  const Index nin = nf * npol * nlos;
1258 
1259  // Initialise an output stream for runtime errors and a flag for errors
1260  ostringstream os;
1261  bool error_found = false;
1262 
1263  // Check that sensor_response variables are consistent in size
1264  if (sensor_response_f.nelem() != nin) {
1265  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
1266  << "grid variables (sensor_response_f_grid etc.).\n";
1267  error_found = true;
1268  }
1269  if (sensor_response.nrows() != nin) {
1270  os << "The sensor block response matrix *sensor_response* does not have\n"
1271  << "right size compared to the sensor grid variables\n"
1272  << "(sensor_response_f_grid etc.).\n";
1273  error_found = true;
1274  }
1275 
1276  // We allow f_backend to be unsorted, but must be inside sensor_response_f_grid
1278  os << "At least one value in *f_backend* (" << min(f_backend)
1279  << ") below range\ncovered by *sensor_response_f_grid* ("
1280  << min(sensor_response_f_grid) << ").\n";
1281  error_found = true;
1282  }
1284  os << "At least one value in *f_backend* (" << max(f_backend)
1285  << ") above range\ncovered by *sensor_response_f_grid* ("
1286  << max(sensor_response_f_grid) << ").\n";
1287  error_found = true;
1288  }
1289 
1290  // Check number of columns in backend_channel_response
1291  //
1292  const Index nrp = backend_channel_response.nelem();
1293  //
1294  if (nrp != 1 && nrp != f_backend.nelem()) {
1295  os << "The WSV *backend_channel_response* must have 1 or n elements,\n"
1296  << "where n is the length of *f_backend*.\n";
1297  error_found = true;
1298  }
1299 
1300  // If errors where found throw runtime_error with the collected error
1301  // message (before error message gets too long).
1302  if (error_found) throw runtime_error(os.str());
1303 
1304  Numeric f_dlow = 0.0;
1305  Numeric f_dhigh = 0.0;
1306 
1307  Index freq_full = nrp > 1;
1308  for (Index i = 0; i < f_backend.nelem(); i++) {
1309  const Index irp = i * freq_full;
1310  ConstVectorView bchr_f_grid =
1311  backend_channel_response[irp].get_numeric_grid(GFIELD1_F_GRID);
1312 
1313  if (bchr_f_grid.nelem() != backend_channel_response[irp].data.nelem()) {
1314  os << "Mismatch in size of grid and data in element " << i
1315  << "\nof *sideband_response*.\n";
1316  error_found = true;
1317  }
1318 
1319  if (!is_increasing(bchr_f_grid)) {
1320  os << "The frequency grid of element " << irp
1321  << " in *backend_channel_response*\nis not strictly increasing.\n";
1322  error_found = true;
1323  }
1324 
1325  // Check if the relative grid added to the channel frequencies expands
1326  // outside sensor_response_f_grid.
1327  //
1328  Numeric f1 = f_backend[i] + bchr_f_grid[0] - min(sensor_response_f_grid);
1329  Numeric f2 =
1330  (max(sensor_response_f_grid) - f_backend[i]) - last(bchr_f_grid);
1331  //
1332  f_dlow = min(f_dlow, f1);
1333  f_dhigh = min(f_dhigh, f2);
1334  }
1335 
1336  if (f_dlow < 0) {
1337  os << "The WSV *sensor_response_f_grid* is too narrow. It should be\n"
1338  << "expanded with " << -f_dlow << " Hz in the lower end. This change\n"
1339  << "should be applied to either *f_grid* or the sensor part in\n"
1340  << "front of *sensor_responseBackend*.\n";
1341  error_found = true;
1342  }
1343  if (f_dhigh < 0) {
1344  os << "The WSV *sensor_response_f_grid* is too narrow. It should be\n"
1345  << "expanded with " << -f_dhigh << " Hz in the higher end. This change\n"
1346  << "should be applied to either *f_grid* or the sensor part in\n"
1347  << "front of *sensor_responseBackend*.\n";
1348  error_found = true;
1349  }
1350 
1351  // If errors where found throw runtime_error with the collected error
1352  // message.
1353  if (error_found) throw runtime_error(os.str());
1354 
1355  // Call the core function
1356  //
1357  Sparse hbackend;
1358  //
1359  spectrometer_matrix(hbackend,
1360  f_backend,
1363  npol,
1364  nlos,
1365  sensor_norm);
1366 
1367  // Here we need a temporary sparse that is copy of the sensor_response
1368  // sparse matrix. We need it since the multiplication function can not
1369  // take the same object as both input and output.
1370  Sparse htmp = sensor_response;
1371  sensor_response.resize(hbackend.nrows(), htmp.ncols());
1372  mult(sensor_response, hbackend, htmp);
1373 
1374  // Some extra output.
1375  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1376  << sensor_response.ncols() << "\n";
1377 
1378  // Update sensor_response_f_grid
1380 
1381  // Set aux variables
1388 }
1389 
1390 
1391 
1392 /* Workspace method: Doxygen documentation will be auto-generated */
1401  const Vector& f_backend,
1403  const Index& sensor_norm,
1404  const Numeric& df1,
1405  const Numeric& df2,
1406  const Verbosity& verbosity) {
1407  // All needed checks are done in sensor_responseBackend
1408 
1410 
1411  // Some needed vectors
1412  Vector f_backend_shifted;
1413  Vector fdummy = sensor_response_f, fdummy_grid = sensor_response_f_grid;
1414 
1415  // Cycle 1
1416  f_backend_shifted = f_backend;
1417  f_backend_shifted += df1;
1418  //
1420  fdummy,
1423  fdummy_grid,
1426  f_backend_shifted,
1428  sensor_norm,
1429  verbosity);
1430  // Cycle 2
1431  f_backend_shifted = f_backend;
1432  f_backend_shifted += df2;
1433  //
1441  f_backend_shifted,
1443  sensor_norm,
1444  verbosity);
1445 
1446  // Total response
1447  sub(sensor_response, H2, H1);
1448 
1449  // sensor_response_f_grid shall be f_backend
1451 
1452  // Set aux variables
1459 }
1460 
1461 
1462 
1463 /* Workspace method: Doxygen documentation will be auto-generated */
1471  const Numeric& w1,
1472  const Numeric& w2,
1473  const Verbosity& verbosity) {
1474  CREATE_OUT3;
1475 
1476  if (sensor_response_dlos_grid.nrows() != 2)
1477  throw runtime_error(
1478  "This method requires that the number of observation directions is 2.");
1479 
1480  if (sensor_response_pol_grid.nelem() != 1)
1481  throw runtime_error(
1482  "This method handles (so far) only single polarisation cases.");
1483 
1484  const Index n = sensor_response_f_grid.nelem();
1485 
1486  // Form H matrix representing beam switching
1487  Sparse Hbswitch(n, 2 * n);
1488  Vector hrow(2 * n, 0.0);
1489  //
1490  for (Index i = 0; i < n; i++) {
1491  hrow[i] = w1;
1492  hrow[i + n] = w2;
1493  //
1494  Hbswitch.insert_row(i, hrow);
1495  //
1496  hrow = 0;
1497  }
1498 
1499  // Here we need a temporary sparse that is copy of the sensor_response
1500  // sparse matrix. We need it since the multiplication function can not
1501  // take the same object as both input and output.
1502  Sparse Htmp = sensor_response;
1503  sensor_response.resize(Hbswitch.nrows(), Htmp.ncols());
1504  mult(sensor_response, Hbswitch, Htmp);
1505 
1506  // Some extra output.
1507  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1508  << sensor_response.ncols() << "\n";
1509 
1510  // Update sensor_response_za_grid
1511  const Vector zaaa = sensor_response_dlos_grid(1, joker);
1512  sensor_response_dlos_grid.resize(1, zaaa.nelem());
1513  sensor_response_dlos_grid(0, joker) = zaaa;
1514 
1515  // Set aux variables
1522 }
1523 
1524 
1525 
1526 /* Workspace method: Doxygen documentation will be auto-generated */
1535  const Verbosity& verbosity) {
1536  CREATE_OUT3;
1537 
1538  if (sensor_response_dlos_grid.nrows() != 1)
1539  throw runtime_error(
1540  "This method requires that the number of observation directions is 1.");
1541 
1542  if (sensor_response_pol_grid.nelem() != 1)
1543  throw runtime_error(
1544  "This method handles (so far) only single polarisation cases.");
1545 
1546  const Index n = sensor_response_f_grid.nelem();
1547  const Index n2 = n / 2;
1548 
1549  if (sensor_response.nrows() != n)
1550  throw runtime_error(
1551  "Assumptions of method are not fulfilled, "
1552  "considering number of rows in *sensor_response* "
1553  "and length of *sensor_response_f_grid*.");
1554 
1555  if (!is_multiple(n, 2))
1556  throw runtime_error(
1557  "There is an odd number of total frequencies, "
1558  "which is not consistent with the assumptions of "
1559  "the method.");
1560 
1561  // Form H matrix representing frequency switching
1562  Sparse Hbswitch(n2, n);
1563  Vector hrow(n, 0.0);
1564  //
1565  for (Index i = 0; i < n2; i++) {
1566  hrow[i] = -1;
1567  hrow[i + n2] = 1;
1568  //
1569  Hbswitch.insert_row(i, hrow);
1570  //
1571  hrow = 0;
1572  }
1573 
1574  // Here we need a temporary sparse that is copy of the sensor_response
1575  // sparse matrix. We need it since the multiplication function can not
1576  // take the same object as both input and output.
1577  Sparse Htmp = sensor_response;
1578  sensor_response.resize(Hbswitch.nrows(), Htmp.ncols());
1579  mult(sensor_response, Hbswitch, Htmp);
1580 
1581  // Some extra output.
1582  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1583  << sensor_response.ncols() << "\n";
1584 
1585  // Update sensor_response_f_grid
1586  const Vector f = sensor_response_f_grid;
1587  sensor_response_f_grid.resize(n2);
1588  sensor_response_f_grid = f[Range(n2, n2)];
1589 
1590  // Set aux variables
1597 }
1598 
1599 
1600 
1601 /* Workspace method: Doxygen documentation will be auto-generated */
1602 void sensor_responseIF2RF( // WS Output:
1605  // WS Input:
1606  const Numeric& lo,
1607  const String& sideband_mode,
1608  const Verbosity&) {
1609  // Check that frequencies are not too high. This might be a floating limit.
1610  // For this we use the variable f_lim, given in Hz.
1611  Numeric f_lim = 30e9;
1612  if (max(sensor_response_f_grid) > f_lim)
1613  throw runtime_error("The frequencies seem to already be given in RF.");
1614 
1615  // Lower band
1616  if (sideband_mode == "lower") {
1617  sensor_response_f *= -1;
1618  sensor_response_f_grid *= -1;
1619  sensor_response_f += lo;
1621  }
1622 
1623  // Upper band
1624  else if (sideband_mode == "upper") {
1625  sensor_response_f += lo;
1627  }
1628 
1629  // Unknown option
1630  else {
1631  throw runtime_error(
1632  "Only allowed options for *sideband _mode* are \"lower\" and \"upper\".");
1633  }
1634 }
1635 
1636 
1637 
1638 /* Workspace method: Doxygen documentation will be auto-generated */
1646  const Index& polyorder,
1647  const Index& nfill,
1648  const Verbosity& verbosity) {
1649  CREATE_OUT3;
1650 
1651  // Some sizes
1652  const Index nf = sensor_response_f_grid.nelem();
1653  const Index npol = sensor_response_pol_grid.nelem();
1654  const Index nlos = sensor_response_dlos_grid.nrows();
1655  const Index nin = nf * npol * nlos;
1656 
1657  // Initialise a output stream for runtime errors and a flag for errors
1658  ostringstream os;
1659  bool error_found = false;
1660 
1661  // Check that sensor_response variables are consistent in size
1662  if (sensor_response_f.nelem() != nin) {
1663  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
1664  << "grid variables (sensor_response_f_grid etc.).\n";
1665  error_found = true;
1666  }
1667  if (sensor_response.nrows() != nin) {
1668  os << "The sensor block response matrix *sensor_response* does not have\n"
1669  << "right size compared to the sensor grid variables\n"
1670  << "(sensor_response_f_grid etc.).\n";
1671  error_found = true;
1672  }
1673 
1674  // Check polyorder and nfill
1675  if (polyorder < 2 || polyorder > 7) {
1676  os << "Accepted range for *polyorder* is [3,7].\n";
1677  error_found = true;
1678  }
1679  if (nfill < 1) {
1680  os << "The argument *nfill* must be > 1.\n";
1681  error_found = true;
1682  }
1683 
1684  // If errors where found throw runtime_error with the collected error
1685  // message.
1686  if (error_found) throw runtime_error(os.str());
1687 
1688  // New frequency grid
1689  //
1690  const Index n1 = nfill + 1;
1691  const Index n2 = nfill + 2;
1692  const Index nnew = (nf - 1) * n1 + 1;
1693  //
1694  Vector fnew(nnew);
1695  //
1696  for (Index i = 0; i < nf - 1; i++) {
1697  Vector fp(n2);
1699  fnew[Range(i * n1, n2)] = fp;
1700  }
1701 
1702  // Find interpolation weights
1703  //
1704  ArrayOfGridPosPoly gp(nnew);
1705  Matrix itw(nnew, polyorder + 1);
1706  //
1707  gridpos_poly(gp, sensor_response_f_grid, fnew, polyorder);
1708  interpweights(itw, gp);
1709 
1710  // Set up H for this part
1711  //
1712  Sparse hpoly(nnew * npol * nlos, nin);
1713  Vector hrow(nin, 0.0);
1714  Index row = 0;
1715  //
1716  for (Index ilos = 0; ilos < nlos; ilos++) {
1717  for (Index iv = 0; iv < nnew; iv++) {
1718  for (Index ip = 0; ip < npol; ip++) {
1719  const Index col0 = ilos * nf * npol;
1720  for (Index i = 0; i < gp[iv].idx.nelem(); i++) {
1721  const Numeric w = gp[iv].w[i];
1722  if (abs(w) > 1e-5) {
1723  hrow[col0 + gp[iv].idx[i] * npol + ip] = w;
1724  }
1725  }
1726  hpoly.insert_row(row, hrow);
1727  for (Index i = 0; i < gp[iv].idx.nelem(); i++) {
1728  hrow[col0 + gp[iv].idx[i] * npol + ip] = 0;
1729  }
1730  row += 1;
1731  }
1732  }
1733  }
1734 
1735  // Here we need a temporary sparse that is copy of the sensor_response
1736  // sparse matrix. We need it since the multiplication function can not
1737  // take the same object as both input and output.
1738  Sparse htmp = sensor_response;
1739  sensor_response.resize(hpoly.nrows(), htmp.ncols());
1740  mult(sensor_response, hpoly, htmp);
1741 
1742  // Some extra output.
1743  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1744  << sensor_response.ncols() << "\n";
1745 
1746  // Update sensor_response_f_grid
1747  sensor_response_f_grid = fnew;
1748 
1749  // Set aux variables
1756 }
1757 
1758 
1759 
1760 /* Workspace method: Doxygen documentation will be auto-generated */
1768  const Vector& f_grid,
1769  const Matrix& mblock_dlos_grid,
1770  const Index& antenna_dim,
1771  const Index& atmosphere_dim,
1772  const Index& stokes_dim,
1773  const Index& sensor_norm,
1774  const Verbosity& verbosity) {
1775  CREATE_OUT2;
1776  CREATE_OUT3;
1777 
1778  // Check input
1779 
1780  // Basic variables
1781  chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
1782  chk_if_in_range("antenna_dim", antenna_dim, 1, 2);
1783  chk_if_bool("sensor_norm", sensor_norm);
1784 
1785  // mblock_dlos_grid
1786  if (mblock_dlos_grid.empty())
1787  throw runtime_error("*mblock_dlos_grid* is empty.");
1788  if (mblock_dlos_grid.ncols() > 2)
1789  throw runtime_error(
1790  "The maximum number of columns in *mblock_dlos_grid* is two.");
1791  if (atmosphere_dim < 3) {
1792  if (mblock_dlos_grid.ncols() != 1)
1793  throw runtime_error(
1794  "For 1D and 2D *mblock_dlos_grid* must have exactly one column.");
1795  if (antenna_dim == 2)
1796  throw runtime_error(
1797  "2D antennas (antenna_dim=2) can only be "
1798  "used with 3D atmospheres.");
1799  }
1800 
1801  // Set grid variables
1804  //
1806  //
1807  for (Index is = 0; is < stokes_dim; is++) {
1808  sensor_response_pol_grid[is] = is + 1;
1809  }
1810 
1811  // Set aux variables
1818 
1819  //Set response matrix to identity matrix
1820  //
1821  const Index n = sensor_response_f.nelem();
1822  //
1823  out2 << " Initialising *sensor_reponse* as a identity matrix.\n";
1824  out3 << " Size of *sensor_response*: " << n << "x" << n << "\n";
1825  //
1826  sensor_response.resize(n, n);
1828 }
1829 
1830 
1831 
1832 /* Workspace method: Doxygen documentation will be auto-generated */
1841  const Index& stokes_dim,
1842  const Vector& f_grid,
1843  const Verbosity& verbosity) {
1844  // Checks are done in sensor_responseInit.
1847 
1848  // Dummy variables (The method is independent of atmosphere_dim.
1849  // atmosphere_dim used below just for some checks when antenna is active, not
1850  // relevant here. ).
1851  const Index sensor_norm = 1, atmosphere_dim = 1;
1852 
1860  f_grid,
1862  antenna_dim,
1864  stokes_dim,
1865  sensor_norm,
1866  verbosity);
1867 }
1868 
1869 
1870 
1871 /* Workspace method: Doxygen documentation will be auto-generated */
1879  const Numeric& lo,
1881  const Index& sensor_norm,
1882  const Verbosity& verbosity) {
1883  CREATE_OUT3;
1884 
1885  // Some sizes
1886  const Index nf = sensor_response_f_grid.nelem();
1887  const Index npol = sensor_response_pol_grid.nelem();
1888  const Index nlos = sensor_response_dlos_grid.nrows();
1889  const Index nin = nf * npol * nlos;
1890 
1891  // Frequency grid of for sideband response specification
1892  ConstVectorView sbresponse_f_grid =
1893  sideband_response.get_numeric_grid(GFIELD1_F_GRID);
1894 
1895  // Initialise a output stream for runtime errors and a flag for errors
1896  ostringstream os;
1897  bool error_found = false;
1898 
1899  // Check that sensor_response variables are consistent in size
1900  if (sensor_response_f.nelem() != nin) {
1901  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
1902  << "grid variables (sensor_response_f_grid etc.).\n";
1903  error_found = true;
1904  }
1905  if (sensor_response.nrows() != nin) {
1906  os << "The sensor block response matrix *sensor_response* does not have\n"
1907  << "right size compared to the sensor grid variables\n"
1908  << "(sensor_response_f_grid etc.).\n";
1909  error_found = true;
1910  }
1911 
1912  // Check that the lo frequency is within the sensor_response_f_grid
1914  os << "The given local oscillator frequency is outside the sensor\n"
1915  << "frequency grid. It must be within the *sensor_response_f_grid*.\n";
1916  error_found = true;
1917  }
1918 
1919  // Checks of sideband_response, partly in combination with lo
1920  if (sbresponse_f_grid.nelem() != sideband_response.data.nelem()) {
1921  os << "Mismatch in size of grid and data in *sideband_response*.\n";
1922  error_found = true;
1923  }
1924  if (sbresponse_f_grid.nelem() < 2) {
1925  os << "At least two data points must be specified in "
1926  << "*sideband_response*.\n";
1927  error_found = true;
1928  }
1929  if (!is_increasing(sbresponse_f_grid)) {
1930  os << "The frequency grid of *sideband_response* must be strictly\n"
1931  << "increasing.\n";
1932  error_found = true;
1933  }
1934  if (fabs(last(sbresponse_f_grid) + sbresponse_f_grid[0]) > 0) {
1935  os << "The end points of the *sideband_response* frequency grid must be\n"
1936  << "symmetrically placed around 0. That is, the grid shall cover a\n"
1937  << "a range that can be written as [-df,df]. \n";
1938  error_found = true;
1939  }
1940 
1941  // Check that response function does not extend outside sensor_response_f_grid
1942  Numeric df_high = lo + last(sbresponse_f_grid) - last(sensor_response_f_grid);
1943  Numeric df_low = sensor_response_f_grid[0] - lo - sbresponse_f_grid[0];
1944  if (df_high > 0 && df_low > 0) {
1945  os << "The *sensor_response_f* grid must be extended by at least\n"
1946  << df_low << " Hz in the lower end and " << df_high << " Hz in the\n"
1947  << "upper end to cover frequency range set by *sideband_response*\n"
1948  << "and *lo*. Or can the frequency grid of *sideband_response* be\n"
1949  << "decreased?";
1950  error_found = true;
1951  } else if (df_high > 0) {
1952  os << "The *sensor_response_f* grid must be extended by at " << df_high
1953  << " Hz\nin the upper end to cover frequency range set by\n"
1954  << "*sideband_response* and *lo*. Or can the frequency grid of\n"
1955  << "*sideband_response* be decreased?";
1956  error_found = true;
1957  } else if (df_low > 0) {
1958  os << "The *sensor_response_f* grid must be extended by at " << df_low
1959  << " Hz\nin the lower end to cover frequency range set by\n"
1960  << "*sideband_response* and *lo*. Or can the frequency grid of\n"
1961  << "*sideband_response* be decreased?";
1962  error_found = true;
1963  }
1964 
1965  // If errors where found throw runtime_error with the collected error
1966  // message.
1967  if (error_found) throw runtime_error(os.str());
1968 
1969  //Call the core function
1970  //
1971  Sparse hmixer;
1972  Vector f_mixer;
1973  //
1974  mixer_matrix(hmixer,
1975  f_mixer,
1976  lo,
1979  npol,
1980  nlos,
1981  sensor_norm);
1982 
1983  // Here we need a temporary sparse that is copy of the sensor_response
1984  // sparse matrix. We need it since the multiplication function can not
1985  // take the same object as both input and output.
1986  Sparse htmp = sensor_response;
1987  sensor_response.resize(hmixer.nrows(), htmp.ncols());
1988  mult(sensor_response, hmixer, htmp);
1989 
1990  // Some extra output.
1991  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
1992  << sensor_response.ncols() << "\n";
1993 
1994  // Update sensor_response_f_grid
1995  sensor_response_f_grid = f_mixer;
1996 
1997  // Set aux variables
2004 }
2005 
2006 
2007 
2008 /* Workspace method: Doxygen documentation will be auto-generated */
2010  // WS Output:
2011  Index& antenna_dim,
2020  Index& sensor_norm,
2021  // WS Input:
2022  const Index& atmosphere_dim,
2023  const Index& stokes_dim,
2024  const Vector& f_grid,
2025  const Vector& f_backend,
2028  const String& iy_unit,
2029  const Matrix& antenna_dlos,
2030  const ArrayOfString& mm_pol, /* met_mm_polarisation */
2031  const Vector& mm_ant, /* met_mm_antenna */
2032  // Control Parameters:
2033  const Index& use_antenna,
2034  const Index& mirror_dza,
2035  const Verbosity& verbosity) {
2036  // Input checks
2037 
2038  chk_if_bool("use_antenna", use_antenna);
2039  chk_if_bool("mirror_dza", mirror_dza);
2040 
2041  if (mm_ant.nelem())
2042  throw std::runtime_error(
2043  "So far inclusion of antenna pattern is NOT supported and\n"
2044  "*met_mm_antenna* must be empty.\n");
2045 
2046  if (!(atmosphere_dim == 1 || atmosphere_dim == 3))
2047  throw std::runtime_error(
2048  "This method only supports 1D and 3D atmospheres.");
2049 
2050  if (antenna_dlos.empty())
2051  throw std::runtime_error("*antenna_dlos* is empty.");
2052 
2053  if (antenna_dlos.ncols() > 2)
2054  throw std::runtime_error(
2055  "The maximum number of columns in *antenna_dlos*"
2056  "is two.");
2057 
2058  // Copy, and possibly extend antenna_dlos
2059  Matrix antenna_dlos_local;
2060 
2061  // Mirror?
2062  if (mirror_dza) {
2063  if (antenna_dlos.ncols() > 1)
2064  throw std::runtime_error(
2065  "With *mirror_dza* set to true, *antenna_dlos*"
2066  "is only allowed to have a single column.");
2067  if (atmosphere_dim != 3)
2068  throw std::runtime_error("*mirror_dza* only makes sense in 3D.");
2069  // We don't want to duplicate zero!
2070  const Index n = antenna_dlos.nrows();
2071  Index nnew = 0;
2072  for (Index i = 0; i < n; i++) {
2073  if (antenna_dlos(i, 0) != 0) {
2074  nnew += 1;
2075  }
2076  }
2077  antenna_dlos_local.resize(n + nnew, 1);
2078  antenna_dlos_local(Range(0, n), 0) = antenna_dlos(joker, 0);
2079  Index pos = n;
2080  for (Index i = n - 1; i >= 0; i--) {
2081  if (antenna_dlos(i, 0) != 0) {
2082  antenna_dlos_local(pos, 0) = -antenna_dlos(i, 0);
2083  pos += 1;
2084  }
2085  }
2086  } else {
2087  antenna_dlos_local = antenna_dlos;
2088  }
2089 
2090  // No normalisation needed here, and set antenna_dim=1 as temporary solution
2091  sensor_norm = 0;
2092  antenna_dim = 1;
2093 
2094  // Create sensor response for mixer+backend, matching one viewing angle
2095  Sparse sensor_response_single;
2096  Matrix mblock_dlos_dummy(1, 1);
2097  mblock_dlos_dummy(0, 0) = 0;
2098  sensor_responseInit(sensor_response_single,
2105  f_grid,
2106  mblock_dlos_dummy,
2107  antenna_dim,
2109  stokes_dim,
2110  sensor_norm,
2111  verbosity);
2112  sensor_responseMixerBackendPrecalcWeights(sensor_response_single,
2119  f_backend,
2122  verbosity);
2123 
2124  // Construct complete sensor_response matrix
2125  const Index num_f = f_grid.nelem();
2126  const Index nchannels = f_backend.nelem();
2127  sensor_response = Sparse(nchannels * antenna_dlos_local.nrows(),
2128  num_f * stokes_dim * antenna_dlos_local.nrows());
2129 
2130  sensor_response_pol_grid.resize(1);
2131  sensor_response_pol_grid[0] = 1;
2132 
2133  if (stokes_dim > 1) {
2134  // With polarisation
2135  if (mm_pol.nelem() != nchannels) {
2136  ostringstream os;
2137  os << "Length of *met_mm_polarisation* (" << mm_pol.nelem()
2138  << ") must match\n"
2139  << "number of channels in *met_mm_backend* (" << nchannels << ").";
2140  throw std::runtime_error(os.str());
2141  }
2142 
2143  Sparse H_pol;
2144  Sparse sensor_response_tmp;
2145 
2146  for (Index iza = 0; iza < antenna_dlos_local.nrows(); iza++) {
2147  sensor_response_tmp = Sparse(nchannels, sensor_response_single.ncols());
2149  H_pol, mm_pol, antenna_dlos_local(iza, 0), stokes_dim, iy_unit);
2150  mult(sensor_response_tmp, H_pol, sensor_response_single);
2151  for (Index r = 0; r < sensor_response_tmp.nrows(); r++)
2152  for (Index c = 0; c < sensor_response_tmp.ncols(); c++) {
2153  const Numeric v = sensor_response_tmp(r, c);
2154 
2155  if (v != 0.)
2156  sensor_response.rw(iza * nchannels + r,
2157  iza * num_f * stokes_dim + c) = v;
2158  }
2159  }
2160  } else {
2161  // No polarisation
2162  for (Index iza = 0; iza < antenna_dlos_local.nrows(); iza++) {
2163  for (Index r = 0; r < sensor_response_single.nrows(); r++)
2164  for (Index c = 0; c < sensor_response_single.ncols(); c++) {
2165  const Numeric v = sensor_response_single(r, c);
2166 
2167  if (v != 0.)
2168  sensor_response.rw(iza * nchannels + r,
2169  iza * num_f * stokes_dim + c) = v;
2170  }
2171  }
2172  }
2173 
2174  antenna_dim = 1;
2175  // Setup antenna
2176  if (use_antenna) {
2177  // FIXME: Do something smart here
2178  throw std::runtime_error("The antenna hasn't arrived yet.");
2179  }
2180 
2181  // mblock angle grids
2182  mblock_dlos_grid = antenna_dlos_local;
2183 
2184  // Set sensor response aux variables
2186 
2187  // Set aux variables
2194 }
2195 
2196 
2197 
2198 /* Workspace method: Doxygen documentation will be auto-generated */
2207  const Vector& f_backend,
2210  const Verbosity& verbosity) {
2211  CREATE_OUT3;
2212 
2213  // Some sizes
2214  const Index nin_f = sensor_response_f_grid.nelem();
2215  const Index nout_f = f_backend.nelem();
2216  const Index npol = sensor_response_pol_grid.nelem();
2217  const Index nlos = sensor_response_dlos_grid.nrows();
2218  const Index nin = nin_f * npol * nlos;
2219  const Index nout = nout_f * npol * nlos;
2220 
2221  // Initialise an output stream for runtime errors and a flag for errors
2222  ostringstream os;
2223  bool error_found = false;
2224 
2225  // Check that sensor_response variables are consistent in size
2226  if (sensor_response_f.nelem() != nin) {
2227  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
2228  << "grid variables (sensor_response_f_grid etc.).\n";
2229  error_found = true;
2230  }
2231  if (sensor_response.nrows() != nin) {
2232  os << "The sensor block response matrix *sensor_response* does not have\n"
2233  << "right size compared to the sensor grid variables\n"
2234  << "(sensor_response_f_grid etc.).\n";
2235  error_found = true;
2236  }
2237 
2238  // We allow f_backend to be unsorted, but must be inside sensor_response_f_grid
2239  if (nout_f == 0) {
2240  os << "*f_backend* is empty !!!\n";
2241  error_found = true;
2242  }
2244  os << "At least one value in *f_backend* (" << min(f_backend)
2245  << ") below range\ncovered by *sensor_response_f_grid* ("
2246  << min(sensor_response_f_grid) << ").\n";
2247  error_found = true;
2248  }
2250  os << "At least one value in *f_backend* (" << max(f_backend)
2251  << ") above range\ncovered by *sensor_response_f_grid* ("
2252  << max(sensor_response_f_grid) << ").\n";
2253  error_found = true;
2254  }
2255 
2256  // frequency index and weights
2257  if (channel2fgrid_indexes.nelem() != nout_f) {
2258  os << "The first size of *channel2fgrid_indexes* an length of *f_backend* "
2259  << "must be equal.\n";
2260  error_found = true;
2261  }
2262  if (channel2fgrid_weights.nelem() != channel2fgrid_indexes.nelem()) {
2263  os << "Leading sizes of *channel2fgrid_indexes* and "
2264  << "*channel2fgrid_weights* differ.\n";
2265  error_found = true;
2266  }
2267  for (Index i = 0; i < nout_f; i++) {
2269  os << "Mismatch in size between *channel2fgrid_indexes* and "
2270  << "*channel2fgrid_weights*, found for array/vector with "
2271  << "index " << i << ".\n";
2272  error_found = true;
2273  }
2274  for (Index j = 0; j < channel2fgrid_indexes[i].nelem(); j++) {
2275  if (channel2fgrid_indexes[i][j] < 0 ||
2276  channel2fgrid_indexes[i][j] >= nin_f) {
2277  os << "At least one value in *channel2fgrid_indexes* is either "
2278  << " < 0 or is too high considering length of "
2279  << "*sensor_response_f_grid*.\n";
2280  error_found = true;
2281  break;
2282  }
2283  }
2284  }
2285 
2286  // If errors where found throw runtime_error with the collected error
2287  if (error_found) throw runtime_error(os.str());
2288 
2289  // Create response matrix
2290  //
2291  Sparse hmb(nout, nin);
2292  {
2293  // Loop output channels
2294  for (Index ifr = 0; ifr < nout_f; ifr++) {
2295  // The summation vector for 1 polarisation and 1 viewing direction
2296  Vector w1(nin_f, 0.0);
2297  for (Index j = 0; j < channel2fgrid_indexes[ifr].nelem(); j++) {
2298  w1[channel2fgrid_indexes[ifr][j]] = channel2fgrid_weights[ifr][j];
2299  }
2300 
2301  // Loop over polarisation and spectra (viewing directions)
2302  // Weights change only with frequency
2303  // (this code is copied from function spectrometer_matrix)
2304  for (Index sp = 0; sp < nlos; sp++) {
2305  for (Index pol = 0; pol < npol; pol++) {
2306  // Distribute the compact weight vector into a complte one
2307  Vector weights_long(nin, 0.0);
2308  weights_long[Range(sp * nin_f * npol + pol, nin_f, npol)] = w1;
2309 
2310  // Insert temp_long into H at the correct row
2311  hmb.insert_row(sp * nout_f * npol + ifr * npol + pol, weights_long);
2312  }
2313  }
2314  }
2315  }
2316 
2317  // Here we need a temporary sparse that is copy of the sensor_response
2318  // sparse matrix. We need it since the multiplication function can not
2319  // take the same object as both input and output.
2320  Sparse htmp = sensor_response;
2321  sensor_response.resize(hmb.nrows(), htmp.ncols());
2322  mult(sensor_response, hmb, htmp);
2323 
2324  // Update sensor_response_f_grid
2326 
2327  // Set aux variables
2334 }
2335 
2336 
2337 
2338 /* Workspace method: Doxygen documentation will be auto-generated */
2347  const Vector& lo_multi,
2352  const Index& sensor_norm,
2353  const Verbosity& verbosity) {
2354  // Some sizes
2355  const Index nf = sensor_response_f_grid.nelem();
2356  const Index npol = sensor_response_pol_grid.nelem();
2357  const Index nlos = sensor_response_dlos_grid.nrows();
2358  const Index nin = nf * npol * nlos;
2359  const Index nlo = lo_multi.nelem();
2360 
2361  // Initialise a output stream for runtime errors and a flag for errors
2362  ostringstream os;
2363  bool error_found = false;
2364 
2365  // Check that sensor_response variables are consistent in size
2366  if (sensor_response_f.nelem() != nin) {
2367  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
2368  << "grid variables (sensor_response_f_grid etc.).\n";
2369  error_found = true;
2370  }
2371  if (sensor_response.nrows() != nin) {
2372  os << "The sensor block response matrix *sensor_response* does not have\n"
2373  << "right size compared to the sensor grid variables\n"
2374  << "(sensor_response_f_grid etc.).\n";
2375  error_found = true;
2376  }
2377 
2378  // Check that response data are consistent with respect to number of
2379  // mixer/reciever chains.
2380  if (sideband_response_multi.nelem() != nlo) {
2381  os << "Inconsistency in length between *lo_mixer* and "
2382  << "*sideband_response_multi*.\n";
2383  error_found = true;
2384  }
2385  if (sideband_mode_multi.nelem() != nlo) {
2386  os << "Inconsistency in length between *lo_mixer* and "
2387  << "*sideband_mode_multi*.\n";
2388  error_found = true;
2389  }
2390  if (f_backend_multi.nelem() != nlo) {
2391  os << "Inconsistency in length between *lo_mixer* and "
2392  << "*f_backend_multi*.\n";
2393  error_found = true;
2394  }
2395  if (backend_channel_response_multi.nelem() != nlo) {
2396  os << "Inconsistency in length between *lo_mixer* and "
2397  << "*backend_channel_response_multi*.\n";
2398  error_found = true;
2399  }
2400 
2401  // If errors where found throw runtime_error with the collected error
2402  // message. Data for each mixer and reciever chain are checked below.
2403  if (error_found) throw runtime_error(os.str());
2404 
2405  // Variables for data to be appended
2406  Array<Sparse> sr;
2407  ArrayOfVector srfgrid;
2408  ArrayOfIndex cumsumf(nlo + 1, 0);
2409 
2410  for (Index ilo = 0; ilo < nlo; ilo++) {
2411  // Copies of variables that will be changed, but must be
2412  // restored for next loop
2413  Sparse sr1 = sensor_response;
2414  Vector srf1 = sensor_response_f;
2416  Matrix srdlos1 = sensor_response_dlos;
2417  Vector srfgrid1 = sensor_response_f_grid;
2418 
2419  // Call single reciever methods. Try/catch for improved error message.
2420  try {
2422  srf1,
2423  srpol1,
2424  srdlos1,
2425  srfgrid1,
2428  lo_multi[ilo],
2430  sensor_norm,
2431  verbosity);
2432 
2434  srf1, srfgrid1, lo_multi[ilo], sideband_mode_multi[ilo], verbosity);
2435 
2437  srf1,
2438  srpol1,
2439  srdlos1,
2440  srfgrid1,
2443  f_backend_multi[ilo],
2445  sensor_norm,
2446  verbosity);
2447  } catch (const runtime_error& e) {
2448  ostringstream os2;
2449  os2 << "Error when dealing with receiver/mixer chain (1-based index) "
2450  << ilo + 1 << ":\n"
2451  << e.what();
2452  throw runtime_error(os2.str());
2453  }
2454 
2455  // Store in temporary arrays
2456  sr.push_back(sr1);
2457  srfgrid.push_back(srfgrid1);
2458  //
2459  cumsumf[ilo + 1] = cumsumf[ilo] + srfgrid1.nelem();
2460  }
2461 
2462  // Append data to create sensor_response_f_grid
2463  //
2464  const Index nfnew = cumsumf[nlo];
2465  sensor_response_f_grid.resize(nfnew);
2466  //
2467  for (Index ilo = 0; ilo < nlo; ilo++) {
2468  for (Index i = 0; i < srfgrid[ilo].nelem(); i++) {
2469  sensor_response_f_grid[cumsumf[ilo] + i] = srfgrid[ilo][i];
2470  }
2471  }
2472 
2473  // Append data to create total sensor_response
2474  //
2475  const Index ncols = sr[0].ncols();
2476  const Index npolnew = sensor_response_pol_grid.nelem();
2477  const Index nfpolnew = nfnew * npolnew;
2478  //
2479  sensor_response.resize(nlos * nfpolnew, ncols);
2480  //
2481  Vector dummy(ncols, 0.0);
2482  //
2483  for (Index ilo = 0; ilo < nlo; ilo++) {
2484  const Index nfpolthis = (cumsumf[ilo + 1] - cumsumf[ilo]) * npolnew;
2485 
2486  assert(sr[ilo].nrows() == nlos * nfpolthis);
2487  assert(sr[ilo].ncols() == ncols);
2488 
2489  for (Index ilos = 0; ilos < nlos; ilos++) {
2490  for (Index i = 0; i < nfpolthis; i++) {
2491  // "Poor mans" transfer of a row from one sparse to another
2492  for (Index ic = 0; ic < ncols; ic++) {
2493  dummy[ic] = sr[ilo](ilos * nfpolthis + i, ic);
2494  }
2495 
2496  sensor_response.insert_row(ilos * nfpolnew + cumsumf[ilo] * npolnew + i,
2497  dummy);
2498  }
2499  }
2500  }
2501 
2502  // Set aux variables
2509 }
2510 
2511 
2512 
2513 /* Workspace method: Doxygen documentation will be auto-generated */
2521  const Index& stokes_dim,
2522  const String& iy_unit,
2524  const Verbosity&) {
2525  // Some sizes
2526  const Index nnew = instrument_pol.nelem();
2527  const Index nf = sensor_response_f_grid.nelem();
2528  const Index npol = sensor_response_pol_grid.nelem();
2529  const Index nlos = sensor_response_dlos_grid.nrows();
2530 
2531  // Initialise an output stream for runtime errors and a flag for errors
2532  ostringstream os;
2533  bool error_found = false;
2534 
2535  Index nfz = nf * nlos;
2536  Index nin = nfz * npol;
2537 
2538  if (sensor_response.nrows() != nin) {
2539  os << "The sensor block response matrix *sensor_response* does not have\n"
2540  << "right size compared to the sensor grid variables\n"
2541  << "(sensor_response_f_grid etc.).\n";
2542  error_found = true;
2543  }
2544 
2545  // Check that sensor_response variables are consistent in size
2546  if (sensor_response_f.nelem() != nin) {
2547  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
2548  << "grid variables (sensor_response_f_grid etc.).\n";
2549  error_found = true;
2550  }
2551  if (npol != stokes_dim) {
2552  os << "Number of input polarisation does not match *stokes_dim*.\n";
2553  error_found = true;
2554  }
2555  if (nnew == 0) {
2556  os << "The WSV *instrument_pol* can not be empty.\n";
2557  error_found = true;
2558  }
2559  // If errors where found throw runtime_error with the collected error
2560  // message (before it gets too long)
2561  if (error_found) throw runtime_error(os.str());
2562 
2563  // Check polarisation data more in detail
2564  for (Index i = 0; i < npol && !error_found; i++) {
2565  if (sensor_response_pol_grid[i] != i + 1) {
2566  os << "The input polarisations must be I, Q, U and V (up to "
2567  << "stokes_dim). It seems that input data are for other "
2568  << "polarisation components.";
2569  error_found = true;
2570  }
2571  }
2572  for (Index i = 0; i < nnew && !error_found; i++) {
2573  if (instrument_pol[i] < 1 || instrument_pol[i] > 10) {
2574  os << "The elements of *instrument_pol* must be inside the range [1,10].\n";
2575  error_found = true;
2576  }
2577  }
2578  // If errors where found throw runtime_error with the collected error
2579  // message (before it gets too long)
2580  if (error_found) throw runtime_error(os.str());
2581 
2582  // If errors where found throw runtime_error with the collected error
2583  // message
2584  if (error_found) throw runtime_error(os.str());
2585 
2586  // Normalisation weight to apply
2587  Numeric w = 0.5;
2588  if (iy_unit == "PlanckBT" || iy_unit == "RJBT") {
2589  w = 1.0;
2590  }
2591 
2592  // Form H matrix representing polarisation response
2593  //
2594  Sparse Hpol(nfz * nnew, nin);
2595  Vector hrow(nin, 0);
2596  Index row = 0;
2597  //
2598  for (Index i = 0; i < nfz; i++) {
2599  Index col = i * npol;
2600  for (Index in = 0; in < nnew; in++) {
2601  /* Old code, matching older version of stokes2pol
2602  Index p = instrument_pol[in] - 1;
2603  //
2604  for( Index iv=0; iv<pv[p].nelem(); iv++ )
2605  { hrow[col+iv] = pv[p][iv]; }
2606  */
2607  stokes2pol(
2608  hrow[Range(col, stokes_dim)], stokes_dim, instrument_pol[in], w);
2609  //
2610  Hpol.insert_row(row, hrow);
2611  //
2612  hrow = 0;
2613  row += 1;
2614  }
2615  }
2616 
2617  // Here we need a temporary sparse that is copy of the sensor_response
2618  // sparse matrix. We need it since the multiplication function can not
2619  // take the same object as both input and output.
2620  Sparse Htmp = sensor_response;
2621  sensor_response.resize(Hpol.nrows(), Htmp.ncols());
2622  mult(sensor_response, Hpol, Htmp);
2623 
2624  // Update sensor_response_pol_grid
2626 
2627  // Set aux variables
2634 }
2635 
2636 
2637 
2638 /* Workspace method: Doxygen documentation will be auto-generated */
2643  const Index& stokes_dim,
2644  const Vector& stokes_rotation,
2645  const Verbosity&) {
2646  // Basic checks
2647  chk_if_in_range("stokes_dim", stokes_dim, 1, 4);
2648 
2649  // Some sizes
2650  const Index nf = sensor_response_f_grid.nelem();
2651  const Index npol = sensor_response_pol_grid.nelem();
2652  const Index nlos = sensor_response_dlos_grid.nrows();
2653  const Index nin = nf * npol * nlos;
2654 
2655  //---------------------------------------------------------------------------
2656  // Initialise a output stream for runtime errors and a flag for errors
2657  ostringstream os;
2658  bool error_found = false;
2659 
2660  // Check that sensor_response variables are consistent in size
2661  if (sensor_response.nrows() != nin) {
2662  os << "The sensor block response matrix *sensor_response* does not have\n"
2663  << "right size compared to the sensor grid variables\n"
2664  << "(sensor_response_f_grid etc.).\n";
2665  error_found = true;
2666  }
2667 
2668  // Check special stuff for this method
2669  if (stokes_dim < 3) {
2670  os << "To perform a rotation of the Stokes coordinate system,\n"
2671  << "*stokes_dim* must be >= 3.\n";
2672  error_found = true;
2673  }
2674  if (stokes_rotation.nelem() != nlos) {
2675  os << "Incorrect number of angles in *stokes_rotation*. The length\n"
2676  << "of this matrix must match *sensor_response_dlos_grid*.\n";
2677  error_found = true;
2678  }
2679  if (npol != stokes_dim) {
2680  os << "Inconsistency detected. The length of *sensor_response_pol_grid*\n"
2681  << "must be equal to *stokes_dim*, and this is not the case.\n";
2682  error_found = true;
2683  }
2684  for (Index is = 0; is < npol; is++) {
2685  if (sensor_response_pol_grid[is] != is + 1) {
2686  os << "For this method, the values in *sensor_response_pol_grid* must\n"
2687  << "be 1,2...stokes_dim. This is not the case, indicating that\n"
2688  << "some previous sensor part has that the data no longer are\n"
2689  << "Stokes vectors.\n";
2690  error_found = true;
2691  break;
2692  }
2693  }
2694 
2695  // If errors where found throw runtime_error with the collected error
2696  // message.
2697  if (error_found) throw runtime_error(os.str());
2698  //---------------------------------------------------------------------------
2699 
2700  // Set up complete the H matrix for applying rotation
2701  //
2702  Sparse H(sensor_response.nrows(), sensor_response.ncols());
2703  {
2704  Sparse Hrot(npol, npol); // Mueller matrix for 1 Stokes vec
2705  Vector row(H.ncols(), 0);
2706  Index irow = 0;
2707  //
2708  for (Index ilos = 0; ilos < nlos; ilos++) {
2709  // Rotation matrix for direction of concern
2710  mueller_rotation(Hrot, npol, stokes_rotation[ilos]);
2711 
2712  for (Index ifr = 0; ifr < nf; ifr++) {
2713  for (Index ip = 0; ip < npol; ip++) {
2714  // Fill relevant part of row with matching (complete) row
2715  // in Hrot, and instert this row in H
2716  for (Index is = 0; is < npol; is++) {
2717  row[irow + is] = Hrot.ro(ip, is);
2718  }
2719  H.insert_row(irow + ip, row);
2720  // Re-zero row.
2721  for (Index is = 0; is < npol; is++) {
2722  row[irow + is] = 0;
2723  }
2724  }
2725  // Update irow, i.e. jump to next frequency
2726  irow += npol;
2727  }
2728  }
2729  }
2730 
2731  // Here we need a temporary sparse that is copy of the sensor_response
2732  // sparse matrix. We need it since the multiplication function can not
2733  // take the same object as both input and output.
2734  Sparse Htmp = sensor_response;
2735  sensor_response.resize(Htmp.nrows(), Htmp.ncols()); //Just in case!
2736  mult(sensor_response, H, Htmp);
2737 }
2738 
2739 
2740 
2741 /* Workspace method: Doxygen documentation will be auto-generated */
2742 void sensor_responseGenericAMSU( // WS Output:
2743  Vector& f_grid,
2744  Index& antenna_dim,
2753  Index& sensor_norm,
2754  // WS Input:
2755  const Index& atmosphere_dim,
2756  const Index& stokes_dim,
2758  // WS Generic Input:
2759  const Numeric& spacing,
2760  const Verbosity& verbosity) {
2761  // Number of instrument channels:
2762  const Index n = sensor_description_amsu.nrows();
2763  const Index m = sensor_description_amsu.ncols();
2764 
2765  // FIXME - Oscar Isoz-090413
2766  // add error checks
2767  //
2768 
2769  // The meaning of the columns in sensor_description_amsu is:
2770  // LO frequency, channel center offset from LO, channel width, second offset from LO (to be extended if needed);
2771  // (All in Hz.)
2772  //
2773  if (5 > sensor_description_amsu.ncols()) {
2774  ostringstream os;
2775  os << "Input variable sensor_description_amsu must have atleast five columns, but it has "
2776  << sensor_description_amsu.ncols() << ".";
2777  throw runtime_error(os.str());
2778  }
2779 
2782  Range(joker), Range(1, m - 2)); // Remember to ignore column 2..
2783  ConstVectorView verbosityVectIn =
2786 
2787  //is there any undefined verbosity values in the vector?
2788  //Set the verbosity to one third of the bandwidth to make sure that the passband flanks does't overlap
2789  const Numeric minRatioVerbosityVsFdiff =
2790  10; // To be used when to passbands are closer then one verbosity value
2791  Index totNumPb = 0;
2792  Vector verbosityVect(n);
2793 
2794  for (Index idx = 0; idx < n; ++idx) {
2795  if ((verbosityVectIn[idx] == 0) || (verbosityVectIn[idx] > width[idx])) {
2796  verbosityVect[idx] = ((Numeric)width[idx]) / 3;
2797  } else {
2798  verbosityVect[idx] = verbosityVectIn[idx];
2799  }
2800  }
2801 
2802  // Create a vector to store the number of passbands (PB) for each channel
2803  ArrayOfIndex numPBpseudo(n); // store values used for calc
2804  ArrayOfIndex numPB(n); // Store the true values
2805  // Find the number of IFs for each channel and calculate the number of passbands
2806  for (Index i = 0; i < n; ++i) {
2807  numPB[i] = 0; // make sure that it is zero
2808  for (Index j = 0; j < (m - 2); ++j) {
2809  if (j != 2) {
2810  if (offset(i, j) > 0) {
2811  numPB[i]++;
2812  }
2813  }
2814  }
2815  numPB[i] = 1 << (int)numPB[i]; // number of passbands= 2^numLO
2816  if (numPB[i] == 1) {
2817  numPBpseudo[i] = 2;
2818  } else {
2819  numPBpseudo[i] = numPB[i];
2820  }
2821 
2822  totNumPb += (int)numPB[i];
2823 
2824  if (numPB[i] > 4) {
2825  ostringstream os;
2826  os << "This function does currently not support more than 4 passbands per channel"
2827  << numPB[i] << ".";
2828  throw runtime_error(os.str());
2829  }
2830  }
2831 
2832  // Find the center frequencies for all sub-channels
2833  // Create one center frequency for each passband
2835  ArrayOfVector f_backend_multi(n); // changed !!!
2836  for (Index i = 0; i < n; ++i) {
2837  // Channel frequencies
2838  Vector& f = f_backend_multi[i];
2839  f.resize(1);
2840  f[0] = lo_multi[i] + 0.0 * width[i]; //(offset(i,0)+offset(i,1));
2841 
2842  //channel response
2843  const Index numVal = 4;
2844  backend_channel_response_multi[i].resize(1);
2846  b_resp.set_name("Backend channel response function");
2847  b_resp.resize(numVal * numPBpseudo[i]);
2848  Vector f_range(numVal * numPBpseudo[i]);
2849  Numeric pbOffset = 0;
2850  b_resp.set_grid_name(0, "Frequency");
2851 
2852  Numeric slope = 0; // 1900;
2853  // To avoid overlapping passbands in the AMSU-A sensor, reduce the passbands of each channel by a few Hz
2854  for (Index pbOffsetIdx = 0; pbOffsetIdx < numPBpseudo[i];
2855  ++pbOffsetIdx) { // Filter response
2856  slope = width[i] / 100;
2857  f_range[pbOffsetIdx * numVal + 0] = -0.5 * width[i] - 0 * slope;
2858  f_range[pbOffsetIdx * numVal + 1] = -0.5 * width[i] + 1 * slope;
2859  f_range[pbOffsetIdx * numVal + 2] = +0.5 * width[i] - 1 * slope;
2860  f_range[pbOffsetIdx * numVal + 3] = +0.5 * width[i] + 0 * slope;
2861 
2862  b_resp.data[pbOffsetIdx * numVal + 0] = 0.0 / (Numeric)numPB[i];
2863  ;
2864  b_resp.data[pbOffsetIdx * numVal + 1] = 1.0 / (Numeric)numPB[i];
2865  b_resp.data[pbOffsetIdx * numVal + 2] = 1.0 / (Numeric)numPB[i];
2866  b_resp.data[pbOffsetIdx * numVal + 3] = 0.0 / (Numeric)numPB[i];
2867 
2868  if (numPB[i] == 1) {
2869  if (pbOffsetIdx == 0) {
2870  pbOffset = -0.0 * width[i];
2871  b_resp.data[pbOffsetIdx * numVal + 0] = 0;
2872  b_resp.data[pbOffsetIdx * numVal + 1] = 1.0 / 1;
2873 
2874  b_resp.data[pbOffsetIdx * numVal + 2] = 1.0 / 1;
2875  b_resp.data[pbOffsetIdx * numVal + 3] = 1.0 / 1;
2876  f_range[pbOffsetIdx * numVal + 0] = -0.5 * width[i] - 2 * slope;
2877  f_range[pbOffsetIdx * numVal + 1] = -0.5 * width[i] - 1 * slope;
2878  f_range[pbOffsetIdx * numVal + 2] = -0.5 * width[i] + 1 * slope;
2879  f_range[pbOffsetIdx * numVal + 3] = -0.5 * width[i] + 2 * slope;
2880  }
2881  if (pbOffsetIdx == 1) {
2882  pbOffset = 0.0 * width[i]; //just a dummy band
2883  b_resp.data[pbOffsetIdx * numVal + 0] = 1.0 / 1;
2884  b_resp.data[pbOffsetIdx * numVal + 1] = 1.0 / 1;
2885  b_resp.data[pbOffsetIdx * numVal + 2] = 1.0 / 1;
2886  b_resp.data[pbOffsetIdx * numVal + 3] = 0;
2887  f_range[pbOffsetIdx * numVal + 0] = +0.5 * width[i] - 3 * slope;
2888  f_range[pbOffsetIdx * numVal + 1] = +0.5 * width[i] - 2 * slope;
2889  f_range[pbOffsetIdx * numVal + 2] = +0.5 * width[i] - 1 * slope;
2890  f_range[pbOffsetIdx * numVal + 3] = +0.5 * width[i] + 0 * slope - 10;
2891  // without the extra '-10' it will fail due to too narrow backend sensor response.
2892  }
2893  } else if (
2894  numPB[i] ==
2895  2) { // move the passband in frequency to the correct frequency
2896  if (pbOffsetIdx == 0) {
2897  pbOffset = -offset(i, 0);
2898  }
2899  if (pbOffsetIdx == 1) {
2900  pbOffset = offset(i, 0);
2901  }
2902  }
2903  if (numPB[i] == 4) {
2904  if (pbOffsetIdx == 0) {
2905  pbOffset = -offset(i, 0) - offset(i, 1);
2906  }
2907  if (pbOffsetIdx == 1) {
2908  pbOffset = -offset(i, 0) + offset(i, 1);
2909  }
2910  if (pbOffsetIdx == 2) {
2911  pbOffset = offset(i, 0) - offset(i, 1);
2912  }
2913  if (pbOffsetIdx == 3) {
2914  pbOffset = offset(i, 0) + offset(i, 1);
2915  }
2916  }
2917  for (Index iii = 0; iii < numVal; ++iii) {
2918  f_range[pbOffsetIdx * numVal + iii] += 1 * pbOffset;
2919  }
2920  }
2921  // Are any passbands overlapping?
2922  for (Index ii = 2; ii < (f_range.nelem() - 2); ++ii) {
2923  if (((b_resp.data[ii - 1] == 1) && (b_resp.data[ii] == 0) &&
2924  (b_resp.data[ii + 1] == 0) && (b_resp.data[ii + 2] == 1))) {
2925  if ((f_range[ii] >= f_range[ii + 1])) // Overlapping passbands
2926  {
2927  if ((f_range[ii + 2] - f_range[ii - 1]) >
2928  verbosityVectIn[i]) // difference in frequency between passbands
2929  {
2930  f_range[ii + 1] = f_range[ii + 2] - verbosityVectIn[i] / 2;
2931  f_range[ii] = f_range[ii - 1] + verbosityVectIn[i] / 2;
2932  } else {
2933  f_range[ii - 1] = (f_range[ii] + f_range[ii + 2]) / 2 -
2934  2 * verbosityVectIn[i] / minRatioVerbosityVsFdiff;
2935  f_range[ii + 1] = (f_range[ii] + f_range[ii + 2]) / 2 +
2936  verbosityVectIn[i] / minRatioVerbosityVsFdiff;
2937  f_range[ii] =
2938  f_range[ii - 1] + verbosityVectIn[i] / minRatioVerbosityVsFdiff;
2939  f_range[ii + 2] =
2940  f_range[ii + 1] + verbosityVectIn[i] / minRatioVerbosityVsFdiff;
2941  }
2942  }
2943  }
2944  }
2945  b_resp.set_grid(0, f_range);
2946  }
2947 
2948  // construct sideband response
2950  for (Index i = 0; i < n; ++i) {
2952  r.set_name("Sideband response function");
2953  r.resize(numPBpseudo[i]);
2954  Vector f(numPBpseudo[i]);
2955  if (numPB[i] == 1) {
2956  r.data[0] = 0.5;
2957  f[0] = -.0 * width[i];
2958  r.data[1] = 0.5;
2959  f[1] = .0 * width[i];
2960  } else if (numPB[i] == 2) {
2961  r.data[0] = 1. / (Numeric)numPB[i];
2962  r.data[1] = 1. / (Numeric)numPB[i];
2963  f[0] = -1 * offset(i, 0) - 0.5 * width[i];
2964  f[1] = +1 * offset(i, 0) + 0.5 * width[i];
2965  } else if (numPB[i] == 4) {
2966  r.data[0] = 1. / (Numeric)numPB[i];
2967  r.data[1] = 1. / (Numeric)numPB[i];
2968  r.data[2] = 1. / (Numeric)numPB[i];
2969  r.data[3] = 1. / (Numeric)numPB[i];
2970  f[0] = -offset(i, 0) - offset(i, 1) - 0.5 * width[i];
2971  ;
2972  f[1] = -offset(i, 0) + offset(i, 1) - 0.5 * width[i];
2973  ;
2974  f[2] = +offset(i, 0) - offset(i, 1) + 0.5 * width[i];
2975  ;
2976  f[3] = +offset(i, 0) + offset(i, 1) + 0.5 * width[i];
2977  ;
2978  }
2979  r.set_grid_name(0, "Frequency");
2980  r.set_grid(0, f);
2981  }
2982 
2983  sensor_norm = 1;
2985  // out
2986  f_grid,
2987  // in
2990  spacing,
2991  verbosityVect,
2992  verbosity);
2993 
2994  // do some final work...
2996 
2998  // out
3006  // in
3007  f_grid,
3009  antenna_dim,
3011  stokes_dim,
3012  sensor_norm,
3013  verbosity);
3014 
3015  Index numLO = lo_multi.nelem();
3016  // Variables for data to be appended
3017  // Based on code from m_sensor->sensor_responseMultiMixerBackend()
3018  Array<Sparse> sr;
3019  ArrayOfVector srfgrid;
3020  ArrayOfIndex cumsumf(numLO + 1, 0);
3021  const Index nlos = sensor_response_dlos_grid.nrows();
3022 
3023  // Do this for all channels ....
3024  for (Index idxLO = 0; idxLO < numLO; idxLO++) {
3025  Sparse sr1 = sensor_response;
3026  Vector srf1 = sensor_response_f;
3028  Matrix srdlos1 = sensor_response_dlos;
3029  Vector srfgrid1 = sensor_response_f_grid;
3030 
3031  sensor_responseBackend( // out
3032  sr1,
3033  srf1,
3034  srpol1,
3035  srdlos1,
3036  srfgrid1,
3037  //in
3040  f_backend_multi[idxLO],
3042  sensor_norm,
3043  verbosity);
3044 
3045  // Store in temporary arrays
3046  sr.push_back(sr1);
3047  srfgrid.push_back(srfgrid1);
3048  //
3049  cumsumf[idxLO + 1] = cumsumf[idxLO] + srfgrid1.nelem();
3050  }
3051 
3052  // Append data to create sensor_response_f_grid
3053  //
3054  const Index nfnew = cumsumf[numLO];
3055  sensor_response_f_grid.resize(nfnew);
3056  //
3057  for (Index ilo = 0; ilo < numLO; ilo++) {
3058  for (Index i = 0; i < srfgrid[ilo].nelem(); i++) {
3059  sensor_response_f_grid[cumsumf[ilo] + i] = srfgrid[ilo][i];
3060  }
3061  }
3062 
3063  const Index ncols = sr[0].ncols();
3064  const Index npolnew = sensor_response_pol_grid.nelem();
3065  const Index nfpolnew = nfnew * npolnew;
3066  //
3067  sensor_response.resize(nlos * nfpolnew, ncols);
3068  //
3069  Vector dummy(ncols, 0.0);
3070  //
3071  for (Index ilo = 0; ilo < numLO; ilo++) {
3072  const Index nfpolthis = (cumsumf[ilo + 1] - cumsumf[ilo]) * npolnew;
3073 
3074  assert(sr[ilo].nrows() == nlos * nfpolthis);
3075  assert(sr[ilo].ncols() == ncols);
3076 
3077  for (Index ilos = 0; ilos < nlos; ilos++) {
3078  for (Index i = 0; i < nfpolthis; i++) {
3079  // "Poor mans" transfer of a row from one sparse to another
3080  for (Index ic = 0; ic < ncols; ic++) {
3081  dummy[ic] = sr[ilo](ilos * nfpolthis + i, ic);
3082  }
3083 
3084  sensor_response.insert_row(ilos * nfpolnew + cumsumf[ilo] * npolnew + i,
3085  dummy);
3086  }
3087  }
3088  }
3089 
3096 }
3097 
3098 
3099 
3100 /* Workspace method: Doxygen documentation will be auto-generated */
3101 void sensor_responseSimpleAMSU( // WS Output:
3102  Vector& f_grid,
3103  Index& antenna_dim,
3112  Index& sensor_norm,
3113  // WS Input:
3114  const Index& atmosphere_dim,
3115  const Index& stokes_dim,
3117  // WS Generic Input:
3118  const Numeric& spacing,
3119  const Verbosity& verbosity) {
3120  // Check that sensor_description_amsu has the right dimension:
3121  if (3 != sensor_description_amsu.ncols()) {
3122  ostringstream os;
3123  os << "Input variable sensor_description_amsu must have three columns, but it has "
3124  << sensor_description_amsu.ncols() << ".";
3125  throw runtime_error(os.str());
3126  }
3127 
3128  // Number of instrument channels:
3129  const Index n = sensor_description_amsu.nrows();
3130 
3131  // The meaning of the columns in sensor_description_amsu is:
3132  // LO frequency, channel center offset from LO, channel width.
3133  // (All in Hz.)
3137 
3138  // Channel frequencies:
3140  for (Index i = 0; i < n; ++i) {
3141  Vector& f = f_backend_multi[i];
3142  f.resize(1);
3143  f[0] = lo_multi[i] + offset[i];
3144  }
3145 
3146  // Construct channel response
3148  for (Index i = 0; i < n; ++i) {
3149  backend_channel_response_multi[i].resize(1);
3151  r.set_name("Backend channel response function");
3152  r.resize(2);
3153 
3154  // Frequency range:
3155  Vector f(2);
3156  f[0] = -0.5 * width[i];
3157  f[1] = +0.5 * width[i];
3158  r.set_grid_name(0, "Frequency");
3159  r.set_grid(0, f);
3160 
3161  // Response:
3162  r.data[0] = 1;
3163  r.data[1] = 1;
3164  }
3165 
3166  // Construct sideband response:
3168  for (Index i = 0; i < n; ++i) {
3170  r.set_name("Sideband response function");
3171  r.resize(2);
3172 
3173  // Frequency range:
3174  Vector f(2);
3175  f[0] = -(offset[i] + 0.5 * width[i]);
3176  f[1] = +(offset[i] + 0.5 * width[i]);
3177  r.set_grid_name(0, "Frequency");
3178  r.set_grid(0, f);
3179 
3180  // Response:
3181  r.data[0] = 0.5;
3182  r.data[1] = 0.5;
3183  }
3184 
3185  // Set sideband mode:
3186  ArrayOfString sideband_mode_multi(n, "upper");
3187 
3188  // We want to automatically normalize the sensor response data, so set sensor_norm to 1:
3189  sensor_norm = 1;
3190 
3191  // Now the rest is just to use some workspace methods:
3192  // ---------------------------------------------------
3193 
3195  lo_multi,
3198  spacing,
3199  verbosity);
3200 
3202 
3210  f_grid,
3212  antenna_dim,
3214  stokes_dim,
3215  sensor_norm,
3216  verbosity);
3217 
3225  lo_multi,
3230  sensor_norm,
3231  verbosity);
3232 }
3233 
3234 
3235 /* Workspace method: Doxygen documentation will be auto-generated */
3236 void WMRFSelectChannels( // WS Output:
3237  Vector& f_grid,
3239  Vector& f_backend,
3240  // WS Input:
3241  const ArrayOfIndex& wmrf_channels,
3242  const Verbosity& verbosity) {
3243  CREATE_OUT2;
3244  CREATE_OUT3;
3245 
3246  // For error messages:
3247  ostringstream os;
3248 
3249  // Some checks of input parameters:
3250 
3251  // wmrf_weights must have same number of rows as f_backend, and same
3252  // number of columns as f_grid.
3253  if ((wmrf_weights.nrows() != f_backend.nelem()) ||
3254  (wmrf_weights.ncols() != f_grid.nelem())) {
3255  os << "The WSV *wmrf_weights* must have same number of rows as\n"
3256  << "*f_backend*, and same number of columns as *f_grid*.\n"
3257  << "wmrf_weights.nrows() = " << wmrf_weights.nrows() << "\n"
3258  << "f_backend.nelem() = " << f_backend.nelem() << "\n"
3259  << "wmrf_weights.ncols() = " << wmrf_weights.ncols() << "\n"
3260  << "f_grid.nelem() = " << f_grid.nelem();
3261  throw runtime_error(os.str());
3262  }
3263 
3264  // wmrf_channels must be strictly increasing (no repetitions).
3265  chk_if_increasing("wmrf_channels", wmrf_channels);
3266 
3267  // All selected channels must be within the original set of
3268  // channels.
3269  if (min(wmrf_channels) < 0) {
3270  os << "Min(wmrf_channels) must be >= 0, but it is " << min(wmrf_channels)
3271  << ".";
3272  }
3273  if (max(wmrf_channels) >= f_backend.nelem()) {
3274  os << "Max(wmrf_channels) must be less than the total number of channels.\n"
3275  << "(We use zero-based indexing!)\n"
3276  << "The actual value you have is " << max(wmrf_channels) << ".";
3277  }
3278 
3279  if (wmrf_channels.nelem() == f_backend.nelem()) {
3280  // No channels to be removed, I can return the original grid.
3281  out2 << " Retaining all channels.\n";
3282  } else {
3283  out2 << " Reducing number of channels from " << f_backend.nelem() << " to "
3284  << wmrf_channels.nelem() << ".\n";
3285  }
3286 
3287  // Now the real work starts:
3288 
3289  // 1. Remove unwanted channels from f_backend:
3291 
3292  // 2. Remove unwanted channels from wmrf_weights. (We also have to
3293  // do something about the frequency dimension of wmrf_weights, but
3294  // we'll do that later.)
3296 
3297  // 3. Identify, which frequencies are still needed, and which are
3298  // now obsolete. We store the still needed frequencies in an
3299  // ArrayOfIndex.
3300 
3301  // Create f_grid_array. We do not store the frequencies themselves,
3302  // but the indices of the frequencies to use.
3303  ArrayOfIndex selection;
3304  // Make sure that selection does not have to be reallocated along
3305  // the way. (This is purely to improve performance a bit.)
3306  selection.reserve(f_grid.nelem());
3307 
3308  // Go through f_grid, and check for each frequency whether it is in
3309  // the set of WMRF frequencies for any of the channels.
3310  assert(wmrf_weights.nrows() == f_backend.nelem());
3311  assert(wmrf_weights.ncols() == f_grid.nelem());
3312  for (Index fi = 0; fi < wmrf_weights.ncols(); ++fi) {
3313  Index i;
3314  for (i = 0; i < wmrf_weights.nrows(); ++i) {
3315  if (wmrf_weights(i, fi) != 0) {
3316  selection.push_back(fi);
3317  break;
3318  }
3319  }
3320  if (i == wmrf_weights.nrows()) {
3321  out3 << " The frequency with index " << fi
3322  << " is not used by any channel.\n";
3323  }
3324  }
3325 
3326  if (selection.nelem() == f_grid.nelem()) {
3327  // No frequencies were removed, I can return the original grid.
3328  out2 << " No unnecessary frequencies, leaving f_grid untouched.\n";
3329  } else if (selection.nelem() == 0) {
3330  throw runtime_error("No frequencies found for any selected channels.\n");
3331  } else {
3332  out2 << " Reducing number of frequency grid points from " << f_grid.nelem()
3333  << " to " << selection.nelem() << ".\n";
3334  }
3335 
3336  // 4. Select the right frequencies in f_grid:
3337  Select(f_grid, f_grid, selection, verbosity);
3338 
3339  // 5. Select the right frequencies in wmrf_weights. This is a bit
3340  // tricky, since Select works on the row dimension. So we have to
3341  // take the transpose.
3342  Sparse wt(wmrf_weights.ncols(), wmrf_weights.nrows());
3343  transpose(wt, wmrf_weights);
3344  Select(wt, wt, selection, verbosity);
3345  wmrf_weights.resize(wt.ncols(), wt.nrows());
3346  transpose(wmrf_weights, wt);
3347 }
3348 
3349 
3350 
3351 /* Workspace method: Doxygen documentation will be auto-generated */
3352 void sensor_responseWMRF( // WS Output:
3358  // WS Input:
3361  const Sparse& wmrf_weights,
3362  const Vector& f_backend,
3363  const Verbosity& verbosity) {
3364  CREATE_OUT3;
3365 
3366  // Some sizes
3367  const Index nf = sensor_response_f_grid.nelem();
3368  const Index npol = sensor_response_pol_grid.nelem();
3369  const Index nlos = sensor_response_dlos_grid.nrows();
3370  const Index nin = nf * npol * nlos;
3371 
3372  // Initialise output stream for runtime errors and a flag for errors
3373  ostringstream os;
3374  bool error_found = false;
3375 
3376  // Check that sensor_response variables are consistent in size
3377  if (sensor_response_f.nelem() != nin) {
3378  os << "Inconsistency in size between *sensor_response_f* and the sensor\n"
3379  << "grid variables (sensor_response_f_grid etc.).\n";
3380  error_found = true;
3381  }
3382  if (sensor_response.nrows() != nin) {
3383  os << "The sensor block response matrix *sensor_response* does not have\n"
3384  << "right size compared to the sensor grid variables\n"
3385  << "(sensor_response_f_grid etc.).\n";
3386  error_found = true;
3387  }
3388 
3389  if (nin == 0) {
3390  os << "One of f_grid, pol_grid, dlos_grid are empty. Sizes are: (" << nf
3391  << ", " << npol << ", " << nlos << ")"
3392  << "\n";
3393  error_found = true;
3394  }
3395 
3396  // Check number of rows in WMRF weight matrix
3397  //
3398  const Index nrw = wmrf_weights.nrows();
3399  //
3400  if (nrw != f_backend.nelem()) {
3401  os << "The WSV *wmrf_weights* must have as many rows\n"
3402  << "as *f_backend* has elements.\n"
3403  << "wmrf_weights.nrows() = " << nrw << "\n"
3404  << "f_backend.nelem() = " << f_backend.nelem() << "\n";
3405  error_found = true;
3406  }
3407 
3408  // Check number of columns in WMRF weight matrix
3409  //
3410  const Index ncw = wmrf_weights.ncols();
3411  //
3412  if (ncw != sensor_response_f_grid.nelem()) {
3413  os << "The WSV *wmrf_weights* must have as many columns\n"
3414  << "as *sensor_response_f_grid* has elements.\n"
3415  << "wmrf_weights.ncols() = " << ncw << "\n"
3416  << "sensor_response_f_grid.nelem() = " << sensor_response_f_grid.nelem()
3417  << "\n";
3418  error_found = true;
3419  }
3420 
3421  // If errors where found throw runtime_error with the collected error
3422  // message (before error message gets too long).
3423  if (error_found) throw runtime_error(os.str());
3424 
3425  // Ok, now the actual work.
3426 
3427  // Here we need a temporary sparse that is copy of the sensor_response
3428  // sparse matrix. We need it since the multiplication function can not
3429  // take the same object as both input and output.
3430  Sparse htmp = sensor_response;
3431  sensor_response.resize(wmrf_weights.nrows(), htmp.ncols());
3433 
3434  // Some extra output.
3435  out3 << " Size of *sensor_response*: " << sensor_response.nrows() << "x"
3436  << sensor_response.ncols() << "\n";
3437 
3438  // Update sensor_response_f_grid
3440 
3441  // Set aux variables
3448 }
3449 
3450 
3451 
3452 /* Workspace method: Doxygen documentation will be auto-generated */
3454  Vector& y_f,
3455  const Matrix& iy,
3456  const Index& stokes_dim,
3457  const Vector& f_grid,
3458  const Numeric& df,
3459  const Verbosity& verbosity) {
3460  // Some dummy values
3461  const Index sensor_norm = 1, atmosphere_dim = 1;
3462 
3463  // Init sensor reponse
3464  //
3470  //
3472 
3480  f_grid,
3482  antenna_dim,
3484  stokes_dim,
3485  sensor_norm,
3486  verbosity);
3487 
3488  // Center position of "channels"
3489  Vector f_backend;
3490  linspace(f_backend, f_grid[0] + df / 2, last(f_grid), df);
3491 
3492  // Create channel response
3495 
3496  // New sensor response
3504  f_backend,
3505  r,
3506  sensor_norm,
3507  verbosity);
3508 
3509  // Some sizes
3510  const Index nf = f_grid.nelem();
3511  const Index n = sensor_response.nrows();
3512 
3513  // Convert iy to a vector
3514  //
3515  Vector iyb(nf * stokes_dim);
3516  //
3517  for (Index is = 0; is < stokes_dim; is++) {
3518  iyb[Range(is, nf, stokes_dim)] = iy(joker, is);
3519  }
3520 
3521  // y and y_f
3522  //
3524  y.resize(n);
3526 }
3527 
3528 
3529 
3530 /* Workspace method: Doxygen documentation will be auto-generated */
3532  Vector& y_f,
3534  Matrix& y_pos,
3535  Matrix& y_los,
3537  Matrix& y_geo,
3538  Matrix& jacobian,
3539  const Index& stokes_dim,
3540  const Index& jacobian_do,
3541  const Matrix& sensor_pos,
3542  const Matrix& sensor_pol,
3543  const Verbosity&) {
3544  // Some sizes
3545  const Index n1 = y.nelem();
3546  const Index nm = sensor_pol.nrows();
3547  const Index nc = sensor_pol.ncols();
3548  const Index n2 = nm * nc;
3549 
3550  // Check consistency of input data
3551  if (y.empty()) throw runtime_error("Input *y* is empty. Use *yCalc*");
3552  if (y_f.nelem() != n1)
3553  throw runtime_error("Lengths of input *y* and *y_f* are inconsistent.");
3554  if (y_pol.nelem() != n1)
3555  throw runtime_error("Lengths of input *y* and *y_pol* are inconsistent.");
3556  if (y_pos.nrows() != n1)
3557  throw runtime_error("Sizes of input *y* and *y_pos* are inconsistent.");
3558  if (y_los.nrows() != n1)
3559  throw runtime_error("Sizes of input *y* and *y_los* are inconsistent.");
3560  if (y_geo.nrows() != n1)
3561  throw runtime_error("Sizes of input *y* and *y_geo* are inconsistent.");
3562  if (jacobian_do) {
3563  if (jacobian.nrows() != n1)
3564  throw runtime_error("Sizes of *y* and *jacobian* are inconsistent.");
3565  }
3566 
3567  // Checks associated with the Stokes vector
3568  if (stokes_dim < 3)
3569  throw runtime_error(
3570  "*stokes_dim* must be >= 3 to correctly apply a "
3571  "polarisation rotation.");
3572  if (n1 < stokes_dim)
3573  throw runtime_error("Length of input *y* smaller than *stokes_dim*.");
3574  for (Index i = 0; i < stokes_dim; i++) {
3575  if (y_pol[i] != i + 1)
3576  throw runtime_error(
3577  "*y* must hold Stokes element values. Data in "
3578  "*y_pol* indicates that this is not the case.");
3579  }
3580 
3581  // Checks of sensor_pos
3582  if (sensor_pos.nrows() != nm)
3583  throw runtime_error(
3584  "Different number of rows in *sensor_pos* and *sensor_pol*.");
3585  if (n2 * stokes_dim != n1)
3586  throw runtime_error(
3587  "Number of columns in *sensor_pol* not consistent with "
3588  "length of *y* and value of *stokes_dim*.");
3589 
3590  // Make copy of all y variables and jacobian
3591  const Vector y1 = y, y_f1 = y_f;
3592  const Matrix y_pos1 = y_pos, y_los1 = y_los, y_geo1 = y_geo;
3593  const ArrayOfIndex y_pol1 = y_pol;
3594  const ArrayOfVector y_aux1 = y_aux;
3595  Matrix jacobian1(0, 0);
3596  if (jacobian_do) {
3597  jacobian1 = jacobian;
3598  }
3599 
3600  // Resize the y variables and jacobian
3601  y.resize(n2);
3602  y_f.resize(n2);
3603  y_pol.resize(n2);
3604  y_pos.resize(n2, y_pos1.ncols());
3605  y_los.resize(n2, y_los1.ncols());
3606  y_geo.resize(n2, y_geo1.ncols());
3607  for (Index a = 0; a < y_aux.nelem(); a++) y_aux[a].resize(n2);
3608  if (jacobian_do) {
3609  jacobian.resize(n2, jacobian1.ncols());
3610  }
3611 
3612  for (Index r = 0; r < nm; r++) {
3613  for (Index c = 0; c < nc; c++) {
3614  const Index iout = r * nc + c;
3615  const Index iin = iout * stokes_dim;
3616 
3617  const Numeric wq = cos(2 * DEG2RAD * sensor_pol(r, c));
3618  const Numeric wu = sin(2 * DEG2RAD * sensor_pol(r, c));
3619 
3620  // Extract radiance for polarisation angle of concern
3621  y[iout] = y1[iin] + wq * y1[iin + 1] + wu * y1[iin + 2];
3622 
3623  // Same operation for jacobian
3624  if (jacobian_do) {
3625  for (Index q = 0; q < jacobian.ncols(); q++)
3626  jacobian(iout, q) = jacobian1(iin, q) + wq * jacobian1(iin + 1, q) +
3627  wu * jacobian1(iin + 2, q);
3628  }
3629 
3630  // Set y_pol
3631  y_pol[iout] = (Index)sensor_pol(r, c);
3632 
3633  // For the rest, copy value matching I of in-data
3634  y_f[iout] = y_f1[iin];
3635  y_pos(iout, joker) = y_pos1(iin, joker);
3636  y_los(iout, joker) = y_los1(iin, joker);
3637  y_geo(iout, joker) = y_geo1(iin, joker);
3638  for (Index a = 0; a < y_aux.nelem(); a++) y_aux[a][iout] = y_aux1[a][iin];
3639  }
3640  }
3641 }
sensor_responseMixer
void sensor_responseMixer(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Numeric &lo, const GriddedField1 &sideband_response, const Index &sensor_norm, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseMixer.
Definition: m_sensor.cc:1872
Matrix
The Matrix class.
Definition: matpackI.h:1193
ARTS::Var::y_pos
Matrix y_pos(Workspace &ws) noexcept
Definition: autoarts.h:7531
ARTS::Var::atmosphere_dim
Index atmosphere_dim(Workspace &ws) noexcept
Definition: autoarts.h:2510
ARTS::Group::Sparse
Sparse Sparse
Definition: autoarts.h:99
ARTS::Var::y_aux
ArrayOfVector y_aux(Workspace &ws) noexcept
Definition: autoarts.h:7423
ARTS::Var::lo_multi
Vector lo_multi(Workspace &ws) noexcept
Definition: autoarts.h:4055
sensor_responseMultiMixerBackend
void sensor_responseMultiMixerBackend(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Vector &lo_multi, const ArrayOfGriddedField1 &sideband_response_multi, const ArrayOfString &sideband_mode_multi, const ArrayOfVector &f_backend_multi, const ArrayOfArrayOfGriddedField1 &backend_channel_response_multi, const Index &sensor_norm, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseMultiMixerBackend.
Definition: m_sensor.cc:2339
backend_channel_responseFlat
void backend_channel_responseFlat(ArrayOfGriddedField1 &r, const Numeric &resolution, const Verbosity &)
WORKSPACE METHOD: backend_channel_responseFlat.
Definition: m_sensor.cc:323
GriddedField::set_grid_name
void set_grid_name(Index i, const String &s)
Set grid name.
Definition: gridded_fields.h:165
f_gridFromSensorAMSUgeneric
void f_gridFromSensorAMSUgeneric(Vector &f_grid, const ArrayOfVector &f_backend_multi, const ArrayOfArrayOfGriddedField1 &backend_channel_response_multi, const Numeric &spacing, const Vector &verbosityVect, const Verbosity &verbosity)
WORKSPACE METHOD: f_gridFromSensorAMSUgeneric.
Definition: m_sensor.cc:516
ARTS::Var::sensor_response_dlos_grid
Matrix sensor_response_dlos_grid(Workspace &ws) noexcept
Definition: autoarts.h:6356
Tensor4::resize
void resize(Index b, Index p, Index r, Index c)
Resize function.
Definition: matpackIV.cc:1064
ySimpleSpectrometer
void ySimpleSpectrometer(Vector &y, Vector &y_f, const Matrix &iy, const Index &stokes_dim, const Vector &f_grid, const Numeric &df, const Verbosity &verbosity)
WORKSPACE METHOD: ySimpleSpectrometer.
Definition: m_sensor.cc:3453
ARTS::Var::antenna_dlos
Matrix antenna_dlos(Workspace &ws) noexcept
Definition: autoarts.h:2354
SPEED_OF_LIGHT
const Numeric SPEED_OF_LIGHT
find_effective_channel_boundaries
void find_effective_channel_boundaries(Vector &fmin, Vector &fmax, const Vector &f_backend, const ArrayOfGriddedField1 &backend_channel_response, const Numeric &delta, const Verbosity &verbosity)
Calculate channel boundaries from instrument response functions.
Definition: sensor.cc:1152
auto_md.h
f_gridFromSensorAMSU
void f_gridFromSensorAMSU(Vector &f_grid, const Vector &lo, const ArrayOfVector &f_backend, const ArrayOfArrayOfGriddedField1 &backend_channel_response, const Numeric &spacing, const Verbosity &verbosity)
WORKSPACE METHOD: f_gridFromSensorAMSU.
Definition: m_sensor.cc:384
gaussian_response_autogrid
void gaussian_response_autogrid(Vector &x, Vector &y, const Numeric &x0, const Numeric &fwhm, const Numeric &xwidth_si, const Numeric &dx_si)
gaussian_response_autogrid
Definition: sensor.cc:606
ARTS::Var::channel2fgrid_weights
ArrayOfVector channel2fgrid_weights(Workspace &ws) noexcept
Definition: autoarts.h:2634
ARTS::Var::instrument_pol
ArrayOfIndex instrument_pol(Workspace &ws) noexcept
Definition: autoarts.h:3616
ARTS::Var::wmrf_weights
Sparse wmrf_weights(Workspace &ws) noexcept
Definition: autoarts.h:7320
ARTS::Var::lo
Numeric lo(Workspace &ws) noexcept
Definition: autoarts.h:4039
Sparse::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackII.cc:69
w
Complex w(Complex z) noexcept
The Faddeeva function.
Definition: linefunctions.cc:42
chk_if_increasing
void chk_if_increasing(const String &x_name, const ArrayOfIndex &x)
chk_if_increasing
Definition: check_input.cc:117
sensor_responseSimpleAMSU
void sensor_responseSimpleAMSU(Vector &f_grid, Index &antenna_dim, Matrix &mblock_dlos_grid, Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, ArrayOfIndex &sensor_response_pol_grid, Matrix &sensor_response_dlos_grid, Index &sensor_norm, const Index &atmosphere_dim, const Index &stokes_dim, const Matrix &sensor_description_amsu, const Numeric &spacing, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseSimpleAMSU.
Definition: m_sensor.cc:3101
gridpos
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
Definition: interpolation.cc:156
mblock_dlos_gridUniformCircular
void mblock_dlos_gridUniformCircular(Matrix &mblock_dlos_grid, const Numeric &spacing, const Numeric &width, const Index &centre, const Verbosity &)
WORKSPACE METHOD: mblock_dlos_gridUniformCircular.
Definition: m_sensor.cc:907
antenna2d_gridded_dlos
void antenna2d_gridded_dlos(Sparse &H, const Index &antenna_dim, ConstMatrixView antenna_dlos, const GriddedField4 &antenna_response, ConstMatrixView mblock_dlos, ConstVectorView f_grid, const Index n_pol)
antenna2d_interp_gridded_dlos
Definition: sensor.cc:199
spectrometer_matrix
void spectrometer_matrix(Sparse &H, ConstVectorView ch_f, const ArrayOfGriddedField1 &ch_response, ConstVectorView sensor_f, const Index &n_pol, const Index &n_sp, const Index &do_norm)
spectrometer_matrix
Definition: sensor.cc:975
sensor_responseMixerBackendPrecalcWeights
void sensor_responseMixerBackendPrecalcWeights(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Vector &f_backend, const ArrayOfArrayOfIndex &channel2fgrid_indexes, const ArrayOfVector &channel2fgrid_weights, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseMixerBackendPrecalcWeights.
Definition: m_sensor.cc:2199
id_mat
void id_mat(MatrixView I)
Identity Matrix.
Definition: lin_alg.cc:2240
chk_if_bool
void chk_if_bool(const String &x_name, const Index &x)
chk_if_bool
Definition: check_input.cc:65
NAT_LOG_2
const Numeric NAT_LOG_2
ARTS::Var::sensor_response
Sparse sensor_response(Workspace &ws) noexcept
Definition: autoarts.h:6300
backend_channel_responseGaussian
void backend_channel_responseGaussian(ArrayOfGriddedField1 &r, const Vector &fwhm, const Vector &xwidth_si, const Vector &dx_si, const Verbosity &)
WORKSPACE METHOD: backend_channel_responseGaussian.
Definition: m_sensor.cc:344
ARTS::Var::backend_channel_response
ArrayOfGriddedField1 backend_channel_response(Workspace &ws) noexcept
Definition: autoarts.h:2548
antenna1d_matrix
void antenna1d_matrix(Sparse &H, const Index &antenna_dim, ConstVectorView antenna_dza, const GriddedField4 &antenna_response, ConstVectorView za_grid, ConstVectorView f_grid, const Index n_pol, const Index do_norm)
antenna1d_matrix
Definition: sensor.cc:58
RAD2DEG
const Numeric RAD2DEG
antenna2d_interp_response
void antenna2d_interp_response(Sparse &H, const Index &antenna_dim, ConstMatrixView antenna_dlos, const GriddedField4 &antenna_response, ConstMatrixView mblock_dlos, ConstVectorView f_grid, const Index n_pol)
antenna2d_interp_response
Definition: sensor.cc:442
joker
const Joker joker
AntennaConstantGaussian1D
void AntennaConstantGaussian1D(Index &antenna_dim, Matrix &mblock_dlos_grid, GriddedField4 &r, Matrix &antenna_dlos, const Index &n_za_grid, const Numeric &fwhm, const Numeric &xwidth_si, const Numeric &dx_si, const Verbosity &verbosity)
WORKSPACE METHOD: AntennaConstantGaussian1D.
Definition: m_sensor.cc:72
GFIELD1_F_GRID
const Index GFIELD1_F_GRID
mblock_dlos_gridUniformRectangular
void mblock_dlos_gridUniformRectangular(Matrix &mblock_dlos_grid, const Numeric &spacing, const Numeric &za_width, const Numeric &aa_width, const Index &centre, const Verbosity &)
WORKSPACE METHOD: mblock_dlos_gridUniformRectangular.
Definition: m_sensor.cc:946
last
Numeric last(ConstVectorView x)
last
Definition: math_funcs.cc:165
interp
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Definition: interpolation.cc:970
ARTS::Var::antenna_response
GriddedField4 antenna_response(Workspace &ws) noexcept
Definition: autoarts.h:2396
sensor_responseFrequencySwitching
void sensor_responseFrequencySwitching(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseFrequencySwitching.
Definition: m_sensor.cc:1527
abs
#define abs(x)
Definition: legacy_continua.cc:20626
ARTS::Var::verbosity
Verbosity verbosity(Workspace &ws) noexcept
Definition: autoarts.h:7112
mult
void mult(ComplexVectorView y, const ConstComplexMatrixView &M, const ConstComplexVectorView &x)
Matrix-Vector Multiplication.
Definition: complex.cc:1579
ARTS::Var::sideband_response
GriddedField1 sideband_response(Workspace &ws) noexcept
Definition: autoarts.h:6497
GFIELD4_ZA_GRID
const Index GFIELD4_ZA_GRID
Sparse
The Sparse class.
Definition: matpackII.h:60
ARTS::Var::f_index
Index f_index(Workspace &ws) noexcept
Definition: autoarts.h:3464
gaussian_response
void gaussian_response(Vector &y, const Vector &x, const Numeric &x0, const Numeric &fwhm)
gaussian_response
Definition: sensor.cc:629
ARTS::Var::jacobian_do
Index jacobian_do(Workspace &ws) noexcept
Definition: autoarts.h:3912
ARTS::Var::y
Vector y(Workspace &ws) noexcept
Definition: autoarts.h:7401
Vector::resize
void resize(Index n)
Resize function.
Definition: matpackI.cc:404
ARTS::Var::y_los
Matrix y_los(Workspace &ws) noexcept
Definition: autoarts.h:7498
ARTS::Var::antenna_dim
Index antenna_dim(Workspace &ws) noexcept
Definition: autoarts.h:2327
data
G0 G2 FVC Y DV Numeric Numeric Numeric Zeeman LowerQuantumNumbers void * data
Definition: arts_api_classes.cc:232
GriddedField1::data
Vector data
Definition: gridded_fields.h:234
ARTS::Var::stokes_dim
Index stokes_dim(Workspace &ws) noexcept
Definition: autoarts.h:6650
ConstMatrixView::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackI.cc:429
Sparse::ro
Numeric ro(Index r, Index c) const
Read only index operator.
Definition: matpackII.cc:130
ARTS::Var::iy_unit
String iy_unit(Workspace &ws) noexcept
Definition: autoarts.h:3856
AntennaMultiBeamsToPencilBeams
void AntennaMultiBeamsToPencilBeams(Matrix &sensor_pos, Matrix &sensor_los, Matrix &antenna_dlos, Index &antenna_dim, Matrix &mblock_dlos_grid, const Index &atmosphere_dim, const Verbosity &verbosity)
WORKSPACE METHOD: AntennaMultiBeamsToPencilBeams.
Definition: m_sensor.cc:118
CREATE_OUT2
#define CREATE_OUT2
Definition: messages.h:206
m_select.h
ARTS::Var::iy
Matrix iy(Workspace &ws) noexcept
Definition: autoarts.h:3690
AntennaOff
void AntennaOff(Index &antenna_dim, Matrix &mblock_dlos_grid, const Verbosity &verbosity)
WORKSPACE METHOD: AntennaOff.
Definition: m_sensor.cc:187
chk_met_mm_backend
void chk_met_mm_backend(const Matrix &mmb)
Check met_mm_backend.
Definition: check_input.cc:1931
ARTS::Var::channel2fgrid_indexes
ArrayOfArrayOfIndex channel2fgrid_indexes(Workspace &ws) noexcept
Definition: autoarts.h:2620
antenna_responseGaussian
void antenna_responseGaussian(GriddedField4 &r, const Numeric &fwhm, const Numeric &xwidth_si, const Numeric &dx_si, const Index &do_2d, const Verbosity &)
WORKSPACE METHOD: antenna_responseGaussian.
Definition: m_sensor.cc:203
sqrt
Numeric sqrt(const Rational r)
Square root.
Definition: rational.h:620
Sparse::insert_row
void insert_row(Index r, Vector v)
Insert row function.
Definition: matpackII.cc:314
sensor_responseFillFgrid
void sensor_responseFillFgrid(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Index &polyorder, const Index &nfill, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseFillFgrid.
Definition: m_sensor.cc:1639
sensor_responseBackendFrequencySwitching
void sensor_responseBackendFrequencySwitching(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Vector &f_backend, const ArrayOfGriddedField1 &backend_channel_response, const Index &sensor_norm, const Numeric &df1, const Numeric &df2, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseBackendFrequencySwitching.
Definition: m_sensor.cc:1393
Array
This can be used to make arrays out of anything.
Definition: array.h:108
pow
Numeric pow(const Rational base, Numeric exp)
Power of.
Definition: rational.h:628
is_decreasing
bool is_decreasing(ConstVectorView x)
Checks if a vector is sorted in reversed order and is strictly decreasing.
Definition: logic.cc:252
Absorption::nelem
Index nelem(const Lines &l)
Number of lines.
Definition: absorptionlines.h:1820
GFIELD4_FIELD_NAMES
const Index GFIELD4_FIELD_NAMES
ARTS::Var::sensor_response_f
Vector sensor_response_f(Workspace &ws) noexcept
Definition: autoarts.h:6376
sensor.h
Sensor modelling functions.
met_mm_polarisation_hmatrix
void met_mm_polarisation_hmatrix(Sparse &H, const ArrayOfString &mm_pol, const Numeric dza, const Index stokes_dim, const String &iy_unit)
Calculate polarisation H-matrix.
Definition: sensor.cc:775
sensor_responseIF2RF
void sensor_responseIF2RF(Vector &sensor_response_f, Vector &sensor_response_f_grid, const Numeric &lo, const String &sideband_mode, const Verbosity &)
WORKSPACE METHOD: sensor_responseIF2RF.
Definition: m_sensor.cc:1602
GriddedField::get_numeric_grid
const Vector & get_numeric_grid(Index i) const
Get a numeric grid.
Definition: gridded_fields.cc:87
min
#define min(a, b)
Definition: legacy_continua.cc:20628
sensorOff
void sensorOff(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, ArrayOfIndex &sensor_response_pol_grid, Matrix &sensor_response_dlos_grid, Matrix &mblock_dlos_grid, const Index &stokes_dim, const Vector &f_grid, const Verbosity &verbosity)
WORKSPACE METHOD: sensorOff.
Definition: m_sensor.cc:1833
CREATE_OUT3
#define CREATE_OUT3
Definition: messages.h:207
messages.h
Declarations having to do with the four output streams.
sensor_responseBackend
void sensor_responseBackend(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Vector &f_backend, const ArrayOfGriddedField1 &backend_channel_response, const Index &sensor_norm, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseBackend.
Definition: m_sensor.cc:1239
ARTS::Var::y_f
Vector y_f(Workspace &ws) noexcept
Definition: autoarts.h:7466
ConstMatrixView::ncols
Index ncols() const
Returns the number of columns.
Definition: matpackI.cc:432
my_basic_string< char >
ARTS::Var::y_pol
ArrayOfIndex y_pol(Workspace &ws) noexcept
Definition: autoarts.h:7515
sensor_responseAntenna
void sensor_responseAntenna(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Matrix &sensor_response_dlos_grid, const Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Index &atmosphere_dim, const Index &antenna_dim, const Matrix &antenna_dlos, const GriddedField4 &antenna_response, const Index &sensor_norm, const String &option_2d, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseAntenna.
Definition: m_sensor.cc:990
ARTS::Var::sensor_pol
Matrix sensor_pol(Workspace &ws) noexcept
Definition: autoarts.h:6253
sub
void sub(Sparse &A, const Sparse &B, const Sparse &C)
Sparse - Sparse subtraction.
Definition: matpackII.cc:667
ConstVectorView::nelem
Index nelem() const
Returns the number of elements.
Definition: matpackI.cc:51
Numeric
NUMERIC Numeric
The type to use for all floating point numbers.
Definition: matpack.h:33
Verbosity
Definition: messages.h:49
ARTS::Var::iyb
Vector iyb(Workspace &ws) noexcept
Definition: autoarts.h:3873
GriddedField4
Definition: gridded_fields.h:340
GriddedField::set_name
void set_name(const String &s)
Set name of this gridded field.
Definition: gridded_fields.h:172
ARTS::Var::f_backend_multi
ArrayOfVector f_backend_multi(Workspace &ws) noexcept
Definition: autoarts.h:3438
math_funcs.h
PI
const Numeric PI
ARTS::Var::f_backend
Vector f_backend(Workspace &ws) noexcept
Definition: autoarts.h:3418
ARTS::Var::f_grid
Vector f_grid(Workspace &ws) noexcept
Definition: autoarts.h:3449
ARTS::Var::sensor_pos
Matrix sensor_pos(Workspace &ws) noexcept
Definition: autoarts.h:6279
antenna_responseVaryingGaussian
void antenna_responseVaryingGaussian(GriddedField4 &r, const Numeric &leff, const Numeric &xwidth_si, const Numeric &dx_si, const Index &nf, const Numeric &fstart, const Numeric &fstop, const Index &do_2d, const Verbosity &verbosity)
WORKSPACE METHOD: antenna_responseVaryingGaussian.
Definition: m_sensor.cc:254
max
#define max(a, b)
Definition: legacy_continua.cc:20629
ARTS::Var::sensor_response_f_grid
Vector sensor_response_f_grid(Workspace &ws) noexcept
Definition: autoarts.h:6392
Matrix::resize
void resize(Index r, Index c)
Resize function.
Definition: matpackI.cc:1056
interpolation_poly.h
Header file for interpolation_poly.cc.
nlinspace
void nlinspace(Vector &x, const Numeric start, const Numeric stop, const Index n)
nlinspace
Definition: math_funcs.cc:231
ConstMatrixView
A constant view of a Matrix.
Definition: matpackI.h:982
ARTS::Var::backend_channel_response_multi
ArrayOfArrayOfGriddedField1 backend_channel_response_multi(Workspace &ws) noexcept
Definition: autoarts.h:2561
DEG2RAD
const Numeric DEG2RAD
ARTS::Var::nrows
Index nrows(Workspace &ws) noexcept
Definition: autoarts.h:4682
oem::Vector
invlib::Vector< ArtsVector > Vector
invlib wrapper type for ARTS vectors.
Definition: oem.h:32
ARTS::Var::mblock_dlos_grid
Matrix mblock_dlos_grid(Workspace &ws) noexcept
Definition: autoarts.h:4247
Range
The range class.
Definition: matpackI.h:160
ARTS::Var::sensor_response_pol
ArrayOfIndex sensor_response_pol(Workspace &ws) noexcept
Definition: autoarts.h:6411
GFIELD4_AA_GRID
const Index GFIELD4_AA_GRID
GriddedField::set_grid
void set_grid(Index i, const Vector &g)
Set a numeric grid.
Definition: gridded_fields.cc:201
ppath.h
Propagation path structure and functions.
mixer_matrix
void mixer_matrix(Sparse &H, Vector &f_mixer, const Numeric &lo, const GriddedField1 &filter, ConstVectorView f_grid, const Index &n_pol, const Index &n_sp, const Index &do_norm)
mixer_matrix
Definition: sensor.cc:645
ARTS::Var::sensor_response_pol_grid
ArrayOfIndex sensor_response_pol_grid(Workspace &ws) noexcept
Definition: autoarts.h:6428
q
#define q
Definition: legacy_continua.cc:21712
is_increasing
bool is_increasing(ConstVectorView x)
Checks if a vector is sorted and strictly increasing.
Definition: logic.cc:215
ARTS::Var::sensor_description_amsu
Matrix sensor_description_amsu(Workspace &ws) noexcept
Definition: autoarts.h:6162
Sparse::nrows
Index nrows() const
Returns the number of rows.
Definition: matpackII.cc:66
ARTS::Var::stokes_rotation
Vector stokes_rotation(Workspace &ws) noexcept
Definition: autoarts.h:6678
sensor_aux_vectors
void sensor_aux_vectors(Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, ConstVectorView sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, ConstMatrixView sensor_response_dlos_grid)
sensor_aux_vectors
Definition: sensor.cc:938
is_multiple
bool is_multiple(const Index &x, const Index &y)
Checks if an integer is a multiple of another integer.
Definition: logic.cc:65
sensor_responsePolarisation
void sensor_responsePolarisation(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, ArrayOfIndex &sensor_response_pol_grid, const Vector &sensor_response_f_grid, const Matrix &sensor_response_dlos_grid, const Index &stokes_dim, const String &iy_unit, const ArrayOfIndex &instrument_pol, const Verbosity &)
WORKSPACE METHOD: sensor_responsePolarisation.
Definition: m_sensor.cc:2514
sensor_responseWMRF
void sensor_responseWMRF(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Sparse &wmrf_weights, const Vector &f_backend, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseWMRF.
Definition: m_sensor.cc:3352
stokes2pol
void stokes2pol(VectorView w, const Index &stokes_dim, const Index &ipol_1based, const Numeric nv)
stokes2pol
Definition: sensor.cc:1041
ARTS::Var::wmrf_channels
ArrayOfIndex wmrf_channels(Workspace &ws) noexcept
Definition: autoarts.h:7306
GriddedField1
Definition: gridded_fields.h:196
VectorNLogSpace
void VectorNLogSpace(Vector &x, const Index &n, const Numeric &start, const Numeric &stop, const Verbosity &verbosity)
WORKSPACE METHOD: VectorNLogSpace.
Definition: m_basic_types.cc:1277
special_interp.h
Header file for special_interp.cc.
GriddedField1::resize
void resize(const GriddedField1 &gf)
Make this GriddedField1 the same size as the given one.
Definition: gridded_fields.h:226
ARTS::Var::sideband_response_multi
ArrayOfGriddedField1 sideband_response_multi(Workspace &ws) noexcept
Definition: autoarts.h:6512
CmpArrayOfNumeric
Helper comparison class to sort an array or vector based on an ArrayOfNumeric.
Definition: array.h:331
sensor_responseBeamSwitching
void sensor_responseBeamSwitching(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Matrix &sensor_response_dlos_grid, const Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Numeric &w1, const Numeric &w2, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseBeamSwitching.
Definition: m_sensor.cc:1464
ARTS::Var::jacobian
Matrix jacobian(Workspace &ws) noexcept
Definition: autoarts.h:3892
ARTS::Var::x
Vector x(Workspace &ws) noexcept
Definition: autoarts.h:7346
sensor_responseStokesRotation
void sensor_responseStokesRotation(Sparse &sensor_response, const Vector &sensor_response_f_grid, const ArrayOfIndex &sensor_response_pol_grid, const Matrix &sensor_response_dlos_grid, const Index &stokes_dim, const Vector &stokes_rotation, const Verbosity &)
WORKSPACE METHOD: sensor_responseStokesRotation.
Definition: m_sensor.cc:2639
Index
INDEX Index
The type to use for all integer numbers and indices.
Definition: matpack.h:39
linspace
std::vector< T > linspace(T s, T e, typename std::vector< T >::size_type count) noexcept
Definition: m_absorptionlines.cc:2539
ARTS::Var::sensor_response_dlos
Matrix sensor_response_dlos(Workspace &ws) noexcept
Definition: autoarts.h:6340
transpose
ConstComplexMatrixView transpose(ConstComplexMatrixView m)
Const version of transpose.
Definition: complex.cc:1509
chk_if_in_range
void chk_if_in_range(const String &x_name, const Index &x, const Index &x_low, const Index &x_high)
chk_if_in_range
Definition: check_input.cc:89
Select
void Select(Array< T > &needles, const Array< T > &haystack, const ArrayOfIndex &needleind, const Verbosity &)
Definition: m_select.h:39
yApplySensorPol
void yApplySensorPol(Vector &y, Vector &y_f, ArrayOfIndex &y_pol, Matrix &y_pos, Matrix &y_los, ArrayOfVector &y_aux, Matrix &y_geo, Matrix &jacobian, const Index &stokes_dim, const Index &jacobian_do, const Matrix &sensor_pos, const Matrix &sensor_pol, const Verbosity &)
WORKSPACE METHOD: yApplySensorPol.
Definition: m_sensor.cc:3531
check_input.h
gridpos_poly
void gridpos_poly(ArrayOfGridPosPoly &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Index order, const Numeric &extpolfac)
The maximum difference from 1 that we allow for a sum check.
Definition: interpolation_poly.cc:120
GriddedField4::data
Tensor4 data
Definition: gridded_fields.h:392
ARTS::Var::ncols
Index ncols(Workspace &ws) noexcept
Definition: autoarts.h:4546
ARTS::Var::sideband_mode_multi
ArrayOfString sideband_mode_multi(Workspace &ws) noexcept
Definition: autoarts.h:6471
Vector
The Vector class.
Definition: matpackI.h:860
sensor_responseMetMM
void sensor_responseMetMM(Index &antenna_dim, Matrix &mblock_dlos_grid, Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, ArrayOfIndex &sensor_response_pol_grid, Matrix &sensor_response_dlos_grid, Index &sensor_norm, const Index &atmosphere_dim, const Index &stokes_dim, const Vector &f_grid, const Vector &f_backend, const ArrayOfArrayOfIndex &channel2fgrid_indexes, const ArrayOfVector &channel2fgrid_weights, const String &iy_unit, const Matrix &antenna_dlos, const ArrayOfString &mm_pol, const Vector &mm_ant, const Index &use_antenna, const Index &mirror_dza, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseMetMM.
Definition: m_sensor.cc:2009
ConstVectorView
A constant view of a Vector.
Definition: matpackI.h:476
f_gridFromSensorHIRS
void f_gridFromSensorHIRS(Vector &f_grid, const Vector &f_backend, const ArrayOfGriddedField1 &backend_channel_response, const Numeric &spacing, const Verbosity &verbosity)
WORKSPACE METHOD: f_gridFromSensorHIRS.
Definition: m_sensor.cc:668
sorting.h
Contains sorting routines.
Array::nelem
Index nelem() const
Number of elements.
Definition: array.h:195
sensor_responseInit
void sensor_responseInit(Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, ArrayOfIndex &sensor_response_pol_grid, Matrix &sensor_response_dlos_grid, const Vector &f_grid, const Matrix &mblock_dlos_grid, const Index &antenna_dim, const Index &atmosphere_dim, const Index &stokes_dim, const Index &sensor_norm, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseInit.
Definition: m_sensor.cc:1761
ARTS::Var::sensor_los
Matrix sensor_los(Workspace &ws) noexcept
Definition: autoarts.h:6201
ARTS::Var::y_geo
Matrix y_geo(Workspace &ws) noexcept
Definition: autoarts.h:7480
GFIELD4_F_GRID
const Index GFIELD4_F_GRID
interpweights
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
Definition: interpolation.cc:741
ARTS::Var::sideband_mode
String sideband_mode(Workspace &ws) noexcept
Definition: autoarts.h:6458
arts.h
The global header file for ARTS.
mueller_rotation
void mueller_rotation(Sparse &H, const Index &stokes_dim, const Numeric &rotangle)
mueller_rotation
Definition: sensor.cc:745
sensor_responseGenericAMSU
void sensor_responseGenericAMSU(Vector &f_grid, Index &antenna_dim, Matrix &mblock_dlos_grid, Sparse &sensor_response, Vector &sensor_response_f, ArrayOfIndex &sensor_response_pol, Matrix &sensor_response_dlos, Vector &sensor_response_f_grid, ArrayOfIndex &sensor_response_pol_grid, Matrix &sensor_response_dlos_grid, Index &sensor_norm, const Index &atmosphere_dim, const Index &stokes_dim, const Matrix &sensor_description_amsu, const Numeric &spacing, const Verbosity &verbosity)
WORKSPACE METHOD: sensor_responseGenericAMSU.
Definition: m_sensor.cc:2742
xml_io.h
This file contains basic functions to handle XML data files.
f_gridMetMM
void f_gridMetMM(Vector &f_grid, Vector &f_backend, ArrayOfArrayOfIndex &channel2fgrid_indexes, ArrayOfVector &channel2fgrid_weights, const Matrix &mm_back, const Vector &freq_spacing, const ArrayOfIndex &freq_number, const Numeric &freq_merge_threshold, const Verbosity &)
WORKSPACE METHOD: f_gridMetMM.
Definition: m_sensor.cc:735
WMRFSelectChannels
void WMRFSelectChannels(Vector &f_grid, Sparse &wmrf_weights, Vector &f_backend, const ArrayOfIndex &wmrf_channels, const Verbosity &verbosity)
WORKSPACE METHOD: WMRFSelectChannels.
Definition: m_sensor.cc:3236
ARTS::Var::sensor_norm
Index sensor_norm(Workspace &ws) noexcept
Definition: autoarts.h:6217