parse_radius.cpp

Go to the documentation of this file.

/* This file is part of Cloudy and is copyright (C)1978-2013 by Gary J. Ferland and
 * others.  For conditions of distribution and use see copyright notice in license.txt */
/*ParseRadius parse the radius command */
#include "cddefines.h"
/*#define       PARSCL  18.489396*/
#include "physconst.h"
#include "optimize.h"
#include "radius.h"
#include "iterations.h"
#include "input.h"
#include "parser.h"
 
void ParseRadius(Parser &p)
{
        bool 
          lgR2set, 
          lgRLog;
        double a, 
          convl,
          r;
 
        DEBUG_ENTRY( "ParseRadius()" );
 
        /* log of inner and outer radii, default second=infinity,
         * if R2<R1 then R2=R1+R2
         * there is an optional keyword, "PARSEC" on the line, to use PC as units */
        convl = p.nMatch("PARS") ? log10( PARSEC ) : 0.;
 
        /* if linear appears on line, then radius is linear, otherwise, log */
        lgRLog = !p.nMatch("LINE");
 
        r = p.FFmtRead();
        if( p.lgEOL() )
                p.NoNumb("radius");
 
        /* option for linear or log radius */
        if( lgRLog )
        {
                r += convl;
        }
        else
        {
                if( r > 0. )
                {
                        r = log10(r) + convl;
                }
                else
                {
                        fprintf(ioQQQ,"The first radius is negative and linear is set - this is impossible.\n");
                        cdEXIT(EXIT_FAILURE);
                }
        }
 
        if( r > 37. || r < -37. )
        {
                fprintf(ioQQQ,"WARNING - the log of the radius is %e - this is too big.\n", r );
                fprintf(ioQQQ," Sorry.\n" );
                cdEXIT(EXIT_FAILURE);
        }
 
        radius.Radius = pow(10.,r);
        radius.lgRadiusKnown = true;
 
        /* check for second number, which indicates thickness or outer radius of model */
        a = p.FFmtRead();
        if( p.lgEOL() )
        {
                /* not set */
                lgR2set = false;
        }
        else
        {
                /* outer radius is set, */
                lgR2set = true;
 
                /* log or linear option is still in place */
                if( lgRLog )
                {
                        a += convl;
                }
                else
                {
                        /* linear radius - convert to log but first make sure that a is > 0 */
                        if( a > 0. )
                        {
                                a = log10(a) + convl;
                        }
                        else
                        {
                                fprintf(ioQQQ,"The second radius is negative and linear is set - this is impossible.\n");
                                cdEXIT(EXIT_FAILURE);
                        }
                }
 
                if( a > 37. || a < -37. )
                {
                        fprintf(ioQQQ,"WARNING - the log of the second radius is %e - this is too big.\n", a );
                        /* flush buffers since we shall soon throw an fpe */
                        fflush( ioQQQ );
                }
                a = pow(10.,a);
                /* check whether it was thickness or outer radius,
                 * we want thickness to be total thickness of modeled region,
                 * NOT outer radius */
                if( a > radius.Radius )
                        radius.StopThickness[0] = a - radius.Radius;
                else
                        radius.StopThickness[0] = a;
 
                for( long int i=1; i < iterations.iter_malloc; i++ )
                {
                        radius.StopThickness[i] = radius.StopThickness[0];
                }
        }
 
        /* vary option */
        if( optimize.lgVarOn )
        {
                /* pointer to where to write */
                optimize.nvfpnt[optimize.nparm] = input.nRead;
 
                /* flag saying second outer radius or thickness was set */
                if( lgR2set )
                {
                        strcpy( optimize.chVarFmt[optimize.nparm], "RADIUS %f depth or outer R %f LOG" );
                        optimize.nvarxt[optimize.nparm] = 2;
                        /* second number is thickness or outer radius */
                        optimize.vparm[1][optimize.nparm] = (realnum)log10(a);
                        fprintf(ioQQQ,
                                " WARNING - outer radius or thickness was set with a variable radius.\n");
                        fprintf(ioQQQ,
                                " The interpretation of the second number can change from radius to depth as radius changes.\n");
                        fprintf(ioQQQ,
                                " Do not use the second parameter unless you are certain that you know what you are doing.\n");
                        fprintf(ioQQQ,
                                " Consider using the STOP THICKNESS or STOP RADIUS command instead.\n");
                }
                else
                {
                        strcpy( optimize.chVarFmt[optimize.nparm], "RADIUS= %f LOG" );
                        optimize.nvarxt[optimize.nparm] = 1;
                }
 
                /* log of radius is first number */
                optimize.vparm[0][optimize.nparm] = (realnum)log10(radius.Radius);
                optimize.vincr[optimize.nparm] = 0.5;
                ++optimize.nparm;
        }
        return;
}