LCOV - code coverage report
Current view: top level - synthesis/MeasurementComponents - Mueller.cc (source / functions) Hit Total Coverage
Test: casacpp_coverage.info Lines: 51 392 13.0 %
Date: 2024-10-04 18:58:15 Functions: 10 48 20.8 %

          Line data    Source code
       1             : //# Mueller.cc: Implementation of Mueller
       2             : //# Copyright (C) 1996,1997,2000,2001,2002,2003
       3             : //# Associated Universities, Inc. Washington DC, USA.
       4             : //#
       5             : //# This library is free software; you can redistribute it and/or modify it
       6             : //# under the terms of the GNU Library General Public License as published by
       7             : //# the Free Software Foundation; either version 2 of the License, or (at your
       8             : //# option) any later version.
       9             : //#
      10             : //# This library is distributed in the hope that it will be useful, but WITHOUT
      11             : //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
      12             : //# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
      13             : //# License for more details.
      14             : //#
      15             : //# You should have received a copy of the GNU Library General Public License
      16             : //# along with this library; if not, write to the Free Software Foundation,
      17             : //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
      18             : //#
      19             : //# Correspondence concerning AIPS++ should be adressed as follows:
      20             : //#        Internet email: casa-feedback@nrao.edu.
      21             : //#        Postal address: AIPS++ Project Office
      22             : //#                        National Radio Astronomy Observatory
      23             : //#                        520 Edgemont Road
      24             : //#                        Charlottesville, VA 22903-2475 USA
      25             : //#
      26             : //#
      27             : 
      28             : #include <synthesis/MeasurementComponents/Mueller.h>
      29             : #include <synthesis/MeasurementComponents/Jones.h>
      30             : #include <synthesis/MeasurementComponents/VisVector.h>
      31             : #include <casacore/casa/aips.h>
      32             : #include <casacore/casa/BasicSL/Complex.h>
      33             : #include <iostream>
      34             : #include <casacore/casa/Exceptions/Error.h>
      35             : 
      36             : #include <casacore/casa/namespace.h>
      37             : 
      38             : using namespace casacore;
      39             : namespace casa { //# NAMESPACE CASA - BEGIN
      40             : 
      41             : // Constructor
      42          16 : Mueller::Mueller() : 
      43          16 :   m0_(NULL), 
      44          16 :   ok0_(NULL),
      45          16 :   m_(NULL), 
      46          16 :   mi_(NULL), 
      47          16 :   ok_(NULL),
      48          16 :   oki_(NULL),
      49          16 :   cOne_(1.0), 
      50          16 :   cZero_(0.0), 
      51          16 :   scalardata_(false),
      52          32 :   vtmp_(VisVector::Four,true)
      53          16 : {}
      54             : 
      55             : // Formation from Jones matrix outer product: General version
      56           0 : void Mueller::fromJones(const Jones& jones1, const Jones& jones2) {
      57             : 
      58           0 :   if (!ok_ || !jones1.ok_ || !jones2.ok_)
      59           0 :     throw(AipsError("Illegal use of Mueller::fromJones."));
      60             : 
      61           0 :   mi_=m_;
      62           0 :   oki_=ok_;
      63             :   Int i,j,io2,jo2,im2,jm2;
      64           0 :   for (i=0;i<4;++i) {
      65           0 :     io2=(i/2)*2;
      66           0 :     im2=(i%2)*2;
      67           0 :     for (j=0;j<4;++j,++mi_,++oki_){
      68           0 :       jo2=(j/2);
      69           0 :       jm2=(j%2);
      70           0 :       (*mi_)  = jones1.j_[io2+jo2]*conj(jones2.j_[im2+jm2]);
      71           0 :       (*oki_) = (jones1.ok_[io2+jo2] && jones2.ok_[im2+jm2]);
      72             :     }
      73             :   }
      74           0 : }
      75             : 
      76             : // In-place multiply onto a VisVector: General version
      77           0 : void Mueller::apply(VisVector& v) {
      78             : 
      79           0 :   switch (v.type()) {
      80           0 :   case VisVector::Four: {
      81           0 :     vtmp_ = v;                          // remember input (necessarily)
      82             :               //  IS VV COPY COMPLETE?
      83           0 :     mi_=m_;                             // iteration within the current matrix
      84             :     Complex *vo,*vi;
      85           0 :     Complex tc;
      86           0 :     vo=v.v_;
      87             :     Int i,j;
      88           0 :     for (j=0;j<4;++j,++vo) {
      89           0 :       vi=vtmp_.v_;
      90           0 :       (*vo) = (*vi);
      91           0 :       (*vo) *= (*mi_);
      92           0 :       ++mi_;
      93           0 :       ++vi;
      94           0 :       for (i=1;i<4;++i,++mi_,++vi) {          // 
      95           0 :         tc = (*vi);
      96           0 :         tc *= (*mi_);
      97           0 :         (*vo) += (tc);   // ((*mi_)*(*vi));
      98             :       }
      99             :     }
     100           0 :     break;
     101             :   }
     102           0 :   default:
     103           0 :     throw(AipsError("Mueller/VisVector type mismatch in apply"));
     104             :     break;
     105             :   }  
     106           0 : }
     107             : 
     108           0 : void Mueller::apply(VisVector& /*v*/, Bool& /*vflag*/) {
     109             : 
     110           0 :   throw(AipsError("Mueller::apply(v,vflag) (general) NYI."));
     111             : 
     112             : }
     113             : 
     114           0 : void Mueller::applyFlag(Bool& /*vflag*/) {
     115           0 :   throw(AipsError("Mueller::applyFlag(vflag) (general) NYI."));
     116             : }
     117             : 
     118           0 : void Mueller::flag(VisVector& /*v*/) {
     119           0 :   throw(AipsError("Mueller::flag(v) (general) NYI."));
     120             : }
     121             : 
     122             : 
     123           0 : void Mueller::invert() {
     124             : 
     125           0 :   throw(AipsError("Invalid attempt to directly invert a general Mueller."));
     126             : 
     127             : }
     128             : 
     129             : // Set matrix elements according to ok flag
     130             : //  (so we don't have to check ok flags atomically in apply)
     131           0 : void Mueller::setMatByOk() {
     132             : 
     133           0 :   throw(AipsError("Invalid attempt to setMatByOk."));
     134             : 
     135             : }
     136             : 
     137             : // Multiply onto a vis VisVector, preserving input (copy then in-place apply)
     138           0 : void Mueller::apply(VisVector& out, const VisVector& in) {
     139           0 :   out=in;
     140           0 :   apply(out);
     141           0 : }
     142             : 
     143           0 : void Mueller::zero() {
     144           0 :   mi_=m_;
     145           0 :   for (Int i=0;i<16;++i,++mi_)
     146           0 :     *mi_=0.0;
     147           0 : }
     148             : 
     149             : // Constructor
     150          16 : MuellerDiag::MuellerDiag() : Mueller() {}
     151             : 
     152             : // Formation from Jones matrix outer product: optimized Diagonal version
     153           0 : void MuellerDiag::fromJones(const Jones& jones1, const Jones& jones2) {
     154             : 
     155           0 :   if (!ok_ || !jones1.ok_ || !jones2.ok_)
     156           0 :     throw(AipsError("Illegal use of MuellerDiag::fromJones."));
     157             : 
     158           0 :   mi_=m_;
     159           0 :   oki_=ok_;
     160           0 :   for (Int i=0;i<4;++i,++mi_) {
     161           0 :     (*mi_)  = jones1.j_[i/2]*conj(jones2.j_[i%2]);
     162           0 :     (*oki_) = (jones1.ok_[i/2]&&jones2.ok_[i%2]);
     163             :   }
     164           0 : }
     165             : 
     166           0 : void MuellerDiag::invert() {
     167             : 
     168           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag::invert."));
     169             : 
     170           0 :   mi_=m_;
     171           0 :   oki_=ok_;
     172           0 :   for (Int i=0;i<4;++i,++mi_,++oki_) 
     173           0 :     if ((*oki_) && abs(*mi_)>0.0)
     174           0 :       (*mi_) = cOne_/(*mi_);
     175             :     else {
     176           0 :       (*mi_)=cZero_;
     177           0 :       (*oki_)=false;
     178             :     }
     179             :   
     180           0 : }
     181             : 
     182             : // Set matrix elements according to ok flag
     183             : //  (so we don't have to check ok flags atomically in apply)
     184           0 : void MuellerDiag::setMatByOk() {
     185             : 
     186             :   // Not needed if ok_ not set
     187           0 :   if (!ok_) return;
     188             : 
     189           0 :   if (!ok_[0]) m_[0]=cOne_;
     190           0 :   if (!ok_[1]) m_[1]=cOne_;
     191           0 :   if (!ok_[2]) m_[2]=cOne_;
     192           0 :   if (!ok_[3]) m_[3]=cOne_;
     193             : 
     194             : }
     195             : 
     196             : 
     197             : 
     198             : // In-place multiply onto a VisVector: optimized Diagonal version
     199           0 : void MuellerDiag::apply(VisVector& v) {
     200             : 
     201             :   // Flag
     202           0 :   if (v.f_) flag(v);
     203             : 
     204           0 :   mi_=m_;
     205           0 :   Complex *vi(v.v_);
     206             : 
     207           0 :   switch (v.type()) {
     208           0 :   case VisVector::Four: {
     209             :     // element-by-element apply of Mueller diagonal to VisVector 4-vector
     210           0 :     for (Int i=0;i<4;++i,++mi_,++vi) (*vi)*=(*mi_);
     211           0 :     break;
     212             :   }
     213           0 :   case VisVector::Two: {
     214             :     // Mueller corner elements apply to VisVector 2-vector
     215           0 :     for (Int i=0;i<2;++i,++vi,mi_+=3) (*vi)*=(*mi_);
     216           0 :     break;
     217             :   }
     218           0 :   case VisVector::One: {
     219             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     220           0 :     (*vi)*=(*mi_);
     221           0 :     break;
     222             :   }
     223             :   }
     224           0 : }
     225             : 
     226           0 : void MuellerDiag::apply(VisVector& v, Bool& vflag) {
     227             : 
     228           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag::apply(v,vflag)."));
     229             : 
     230           0 :   applyFlag(vflag);
     231           0 :   apply(v);
     232             : 
     233           0 : }
     234             : 
     235           0 : void MuellerDiag::applyFlag(Bool& vflag) {
     236             : 
     237           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag::applyFlag(vflag)."));
     238             : 
     239           0 :   if (scalardata_)
     240           0 :     vflag|=(!ok_[0]);
     241             :   else
     242           0 :     vflag|=(!(ok_[0]&&ok_[1]&&ok_[2]&&ok_[3]));
     243           0 : }
     244             : 
     245             : // Flag
     246           0 : void MuellerDiag::flag(VisVector& v) {
     247             : 
     248           0 :   oki_=ok_;
     249           0 :   Bool *fi(v.f_);
     250             : 
     251           0 :   switch (v.type()) {
     252           0 :   case VisVector::Four: {
     253             :     // element-by-element flag of Mueller diagonal to VisVector 4-vector
     254           0 :     for (Int i=0;i<4;++i,++oki_,++fi) (*fi)|=(!(*oki_));
     255           0 :     break;
     256             :   }
     257           0 :   case VisVector::Two: {
     258             :     // Mueller corner elements apply to VisVector 2-vector
     259           0 :     for (Int i=0;i<2;++i,++fi,oki_+=3) (*fi)|=(!(*oki_));
     260           0 :     break;
     261             :   }
     262           0 :   case VisVector::One: {
     263             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     264           0 :     (*fi)|=(!(*oki_));
     265           0 :     break;
     266             :   }
     267             :   }
     268           0 : }
     269             : 
     270             : 
     271           0 : void MuellerDiag::zero() {
     272           0 :   mi_=m_;
     273           0 :   for (Int i=0;i<4;++i,++mi_)
     274           0 :     *mi_=0.0;
     275           0 : }
     276             : 
     277             : // Constructor
     278          16 : MuellerDiag2::MuellerDiag2() : MuellerDiag() {}
     279             :  
     280             : // Formation from Jones matrix outer product: optimized Diag2 version
     281           0 : void MuellerDiag2::fromJones(const Jones& jones1, const Jones& jones2) {
     282             : 
     283           0 :   if (!ok_ || !jones1.ok_ || !jones2.ok_)
     284           0 :     throw(AipsError("Illegal use of MuellerDiag2::fromJones."));
     285             : 
     286           0 :   mi_=m_;
     287           0 :   oki_=ok_;
     288           0 :   for (Int i=0;i<2;++i,++mi_,++oki_) {
     289           0 :     (*mi_) = jones1.j_[i]*conj(jones2.j_[i]);
     290           0 :     (*oki_) = (jones1.ok_[i]&&jones2.ok_[i]);
     291             :   }
     292           0 : }
     293             : 
     294           0 : void MuellerDiag2::invert() {
     295             : 
     296           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag2::invert."));
     297             : 
     298           0 :   mi_=m_;
     299           0 :   oki_=ok_;
     300           0 :   for (Int i=0;i<2;++i,++mi_,++oki_) 
     301           0 :     if ((*oki_)&& abs(*mi_)>0.0)
     302           0 :       (*mi_) = cOne_/(*mi_);
     303             :     else {
     304           0 :       (*mi_) = cZero_;
     305           0 :       (*oki_) = false;
     306             :     }
     307             :   
     308           0 : }
     309             : 
     310             : // Set matrix elements according to ok flag
     311             : //  (so we don't have to check ok flags atomically in apply)
     312           0 : void MuellerDiag2::setMatByOk() {
     313             : 
     314             :   // Not needed if ok_ not set
     315           0 :   if (!ok_) return;
     316             : 
     317           0 :   if (!ok_[0]) m_[0]=cOne_;
     318           0 :   if (!ok_[1]) m_[1]=cOne_;
     319             : 
     320             : }
     321             : 
     322             : 
     323             : // In-place multiply onto a VisVector: optimized Diag2 version
     324           0 : void MuellerDiag2::apply(VisVector& v) {
     325             : 
     326             :   // Flag
     327           0 :   if (v.f_) flag(v);
     328             : 
     329           0 :   switch (v.type()) {
     330           0 :   case VisVector::Four: 
     331             :     // Apply of "corners-only" Mueller to VisVector 4-vector (x-hands zeroed)
     332           0 :     v.v_[0]*=m_[0];
     333           0 :     v.v_[1]*=0.0;
     334           0 :     v.v_[2]*=0.0;
     335           0 :     v.v_[3]*=m_[1];
     336           0 :     break;
     337           0 :   case VisVector::Two: 
     338             :     // Element-by-element apply of "corners-only" Mueller to VisVector 2-vector
     339           0 :     for (Int i=0;i<2;++i) v.v_[i]*=m_[i];
     340           0 :     break;
     341           0 :   case VisVector::One: {
     342             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     343           0 :     v.v_[0]*=(*m_);
     344           0 :     break;
     345             :   }
     346             :   }
     347             :  
     348           0 : }
     349             : 
     350       13975 : void MuellerDiag2::apply(VisVector& v, Bool& vflag) {
     351             : 
     352       13975 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag2::apply(v,vflag)."));
     353             : 
     354       13975 :   applyFlag(vflag);
     355       13975 :   apply(v);
     356             : 
     357       13975 : }
     358             : 
     359       13975 : void MuellerDiag2::applyFlag(Bool& vflag) {
     360             : 
     361       13975 :   if (!ok_) throw(AipsError("Illegal use of MuellerDiag2::applyFlag(vflag)."));
     362             : 
     363       13975 :   if (scalardata_)
     364           0 :     vflag|=(!ok_[0]);
     365             :   else
     366       13975 :     vflag|=(!(ok_[0]&&ok_[1]));
     367       13975 : }
     368             : 
     369             : // Flag a VisVector: optimized Diag2 version
     370           0 : void MuellerDiag2::flag(VisVector& v) {
     371             : 
     372           0 :   switch (v.type()) {
     373           0 :   case VisVector::Four: 
     374             :     // Apply of "corners-only" Mueller to VisVector 4-vector
     375           0 :     v.f_[0]|=(!ok_[0]);
     376           0 :     v.f_[3]|=(!ok_[1]);
     377           0 :     break;
     378           0 :   case VisVector::Two: 
     379             :     // Element-by-element apply of "corners-only" Mueller to VisVector 2-vector
     380           0 :     for (Int i=0;i<2;++i) v.f_[i]|=(!ok_[i]);
     381           0 :     break;
     382           0 :   case VisVector::One: {
     383             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     384           0 :     v.f_[0]|=(!ok_[0]);
     385           0 :     break;
     386             :   }
     387             :   }
     388             :  
     389           0 : }
     390             : 
     391             : 
     392             : 
     393             : 
     394           0 : void MuellerDiag2::zero() {
     395           0 :   mi_=m_;
     396           0 :   for (Int i=0;i<2;++i,++mi_)
     397           0 :     *mi_=0.0;
     398           0 : }
     399             : 
     400             : // Constructor
     401           0 : MuellerScal::MuellerScal() : MuellerDiag() {}
     402             : 
     403             : // Formation from Jones matrix outer product: optimized Scalar version
     404           0 : void MuellerScal::fromJones(const Jones& jones1, const Jones& jones2) {
     405             : 
     406           0 :   if (!ok_ || !jones1.ok_ || !jones2.ok_)
     407           0 :     throw(AipsError("Illegal use of MuellerScal::fromJones."));
     408             : 
     409           0 :   (*m_) = (*jones1.j_)*conj(*jones2.j_);
     410           0 :   (*ok_) = ((*jones1.ok_)&&(*jones2.ok_));
     411           0 : }
     412             : 
     413           0 : void MuellerScal::invert() {
     414             : 
     415           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerScal::invert."));
     416             : 
     417           0 :   if ((*ok_) && abs(*m_)>0.0)
     418           0 :     (*m_) = cOne_/(*m_);
     419             :   else {
     420           0 :     zero();
     421           0 :     (*ok_)=false;
     422             :   }
     423             : 
     424           0 : }
     425             : 
     426             : // Set matrix elements according to ok flag
     427             : //  (so we don't have to check ok flags atomically in apply)
     428           0 : void MuellerScal::setMatByOk() {
     429             : 
     430             :   // Not needed if ok_ not set
     431           0 :   if (!ok_) return;
     432             : 
     433           0 :   if (!ok_[0]) m_[0]=cOne_;
     434             : 
     435             : }
     436             : 
     437             : // In-place multiply onto a VisVector: optimized Scalar version
     438           0 : void MuellerScal::apply(VisVector& v) {
     439             : 
     440             :   // Flag
     441           0 :   if (v.f_) flag(v);
     442             : 
     443             :   // Apply single value to all vector elements
     444           0 :   for (Int i=0;i<v.vistype_;i++) v.v_[i]*=(*m_);
     445           0 : }
     446             : 
     447           0 : void MuellerScal::apply(VisVector& v, Bool& vflag) {
     448             : 
     449           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerScal::apply(v,vflag)."));
     450             : 
     451           0 :   applyFlag(vflag);
     452           0 :   apply(v);
     453             : 
     454           0 : }
     455             : 
     456           0 : void MuellerScal::applyFlag(Bool& vflag) {
     457             : 
     458           0 :   if (!ok_) throw(AipsError("Illegal use of MuellerScal::applyFlag(vflag)."));
     459             : 
     460           0 :   vflag|=(!*ok_);
     461           0 : }
     462             : 
     463             : // Flag a VisVector: optimized Scalar version
     464           0 : void MuellerScal::flag(VisVector& v) {
     465             :   // Apply single value to all vector elements
     466           0 :   Bool f=(!ok_[0]);
     467           0 :   for (Int i=0;i<v.vistype_;i++) v.f_[i]|=f;
     468           0 : }
     469             : 
     470             : 
     471           0 : void MuellerScal::zero() {
     472           0 :     *m_=0.0;
     473           0 : }
     474             : 
     475             : 
     476             : // Constructor
     477          16 : AddMuellerDiag2::AddMuellerDiag2() : MuellerDiag2() {}
     478             :  
     479       14550 : void AddMuellerDiag2::invert() {
     480             : 
     481             :   // "invert" means "negate" for additive term
     482             : 
     483       14550 :   if (!ok_) throw(AipsError("Illegal use of AddMuellerDiag2::invert."));
     484             : 
     485       14550 :   mi_=m_;
     486       14550 :   oki_=ok_;
     487       43650 :   for (Int i=0;i<2;++i,++mi_,++oki_) 
     488       29100 :     if ((*oki_))
     489       29100 :       (*mi_)*=-1.0;  // negate
     490             :     else {
     491           0 :       (*mi_) = cZero_;
     492           0 :       (*oki_) = false;
     493             :     }
     494             :   
     495       14550 : }
     496             : 
     497             : // Set matrix elements according to ok flag
     498             : //  (so we don't have to check ok flags atomically in apply)
     499       14550 : void AddMuellerDiag2::setMatByOk() {
     500             : 
     501             :   // Not needed if ok_ not set
     502       14550 :   if (!ok_) return;
     503             : 
     504       14550 :   if (!ok_[0]) m_[0]=cZero_;
     505       14550 :   if (!ok_[1]) m_[1]=cZero_;
     506             : 
     507             : }
     508             : 
     509             : // In-place multiply onto a VisVector: optimized Diag2 version
     510       13975 : void AddMuellerDiag2::apply(VisVector& v) {
     511             : 
     512             :   // Flag
     513       13975 :   if (v.f_) flag(v);
     514             : 
     515             :   // "apply" means "add" for additive term
     516             : 
     517       13975 :   switch (v.type()) {
     518           0 :   case VisVector::Four: 
     519             :     // Apply of "corners-only" Mueller to VisVector 4-vector (x-hands zeroed)
     520           0 :     v.v_[0]+=m_[0];   // add
     521           0 :     v.v_[1]*=0.0;
     522           0 :     v.v_[2]*=0.0;
     523           0 :     v.v_[3]+=m_[1];   // add
     524           0 :     break;
     525       13975 :   case VisVector::Two: 
     526             :     // Element-by-element apply of "corners-only" Mueller to VisVector 2-vector
     527       41925 :     for (Int i=0;i<2;++i) v.v_[i]+=m_[i]; // add
     528       13975 :     break;
     529           0 :   case VisVector::One: {
     530             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     531           0 :     v.v_[0]+=(*m_);  // add
     532           0 :     break;
     533             :   }
     534             :   }
     535             :  
     536       13975 : }
     537             : 
     538             : // Constructor
     539           0 : AddMuellerDiag::AddMuellerDiag() : MuellerDiag() {}
     540             :  
     541           0 : void AddMuellerDiag::invert() {
     542             : 
     543             :   // "invert" means "negate" for additive term
     544             : 
     545           0 :   if (!ok_) throw(AipsError("Illegal use of AddMuellerDiag2::invert."));
     546             : 
     547           0 :   mi_=m_;
     548           0 :   oki_=ok_;
     549           0 :   for (Int i=0;i<4;++i,++mi_,++oki_) 
     550           0 :     if ((*oki_))
     551           0 :       (*mi_)*=-1.0;  // negate
     552             :     else {
     553           0 :       (*mi_) = cZero_;
     554           0 :       (*oki_) = false;
     555             :     }
     556             :   
     557           0 : }
     558             : 
     559             : // Set matrix elements according to ok flag
     560             : //  (so we don't have to check ok flags atomically in apply)
     561           0 : void AddMuellerDiag::setMatByOk() {
     562             : 
     563             :   // Not needed if ok_ not set
     564           0 :   if (!ok_) return;
     565             : 
     566           0 :   if (!ok_[0]) m_[0]=cZero_;
     567           0 :   if (!ok_[1]) m_[1]=cZero_;
     568           0 :   if (!ok_[2]) m_[2]=cZero_;
     569           0 :   if (!ok_[3]) m_[3]=cZero_;
     570             : 
     571             : }
     572             : 
     573             : 
     574             : // In-place add onto a VisVector: optimized Diag2 version
     575           0 : void AddMuellerDiag::apply(VisVector& v) {
     576             : 
     577             :   // Flag
     578           0 :   if (v.f_) flag(v);
     579             : 
     580             :   // "apply" means "add" for additive term
     581             : 
     582           0 :   switch (v.type()) {
     583           0 :   case VisVector::Four: 
     584             :     // Direct apply of Mueller to VisVector 4-vector 
     585           0 :     v.v_[0]+=m_[0];   // add
     586           0 :     v.v_[1]+=m_[1];
     587           0 :     v.v_[2]+=m_[2];
     588           0 :     v.v_[3]+=m_[3];
     589           0 :     break;
     590           0 :   case VisVector::Two: 
     591             :     // Element-by-element apply of Mueller corners to VisVector 2-vector
     592           0 :     v.v_[0]+=m_[0]; // add
     593           0 :     v.v_[1]+=m_[3]; // add
     594           0 :     break;
     595           0 :   case VisVector::One: {
     596             :     // Mueller corner element used as scalar (pol-sensitivity TBD)
     597           0 :     v.v_[0]+=m_[0];  // add
     598           0 :     break;
     599             :   }
     600             :   }
     601             :  
     602           0 : }
     603             : 
     604             : 
     605             : 
     606             : 
     607             : // GLOBALS---------------------------
     608             : 
     609          16 : Mueller* createMueller(const Mueller::MuellerType& mtype) {
     610          16 :   Mueller *m(NULL);
     611          16 :   switch (mtype) {
     612           0 :   case Mueller::General:
     613           0 :     m = new Mueller();
     614           0 :     break;
     615           0 :   case Mueller::Diagonal:
     616           0 :     m = new MuellerDiag();
     617           0 :     break;
     618           0 :   case Mueller::Diag2:
     619           0 :     m = new MuellerDiag2();
     620           0 :     break;
     621          16 :   case Mueller::AddDiag2:
     622          16 :     m = new AddMuellerDiag2();
     623          16 :     break;
     624           0 :   case Mueller::AddDiag:
     625           0 :     m = new AddMuellerDiag();
     626           0 :     break;
     627           0 :   case Mueller::Scalar:
     628           0 :     m = new MuellerScal();
     629           0 :     break;
     630             :   }
     631          16 :   return m;
     632             : }
     633             : 
     634             : // Return the enum for from an int
     635             : /* 
     636             : Mueller::MuellerType muellerType(const Int& n) {
     637             :   switch (n) {
     638             :   case 1:
     639             :     return Mueller::Scalar;
     640             :   case 2:
     641             :     return Mueller::Diag2;
     642             :   case 4:
     643             :     return Mueller::Diagonal;
     644             :   case 16:
     645             :     return Mueller::General;
     646             :   default:
     647             :     throw(AipsError("Bad MuellerType."));
     648             :   }
     649             : }
     650             : */
     651             : 
     652           0 : Mueller::MuellerType muellerType(const Jones::JonesType& jtype,const VisVector::VisType& vtype) {
     653           0 :   switch(jtype) {
     654           0 :   case Jones::General:
     655           0 :     switch(vtype) {
     656           0 :     case VisVector::Four:
     657           0 :       return Mueller::General;
     658             :       break;
     659           0 :     default:
     660           0 :       throw(AipsError("Cannot apply General calibration matrices to incomplete visibilities."));
     661             :       // TBD: Must cater for the case of some spws with differing VisVector VisTypes
     662             :     }
     663             :     break;
     664           0 :   case Jones::Diagonal:
     665           0 :     switch(vtype) {
     666           0 :     case VisVector::Four:
     667           0 :       return Mueller::Diagonal;
     668             :       break;
     669           0 :     default:
     670           0 :       return Mueller::Diag2;
     671             :       break;
     672             :     }
     673             :     break;
     674           0 :   case Jones::Scalar:
     675           0 :     return Mueller::Scalar;
     676             :     break;
     677           0 :   default:
     678             :     // can't reach here
     679           0 :     throw(AipsError("Cannot figure out Mueller algebra from Jones algebra"));
     680             : 
     681             :   }
     682             :   // Signature consistency (can't reach here)
     683             :   return Mueller::General;
     684             : }
     685             : 
     686             : 
     687             : // print out a Mueller matrix (according to type)
     688           0 : ostream& operator<<(ostream& os, const Mueller& mat) {
     689             : 
     690             :   Complex *mi;
     691           0 :   mi=mat.m_;
     692             : 
     693           0 :   switch (mat.type()) {
     694           0 :   case Mueller::General: {
     695           0 :     cout << "General Mueller: " << endl;;
     696           0 :     for (Int i=0;i<4;++i) {
     697           0 :       cout << " [" << *mi++;
     698           0 :       for (Int j=1;j<4;++j,++mi) cout << ", " << *mi;
     699           0 :       cout << "]";
     700           0 :       if (i<3) cout << endl;
     701             :     }
     702           0 :     break;
     703             :   }
     704           0 :   case Mueller::Diagonal: {
     705           0 :     cout << "Diagonal Mueller: " << endl;
     706           0 :     cout << " [";
     707           0 :     cout << *mi;
     708           0 :     ++mi;
     709           0 :     for (Int i=1;i<4;++i,++mi) cout << ", " << *mi;
     710           0 :     cout << "]";
     711           0 :     break;
     712             :   }
     713           0 :   case Mueller::Diag2: {
     714           0 :     cout << "Diag2 Mueller: " << endl;
     715           0 :     cout << " [" << *mi << ", ";
     716           0 :     ++mi;
     717           0 :     cout << *mi << "]";
     718           0 :     break;
     719             :   }
     720           0 :   case Mueller::AddDiag: {
     721           0 :     cout << "AddDiag Mueller: " << endl;
     722           0 :     cout << " [";
     723           0 :     cout << *mi;
     724           0 :     ++mi;
     725           0 :     for (Int i=1;i<4;++i,++mi) cout << ", " << *mi;
     726           0 :     cout << "]";
     727           0 :     break;
     728             :   }
     729           0 :   case Mueller::AddDiag2: {
     730           0 :     cout << "AddDiag2 Mueller: " << endl;
     731           0 :     cout << " [" << *mi << ", ";
     732           0 :     ++mi;
     733           0 :     cout << *mi << "]";
     734           0 :     break;
     735             :   }
     736           0 :   case Mueller::Scalar: {
     737           0 :     cout << "Scalar Mueller: " << endl;
     738           0 :     cout << " [ " << *mi << " ]";
     739           0 :     break;
     740             :   }
     741             :   }
     742             : 
     743           0 :   return os;
     744             : }
     745             : 
     746             : } //# NAMESPACE CASA - END
     747             : 
     748             : 

Generated by: LCOV version 1.16