coco

--------------------------------------------------------
 coco function      General coordinate convertor.
                  Convert/precess coordinate from/to
                  Equatorial/galactic/Ecliptic.
 Input  : - Matrix of input coordinates.
            First column for Long/RA and second column
            for lat/Dec.
          - Type of input coordinates.
            'j####.#' - equatorial, mean of date (default 'j2000.0'). [Julian].
            'J####.#' - equatorial, true of date (default 'j2000.0'). [Julian].
            'g' - J2000.0 galactic.
            'e' - Ecliptic with J2000.0 equinox.
          - Type of outpt coordinates.
            'j####.#' - equatorial, mean of date (default 'j2000.0'). [Julian].
            'J####.#' - equatorial, true of date (default 'j2000.0'). [Julian].
            'g' - J2000.0 galactic. (default)
            'e' - Ecliptic with J2000.0 equinox.
          - Units for input coordinates.
            'r' - radians. (default)
            'd' - degrees.
            'h' - hours/deg.
          - Units for outpu coordinates.
            'r' - radians. (default)
            'd' - degrees.
            'h' - hours/deg.
 Output : - Matrix of output coordinates.
          - Total rotation matrix.
 Tested : Matlab 5.3
     By : Eran O. Ofek             Febuary 2000 / June 2000 / December 2002
    URL : http://wise-obs.tau.ac.il/~eran/matlab.html 
--------------------------------------------------------

