Changeset 25143


Ignore:
Timestamp:
06/24/20 19:13:14 (5 years ago)
Author:
Eric.Larour
Message:

CHG: changing name.

File:
1 copied

Legend:

Unmodified
Added
Removed
  • TabularUnified issm/trunk-jpl/src/m/boundaryconditions/getlovenumbers.m

    r21307 r25143  
    1 function series=love_numbers(type,reference_frame)
     1function series=getlovenumbers(varargin)
    22%LOVE_NUMBERS: provide love numbers
    33%       retrieved from: http://www.srosat.com/iag-jsg/loveNb.php
     
    55
    66%       Usage:
    7 %       series = love_numbers(type)
    8 %       series = love_numbers(type,reference_frame)
     7%       series = love_numbers('type','loadingverticaldisplacement','referenceframe','CM','maxdeg',1001)
    98%
    10 %       - type = one of 'h','k','l','gamma' and 'lambda'.
     9%       - type = one of 'loadingverticaldisplacement', 'loadinggravitationalpotential', 'loadinghorizontaldisplacement',  ...
     10%                    'tidalverticaldisplacement', 'tidalgravitationalpotential', 'tidalhorizontaldisplacement'
    1111%       - reference_frame = one of 'CM' (default) and 'CF'.
     12%   - maxdeg = default 1001
    1213%
    1314%       Example: 
    14 %       love_k = love_numbers('k');
    15 %       love_k = love_numbers('k','CF');
    16 %
     15%       h=love_number('type','loadingverticaldisplacement','referenceframe','CM','maxdeg',maxdeg);
     16%       k=love_number('type','loadinggravitationalpotential','referenceframe','CM','maxdeg',maxdeg);
     17%       l=love_number('type','loadinghorizontaldisplacement','referenceframe','CM','maxdeg',maxdeg);
     18%       th=love_number('type','tidalverticaldisplacement','referenceframe','CM','maxdeg',maxdeg);
     19%       tk=love_number('type','tidalgravitationalpotential','referenceframe','CM','maxdeg',maxdeg);
     20%       tl=love_number('type','tidalhorizontaldisplacement','referenceframe','CM','maxdeg',maxdeg);
    1721
    18 % some checks:
    19 if nargin==1,
    20         frame='CM';
    21         disp('Info: computation is done in Center of Mass (CM) reference frame');
    22 elseif nargin==2,
    23         frame=reference_frame;
    24         if ~( strcmpi(reference_frame,'CM') | strcmpi(reference_frame,'CF')),
    25                 error('reference_frame should be one of ''CM'' or ''CF''');
    26         end
    27 else
    28         help love_numbers
    29         error('bad usage');
     22%recover options:
     23options=pairoptions(varargin{:});
     24type=getfieldvalue(options,'type');
     25frame=getfieldvalue(options,'referenceframe','CM');
     26maxdeg=getfieldvalue(options,'maxdeg',1000);
     27
     28if(maxdeg>10000),
     29        error('PREM love numbers computed only for deg < 10,000. Request lower maxdeg' );
    3030end
    3131
    32 if ~( strcmpi(type,'h') | strcmpi(type,'k') | strcmpi(type,'l') | strcmpi(type,'gamma') | strcmpi(type,'lambda') ),
    33         error('type should be one of ''h'',''k'',''l'',''gamma'' and ''lambda''');
     32if ~( strcmpi(type,'loadingverticaldisplacement') | strcmpi(type,'loadinggravitationalpotential') | strcmpi(type,'loadinghorizontaldisplacement') ...
     33                | strcmpi(type,'tidalverticaldisplacement') | strcmpi(type,'tidalgravitationalpotential') | strcmpi(type,'tidalhorizontaldisplacement')),
     34        error('type should be one of ''loadingverticaldisplacement'', ''loadinggravitationalpotential'', ''loadinghorizontaldisplacement'',''tidalverticaldisplacement'', ''tidalgravitationalpotential'', ''tidalhorizontaldisplacement''');
    3435end
    3536
    36                 love_numbers=[...
     37        love_numbers=[...
    3738     0           0          0          0          0          0          0
    3839        -1.28740059     -1.00000000     -0.89858519 1.28740059 0.42519882 0.89858519 0.00000000
     
    1003710038        -6.27342778 -0.00030945 0.00018956 7.27311833 0.99905480 0.99950099 0.49327194];
    1003810039
    10039         if type=='h',
     10040        %cut off series at maxdeg:
     10041        love_numbers(maxdeg+2:end,:)=[];
     10042
     10043        %retrive right type:
     10044        if strcmpi(type,'loadingverticaldisplacement'),
    1004010045                series=love_numbers(:,1);
    10041         elseif type=='k',
     10046        elseif strcmpi(type,'loadinggravitationalpotential'),
    1004210047                series=love_numbers(:,2);
    10043         elseif type=='l',
     10048        elseif strcmpi(type,'loadinghorizontaldisplacement'),
    1004410049                series=love_numbers(:,3);
    10045         elseif type=='gamma',
     10050        elseif strcmpi(type,'tidalverticaldisplacement'),
    1004610051                series=love_numbers(:,4);
    10047         elseif type=='lambda',
     10052        elseif strcmpi(type,'tidalgravitationalpotential'),
    1004810053                series=love_numbers(:,5);
     10054        elseif strcmpi(type,'tidalhorizontaldisplacement'),
     10055                series=love_numbers(:,6);
    1004910056        else
    1005010057                error(['love_numbers error message: unknow type:' type]);
     
    1005210059
    1005310060        % choose degree 1 term for CF reference system
    10054         if frame=='CM',
     10061        if strcmpi(frame,'CM'),
    1005510062                return;
    10056         elseif frame=='CF', % from Blewitt, 2003, JGR
    10057                 if type=='h',
     10063        elseif strcmpi(frame,'CF'), % from Blewitt, 2003, JGR
     10064                if strcmpi(type,'loadingverticaldisplacement'),
    1005810065                        series(2,1) = -0.269;
    10059                 elseif type=='k',
     10066                elseif strcmpi(type,'loadinggravitationalpotential'),
    1006010067                        series(2,1) = 0.021; 
    10061                 elseif type=='l',
     10068                elseif strcmpi(type,'loadinghorizontaldisplacement'),
    1006210069                        series(2,1) = 0.134;
    1006310070                end
     
    1006610073        end
    1006710074
     10075
     10076
Note: See TracChangeset for help on using the changeset viewer.