0001 function [OutList,TotRot]=coco(InList,InCooType,OutCooType,InUnits,OutUnits);
0002 %--------------------------------------------------------
0003 % coco function      General coordinate convertor.
0004 %                  Convert/precess coordinate from/to
0005 %                  Equatorial/galactic/Ecliptic.
0006 % Input  : - Matrix of input coordinates.
0007 %            First column for Long/RA and second column
0008 %            for lat/Dec.
0009 %          - Type of input coordinates.
0010 %            'j####.#' - equatorial, mean of date (default 'j2000.0'). [Julian].
0011 %            'J####.#' - equatorial, true of date (default 'j2000.0'). [Julian].
0012 %            'g' - J2000.0 galactic.
0013 %            'e' - Ecliptic with J2000.0 equinox.
0014 %          - Type of outpt coordinates.
0015 %            'j####.#' - equatorial, mean of date (default 'j2000.0'). [Julian].
0016 %            'J####.#' - equatorial, true of date (default 'j2000.0'). [Julian].
0017 %            'g' - J2000.0 galactic. (default)
0018 %            'e' - Ecliptic with J2000.0 equinox.
0019 %          - Units for input coordinates.
0020 %            'r' - radians. (default)
0021 %            'd' - degrees.
0022 %            'h' - hours/deg.
0023 %          - Units for outpu coordinates.
0024 %            'r' - radians. (default)
0025 %            'd' - degrees.
0026 %            'h' - hours/deg.
0027 % Output : - Matrix of output coordinates.
0028 %          - Total rotation matrix.
0029 % Tested : Matlab 5.3
0030 %     By : Eran O. Ofek             Febuary 2000 / June 2000 / December 2002
0031 %    URL : http://wise-obs.tau.ac.il/~eran/matlab.html
0032 %--------------------------------------------------------
0033 RADIAN = 180./pi;
0034 
0035 if (nargin==1),
0036    InCooType  = 'j2000.0';
0037    OutCooType = 'g';
0038    InUnits    = 'r';
0039    OutUnits   = 'r';
0040 elseif (nargin==2),
0041    OutCooType = 'g';
0042    InUnits    = 'r';
0043    OutUnits   = 'r';
0044 elseif (nargin==3),
0045    InUnits    = 'r';
0046    OutUnits   = 'r';
0047 elseif (nargin==4),
0048    OutUnits   = 'r';
0049 elseif (nargin==5),
0050    % do nothing
0051 else
0052    error('Illigal number of input arguments');
0053 end
0054 
0055 LenInType  = length(InCooType);
0056 LenOutType = length(OutCooType);
0057 
0058 if (LenInType>1),
0059    InEquinox  = str2num(InCooType(2:LenInType));
0060    InEquinoxJD = 2451545.5 + 365.25.*(InEquinox - 2000);
0061    InCooType  = InCooType(1);
0062 end
0063 
0064 if (LenOutType>1),
0065    OutEquinox = str2num(OutCooType(2:LenOutType));
0066    OutEquinoxJD = 2451545.5 + 365.25.*(OutEquinox - 2000);
0067    OutCooType  = OutCooType(1);
0068 end
0069 
0070 InCoo   = zeros(size(InList));
0071 OutCoo  = zeros(size(InList));
0072 OutList = zeros(size(InList));
0073 
0074 switch InUnits
0075  case {'r'}
0076     InCoo = InList;
0077  case {'d'}
0078     % convert deg. to radians
0079     InCoo = InList./RADIAN;
0080  case {'h'}
0081     % convert h/d to radians
0082     InCoo(:,1) = InList(:,1).*15./RADIAN;
0083     InCoo(:,2) = InList(:,2)./RADIAN;
0084  otherwise
0085     error('Unknown type of input units');
0086 end
0087 
0088 % convert coordinates to direction cosines
0089 InCosDir = cosined(InCoo);
0090 
0091 
0092 RotM1 = diag([1 1 1]);
0093 
0094 % calculate the first rotation matrix
0095 switch InCooType
0096  case {'j','J'}
0097     if (InEquinox~=2000.0),
0098        % precess coordinates to J2000.0
0099        switch InCooType
0100        case {'j'}
0101            % mean equinox ...
0102            RotM1 = rotm_coo('p',InEquinoxJD);
0103         case {'J'}
0104            % true equinox ...
0105            RotM1 = rotm_coo('pd',InEquinoxJD);
0106         otherwise
0107            error('Illegal InCooType');
0108        end
0109     end
0110  case {'g'}
0111     % convert to Equatorial J2000.0
0112     RotM1 = rotm_coo('G',2451545.5);
0113  case {'e'}
0114     % convert to Equatorial J2000.0
0115     RotM1 = rotm_coo('E',2451545.5);
0116  otherwise
0117     error('Unknown input coordinaytes type');
0118 end
0119 
0120 RotM2 = diag([1 1 1]);
0121 % calculate the second rotation matrix
0122 switch OutCooType
0123  case {'j','J'}
0124     if (OutEquinox~=2000.0),
0125        % precess coordinates from J2000.0
0126        switch OutCooType
0127         case {'j'}   
0128            % mean equinox ...
0129            RotM2 = rotm_coo('P',OutEquinoxJD);
0130         case {'J'}   
0131            % true equinox ...
0132            RotM2 = rotm_coo('Pd',OutEquinoxJD);
0133         otherwise
0134            error('Illegal OutCooType');
0135        end           
0136     end
0137  case {'g'}
0138     % convert to galactic
0139     RotM2 = rotm_coo('g',2451545.5);
0140  case {'e'}
0141     % convert to ecliptic
0142     RotM2 = rotm_coo('e',2451545.5);
0143  otherwise
0144     error('Unknown output coordinaytes type');
0145 end
0146 
0147 % rotate coordinates
0148 TotRot = RotM2*RotM1;
0149 OutCosDir = TotRot*[InCosDir.'];
0150 
0151 
0152 % convert coordinates from direction cosines
0153 OutCoo = cosined([OutCosDir.']);
0154 
0155 
0156 switch OutUnits
0157  case {'r'}
0158     OutList = OutCoo;
0159  case {'d'}
0160     % convert radians to deg.
0161     OutList = OutCoo.*RADIAN;
0162  case {'h'}
0163     % convert radians to h/d
0164     OutList(:,1) = OutCoo(:,1).*RADIAN./15;
0165     OutList(:,2) = OutCoo(:,2).*RADIAN;
0166  otherwise
0167     error('Unknown type of output units');
0168 end
0169 
0170 
0171 
0172 
0173 
0174