/* NIGHTFALL Light Curve Synthesis Program */ /* Copyright (C) 1998 Rainer Wichmann */ /* */ /* This program is free software; you can redistribute it */ /* and/or modify */ /* it under the terms of the GNU General Public License as */ /* published by */ /* the Free Software Foundation; either version 2 of the License, or */ /* (at your option) any later version. */ /* */ /* This program is distributed in the hope that it will be useful, */ /* but WITHOUT ANY WARRANTY; without even the implied warranty of */ /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */ /* GNU General Public License for more details. */ /* */ /* You should have received a copy of the GNU General Public License */ /* along with this program; if not, write to the Free Software */ /* Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */ #include #include #include #include #include #include #include "Light.h" /****************************************************************** @package nightfall @author Rainer Wichmann (rwichman@lsw.uni-heidelberg.de) @version 1.0 @short Parse the command line @tip This is the most awkward function in the whole program ... @param (int) argc The number of arguments @param (char *[]) argv The arguments @return (void) @heading Main *******************************************************************/ /* Get The Command Line Arguments */ void GetArguments(int *argc, char *argv[]) { int i; /* loop variable */ int Ichar; int CheckIn, mandatory = 0; char infile[256]; char ErrMsg[256]; char *debug; double G_SI = 6.6726e-11; /* gravitational constant */ int FitDim = 0; char Fitopt[] = "0123456789ABCDEFGHIJKLMNOPabcdefghijkl"; /* these are globals, we set the default values */ if (PHASESTEPS > 80) PhaseSteps = 80; else PhaseSteps = PHASESTEPS; StepsPerPi = STEPS_PER_PI; MaximumElements = MAXELEMENTS; /* optional arguments */ while ((*argc) > 1) { if ((argv[1][0] == '-') && (argv[1][1] != '-')) { switch (argv[1][1]) { case 'h': UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit (0); break; case 'D': debug = &argv[1][2]; while (debug != NULL && isalpha (*debug)) { if ( ! isalpha (*debug) ) { WARNERR(_("Invalid character in debug string encountered, skipped"), argv[1]); } else { if ( (*debug) == 'w') Flags.debug[WARN] = 1; else if ( (*debug) == 'v') Flags.debug[VERBOSE] = 1; else if ( (*debug) == 'b') Flags.debug[BUSY] = 1; else if ( (*debug) == 'g') Flags.debug[GEO] = 1; else if ( (*debug) == 'i') Flags.debug[INFILE] = 1; else if ( (*debug) == 'm') Flags.debug[MINFIND] = 1; else if ( (*debug) == 'q') Flags.debug[CHISQR] = 1; else if ( (*debug) == 'o') Flags.debug[OPTIMIZE] = 1; else if ( (*debug) == 's') Flags.debug[SURFACE] = 1; } debug++; } break; case 'U': Flags.interactive = ON; break; case 'A': Flags.animate = ON; break; #ifdef _WITH_OPENGL case 'g': Flags.use_opengl = OFF; break; #endif case 'V': if (argv[1][2] == 'i') Flags.visualize = 2; else if (argv[1][2] == 'c') Flags.visualize = 3; else if (argv[1][2] == 'a') Flags.visualize = 4; else Flags.visualize = ON; break; case 'G': if (argv[1][2] == 'P') Flags.plot = 3; else if (argv[1][2] == '2') Flags.plot = 2; else if (argv[1][2] == 'S') Flags.plot = 4; else Flags.plot = ON; break; case 'H': Flags.eps = ON; break; case 'N': CheckIn = sscanf(&argv[1][2], "%d", &PhaseSteps); if (CheckIn != 1) { WARNERR(_("Steps for Lightcurve (maybe a 'blank' after 'N'), Default used"), argv[1] ); PhaseSteps = IMIN(PHASESTEPS, 80); } else if (PhaseSteps > PHASESTEPS) { sprintf(ErrMsg, _("Lightcurve Steps > Maximum (%d), Default used"), PHASESTEPS); WARNERR(ErrMsg, argv[1] ); PhaseSteps = IMIN(PHASESTEPS, 80); } else if (PhaseSteps < 0) { WARNERR(_("Lightcurve Steps < Zero, Default used"), argv[1] ); PhaseSteps = IMIN(PHASESTEPS, 80); } break; case 'P': CheckIn = sscanf(&argv[1][2], "%lf", &Orbit.LambdaZero); if (CheckIn != 1) { WARNERR(_("Line Wavelength missing (maybe a 'blank' after 'P'), default used"), argv[1] ); Orbit.LambdaZero = LAMDAZERO_DEFAULT; } else if ( (Orbit.LambdaZero <= (LIM_PRF_L-FLT_EPSILON)) || (Orbit.LambdaZero >= (LIM_PRF_H+FLT_EPSILON*LIM_PRF_H))){ sprintf(ErrMsg, _("Line Profile Rest Wavelength must be in [%f, %f]"), LIM_PRF_L, LIM_PRF_H ); WARNERR(ErrMsg, argv[1]); Orbit.LambdaZero = LAMDAZERO_DEFAULT; } Flags.lineprofile = ON; break; case 'X': debug = &argv[1][2]; if ( ! isalnum (*debug) ) { sprintf(ErrMsg, _("Bad Option for Simplex Fit %c, skipped"), *debug); WARNERR(ErrMsg, argv[1]); } while ( *debug != '\0' ) { for (Ichar = 0; Ichar < (int) strlen(Fitopt); ++Ichar) { if (*debug == Fitopt[Ichar]) Flags.simplex[Ichar] = ON; } debug++; } if (argv[2] != NULL) CheckIn = sscanf(argv[2],"%f", &Flags.SimplexTol); else CheckIn = 0; if (CheckIn != 1) { WARNERR(_("Simplex Fit: No Tolerance specified, default used"), argv[1]); Flags.SimplexTol = SIMPLEX_FIT_TOLERANCE; } else { ++argv; --(*argc); } Flags.WantFit = ON; if (fabs (Flags.SimplexTol - 0.001) <= FLT_EPSILON) Flags.anneal = ON; break; case 'Y': debug = &argv[1][2]; if ( ! isalnum (*debug) ) { sprintf(ErrMsg,_("Bad Option for ChiSquare Map %c, skipped"), *debug); WARNERR(ErrMsg, argv[1]); } while ( *debug != '\0' ) { for (Ichar = 0; Ichar < (int) strlen(Fitopt); ++Ichar) { if (*debug == Fitopt[Ichar]) Flags.simplex[Ichar] = ON; } debug++; } ++argv; --(*argc); if (argv[1] != NULL) CheckIn = sscanf(argv[1],"%f", &Flags.Step[0]); else CheckIn = 0; if (CheckIn != 1) INERR(_("Chi Square Map: No Step1 specified"), argv[1]) ++argv; --(*argc); if (argv[1] != NULL) CheckIn = sscanf(argv[1],"%f", &Flags.Step[1]); else CheckIn = 0; if (CheckIn != 1) INERR(_("Chi Square Map: No Step2 specified"), argv[1]) Flags.WantMap = ON; break; case 'M': Flags.blackbody = OFF; break; case 'F': Flags.fractional = ON; break; case 'L': (void) sscanf(&argv[1][2], "%d", &Flags.limb); if (Flags.limb > LD_LAW_MAX || Flags.limb < 0) { WARNERR(_("Invalid limb darkening option (maybe 'blank' after 'L'), linear used"), argv[1]); Flags.limb = 0; } break; case 'R': (void) sscanf(&argv[1][2], "%d", &Flags.reflect); if (Flags.reflect > 9 || Flags.reflect < 1) { WARNERR(_("Invalid Reflection Iteration Count (maybe 'blank' after 'R'), 3 used"), argv[1]); Flags.reflect = 3; } break; case 'B': switch (argv[1][2]) { case 'U': Flags.PlotBand = Umag; break; case 'B': Flags.PlotBand = Bmag; break; case 'V': Flags.PlotBand = Vmag; break; case 'R': Flags.PlotBand = Rmag; break; case 'I': Flags.PlotBand = Imag; break; case 'u': Flags.PlotBand = umag; break; case 'v': Flags.PlotBand = vmag; break; case 'b': Flags.PlotBand = bmag; break; case 'y': Flags.PlotBand = ymag; break; case 'J': Flags.PlotBand = Jmag; break; case 'H': Flags.PlotBand = Hmag; break; case 'K': Flags.PlotBand = Kmag; break; case '1': Flags.PlotBand = mag1; break; case '2': Flags.PlotBand = mag2; break; default: WARNERR(_("Bad Option For Bandpass to Graph, V Band used"), argv[1]) break; } break; case '3': switch (argv[1][2]) { case 'U': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Umag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'B': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Bmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'V': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Vmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'R': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Rmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'I': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Imag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'u': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[umag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'v': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[vmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'b': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[bmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'y': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[ymag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'J': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Jmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'H': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Hmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; case 'K': CheckIn = sscanf(&argv[1][3], "%lf", &Orbit.Third[Kmag]); if (CheckIn != 1) WARNERR(_("No Third Light given (misplaced 'blank' ?)"), argv[1]) break; default: WARNERR(_("Bad Option For Third Light"), argv[1]) break; } break; case 'C': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } else if (strlen(argv[2]) < (sizeof(infile)-1)) { (void) sscanf(argv[2], "%s", infile); ParseConfig(infile, &mandatory); ++argv; --(*argc); } else { WARNERR(_("Path to Config File exceeds Maximum Length (255), skipped"), argv[2]); ++argv; --(*argc); } break; case 'I': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } else if (strlen(argv[2]) < (sizeof(infile)-1)) { (void) sscanf(argv[2], "%s", infile); Read(infile); ++argv; --(*argc); } else { WARNING(_("Path to Input File exceeds Maximum Length (255), skipped")); ++argv; --(*argc); } break; case 's': if ((*argc) < 5) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } switch (argv[1][2]) { case 'P': ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Primary][Flags.Spots1].longitude); if (CheckIn != 1) INERR(_("No Longitude for Spot"), argv[1] ); if( (Spot[Primary][Flags.Spots1].longitude-LIM_SLO_L) <=(-FLT_EPSILON) || (Spot[Primary][Flags.Spots1].longitude-LIM_SLO_H)>=FLT_EPSILON) INERR(_("Spot longitude must be in [-360,360]"), argv[1]); ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Primary][Flags.Spots1].latitude); if (CheckIn != 1) INERR(_("No Latitude for Spot"), argv[1]); if((Spot[Primary][Flags.Spots1].latitude-LIM_SLA_L) <=(-FLT_EPSILON) || (Spot[Primary][Flags.Spots1].latitude-LIM_SLA_H)>=FLT_EPSILON) INERR(_("Spot latitude must be in [-90,90]"), argv[1]); ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf",&Spot[Primary][Flags.Spots1].radius); if (CheckIn != 1) INERR(_("No Radius for Spot"), argv[1]); if((Spot[Primary][Flags.Spots1].radius-LIM_SRA_L)<=(-FLT_EPSILON)|| (Spot[Primary][Flags.Spots1].radius-LIM_SRA_H)>=FLT_EPSILON) INERR(_("Spot radius must be in [1,180]"), argv[1] ); ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Primary][Flags.Spots1].dimfactor); if (CheckIn != 1) INERR(_("No Dim Factor for Spot"), argv[1]); if((Spot[Primary][Flags.Spots1].dimfactor-LIM_SDF_L) <=(-FLT_EPSILON)) INERR(_("Spot Dim Factor must be > 0.5"), argv[1]); if((Spot[Primary][Flags.Spots1].dimfactor-LIM_SDF_H)>=FLT_EPSILON) INERR(_("Spot Dim Factor must be < 2. "), argv[1]); ++Flags.Spots1; break; case 'S': ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Secondary][Flags.Spots2].longitude); if (CheckIn != 1) INERR(_("No Longitude for Spot"), argv[1] ); if((Spot[Secondary][Flags.Spots2].longitude-LIM_SLO_L) <=(-FLT_EPSILON) || (Spot[Secondary][Flags.Spots2].longitude-LIM_SLO_H)>=FLT_EPSILON) INERR(_("Spot longitude must be in [-360,360]"), argv[1]); ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Secondary][Flags.Spots2].latitude); if (CheckIn != 1) INERR(_("No Latitude for Spot"), argv[1]); if((Spot[Secondary][Flags.Spots2].latitude-LIM_SLA_L) <=(-FLT_EPSILON) || (Spot[Secondary][Flags.Spots2].latitude-LIM_SLA_H)>=FLT_EPSILON) INERR(_("Spot latitude must be in [-90,90]"), argv[1]); ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Secondary][Flags.Spots2].radius); if (CheckIn != 1) INERR(_("No Radius for Spot"), argv[1]); if((Spot[Secondary][Flags.Spots2].radius-LIM_SRA_L) <=(-FLT_EPSILON) || (Spot[Secondary][Flags.Spots2].radius-LIM_SRA_H)>=FLT_EPSILON) INERR(_("Spot radius must be in [1,180]"), argv[1] ); ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Spot[Secondary][Flags.Spots2].dimfactor); if (CheckIn != 1) INERR(_("No Dim Factor for Spot"), argv[1]); if((Spot[Secondary][Flags.Spots2].dimfactor-LIM_SDF_L) <=(-FLT_EPSILON)) INERR(_("Spot Dim Factor must be > 0.5"), argv[1]); if((Spot[Secondary][Flags.Spots2].dimfactor-LIM_SDF_H) >=FLT_EPSILON) INERR(_("Spot Dim Factor must be < 2."), argv[1]); ++Flags.Spots2; break; default: WARNERR(_("Bad option for spot, skipped"), argv[1]); break; } break; case 't': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } switch (argv[1][2]) { case 'P': ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Orbit.TruePeriod); if (CheckIn != 1) WARNERR(_("No Period given (maybe a 'blank' after 'P')"), argv[1]); if (Orbit.TruePeriod <= LIM_PERI_L) INERR(_("Absolute Period"), argv[1]); /* convert days to seconds */ Orbit.TruePeriod = Orbit.TruePeriod * 86400.; break; case 'M': ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Orbit.TrueMass); if (CheckIn != 1) WARNERR(_("No Mass given (maybe a 'blank' after 'M')"), argv[1]); if (Orbit.TrueMass <= LIM_MASS_L) INERR(_("Absolute Mass"), argv[1]); /* convert solar masses to KG */ Orbit.TrueMass = Orbit.TrueMass * 1.989E30; break; case 'D': ++argv; --(*argc); CheckIn = sscanf(argv[1], "%lf", &Orbit.TrueDistance); if (CheckIn != 1) WARNERR(_("No Separation given (maybe a 'blank' after 'D')"), argv[1]); if (Orbit.TrueDistance <= LIM_DIST_L) INERR(_("Absolute Separation"), argv[1]); /* convert solar radius to meter */ Orbit.TrueDistance = Orbit.TrueDistance * 6.960E8; break; default: WARNERR(_("Bad option for absolute values, default used"), argv[1]); break; } break; case 'e': if ((*argc) < 4) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Orbit.Excentricity); if (CheckIn != 1) INERR(_("No Eccentricity found"), argv[1]); if ( (Orbit.Excentricity - LIM_EX_L) <= (-FLT_EPSILON) || (Orbit.Excentricity - LIM_EX_H) >= FLT_EPSILON ) INERR(_("Eccentricity must be in (0., 0.95)"), argv[1]); ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Orbit.Omega); if (CheckIn != 1) INERR(_("No Periastron Length found"), argv[1]); if ((Orbit.Omega - LIM_PA_L) <= (-FLT_EPSILON) || (Orbit.Omega - LIM_PA_H) >= FLT_EPSILON ) INERR(_("Periastron Length must be in [-360, 360]"), argv[1]); Flags.elliptic = ON; break; case 'f': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } switch (argv[1][2]) { case 'P': ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Binary[Primary].Fratio); if (CheckIn != 1) INERR(_("No Rotation Rate found"), argv[1]); if( (Binary[Primary].Fratio - LIM_FR_L) <= (-FLT_EPSILON) || (Binary[Primary].Fratio - LIM_FR_H) >= FLT_EPSILON ) INERR(_("Rotation Rate must be in (0., 10.]"), argv[1]); Flags.asynchron1 = ON; break; case 'S': ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Binary[Secondary].Fratio); if (CheckIn != 1) INERR(_("No Rotation Rate found"), argv[1]); if( (Binary[Secondary].Fratio - LIM_FR_L) <= (-FLT_EPSILON) || (Binary[Primary].Fratio - LIM_FR_H) >= FLT_EPSILON ) INERR(_("Rotation Rate must be in (0., 10.]"), argv[1]); Flags.asynchron2 = ON; break; default: WARNERR(_("Bad Option (use 'nightfall to get a list of options): "), argv[1]); break; } break; case 'O': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } switch (argv[1][2]) { case 'P': ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Binary[Primary].log_g); if (CheckIn != 1) INERR(_("No log g found"), argv[1]); if( (Binary[Primary].log_g - LIM_LOGG_L) <= (-FLT_EPSILON) || (Binary[Primary].log_g - LIM_LOGG_H) >= FLT_EPSILON ) INERR(_("Log g must be in [3.5, 5.0]"), argv[1]); Binary[Primary].log_g = 0.5 * ROUND(2.0 * Binary[Primary].log_g); break; case 'S': ++argv; --(*argc); CheckIn = sscanf(argv[1],"%lf", &Binary[Secondary].log_g); if (CheckIn != 1) INERR(_("No log_g found"), argv[1]); if( (Binary[Secondary].log_g - LIM_LOGG_L) <= (-FLT_EPSILON) || (Binary[Primary].log_g - LIM_LOGG_H) >= FLT_EPSILON ) INERR(_("Log g must be in [3.5, 5.0]"), argv[1]); Binary[Secondary].log_g = 0.5 * ROUND(2.0 * Binary[Secondary].log_g); break; default: WARNERR(_("Bad Option (use 'nightfall to get a list of options): "), argv[1]); break; } break; case 'a': if ((*argc) < 3) { UsageInfo(); MPI_Abort(MPI_COMM_WORLD, EXIT_FAILURE); exit(EXIT_FAILURE); } switch (argv[1][2]) { case 'P': ++argv; --(*argc); SetAlbedo(Primary, argv[1]); break; case 'S': ++argv; --(*argc); SetAlbedo(Secondary, argv[1]); break; default: WARNERR(_("Bad Option (use 'nightfall to get a list of options): "), argv[1]); break; } break; default: WARNERR(_("Bad Option (use 'nightfall to get a list of options): "), argv[1]); } ++argv; --(*argc); } else if ((argv[1][0] == '-') && (argv[1][1] == '-')) { /* --------------- ignore ----------------------------------------- */ ++argv; --(*argc); if ( (*argc) > 1 ) { if (argv[1][0] != '-') { ++argv; --(*argc); } } } else { /* -------------- Argument did not start with a '-' or '--' ----- */ break; } } /* >>>>>>>>>>>>>>>>>>>>>>>> mandatory arguments <<<<<<<<<<<<<<<<<<< */ if ( ((*argc) < 7) && (mandatory < 6) && (Flags.interactive == OFF)) nf_error(_("One of the mandatory Arguments is missing. Type 'nightfall' to get more info. (One possible reason might be a missing argument for a preceding option)")); /* + now we scan the rest of the CL. if there are gtk arguments at end, + they will be passed. if they are in front of the mandatory arguments, + a read error will occur (hopefully) */ if ((*argc) >= 7 ) { (void) sscanf(argv[1], "%lf", &Binary[Primary].Mq); if((Binary[Primary].Mq - LIM_MQ_L) <= (-FLT_EPSILON) || (Binary[Primary].Mq - LIM_MQ_H) >= FLT_EPSILON ) INERR(_("mass ratio must be in [0.001,50.]"), argv[1] ); Binary[Secondary].Mq = 1.0/Binary[Primary].Mq; ++argv; --(*argc); (void) sscanf(argv[1], "%lf", &Orbit.Inclination); if((Orbit.Inclination-LIM_IN_L) <= (-FLT_EPSILON) || (Orbit.Inclination - RTOD*LIM_IN_H) >= FLT_EPSILON) INERR(_("Inclination must be in [0,90]"), argv[1]); Orbit.Inclination = DTOR * Orbit.Inclination; ++argv; --(*argc); (void) sscanf(argv[1], "%lf", &Binary[Primary].RocheFill); if((Binary[Primary].RocheFill - LIM_RF_L) <= (-FLT_EPSILON) || (Binary[Primary].RocheFill - LIM_RO_H) >= FLT_EPSILON) INERR(_("Roche lobe filling must be in [0.001,1.3]"),argv[1] ); ++argv; --(*argc); (void) sscanf(argv[1], "%lf", &Binary[Secondary].RocheFill); if((Binary[Secondary].RocheFill - LIM_RF_L) <= (-FLT_EPSILON) || (Binary[Secondary].RocheFill - LIM_RO_H) >= FLT_EPSILON) INERR(_("Roche lobe filling must be in [0.001,1.3]"),argv[1] ); ++argv; --(*argc); (void) sscanf(argv[1], "%lf", &Binary[Primary].Temperature); if(Flags.blackbody == OFF) { if((Binary[Primary].Temperature - LIM_TM_L ) <= (-FLT_EPSILON) || (Binary[Primary].Temperature - LIM_TM_H(Binary[Primary].log_g)) >= FLT_EPSILON) INERR(_("Temperature must be in [3000,35000]"), argv[1]); } else { if((Binary[Primary].Temperature - LIM_TB_L) <= (-FLT_EPSILON) || (Binary[Primary].Temperature - LIM_TB_H )>= FLT_EPSILON) INERR(_("Temperature must be in [350,350000]"), argv[1]); } ++argv; --(*argc); (void) sscanf(argv[1], "%lf", &Binary[Secondary].Temperature); if(Flags.blackbody == OFF) { if((Binary[Secondary].Temperature - LIM_TM_L) <= (-FLT_EPSILON) || (Binary[Secondary].Temperature - LIM_TM_H(Binary[Secondary].log_g)) >= FLT_EPSILON) INERR(_("Temperature must be in [3000,35000]"), argv[4]); } else { if((Binary[Secondary].Temperature - LIM_TB_L) <= (-FLT_EPSILON) || (Binary[Secondary].Temperature - LIM_TB_H) >= FLT_EPSILON) INERR(_("Temperature must be in [350,350000]"), argv[1]); } ++argv; --(*argc); } /* >>>>>>>>>>>>>>>>>>> ADJUST INPUT <<<<<<<<<<<<<<<<<<<<<<<< */ for (i = 0; i < 32; ++i) { if (Flags.simplex[i] == ON) ++FitDim; } if (FitDim < 2 && Flags.WantFit == ON) { Flags.WantFit = OFF; WARNING(_("Less than two Parameters to fit, discard Fit Option")); } /* no overcontact if elliptic */ if ( (Flags.elliptic == ON) && (Binary[Primary].RocheFill >= (1+FLT_EPSILON) ) ) { Binary[Primary].RocheFill = 1.0; WARNING(_("Eccentric Orbit --> Decrease Roche Fill Factor to One")); } if ( (Flags.elliptic == ON) && (Binary[Secondary].RocheFill >= (1+FLT_EPSILON) ) ) { Binary[Secondary].RocheFill = 1.0; WARNING(_("Eccentric Orbit --> Decrease Roche Fill Factor to One")); } /* no overcontact if disk */ /* no disk if elliptic */ #ifdef HAVE_DISK if ( (Flags.elliptic == ON) && (Flags.disk == ON) ) { Flags.disk = OFF; WARNING(_("Eccentric Orbit --> No accretion disk possible")); } if ( (Flags.disk == ON) && (Binary[Primary].RocheFill >= (1+FLT_EPSILON) ) ) { Binary[Primary].RocheFill = 1.0; WARNING(_("Eccentric Orbit --> Decrease Roche Fill Factor to One")); } if ( (Flags.disk == ON) && (Binary[Secondary].RocheFill >= (1+FLT_EPSILON) ) ) { Binary[Secondary].RocheFill = 1.0; WARNING(_("Eccentric Orbit --> Decrease Roche Fill Factor to One")); } #endif /* common fillfactor if overcontact */ if ( Binary[Primary].RocheFill >= (1+FLT_EPSILON) || Binary[Secondary].RocheFill >= (1+FLT_EPSILON) ){ if ( fabs(Binary[Primary].RocheFill-Binary[Secondary].RocheFill) >= FLT_EPSILON) WARNING(_("Overcontact --> Set Common Roche Fill Factor to Maximum of Both")); Binary[Primary].RocheFill = MAX(Binary[Primary].RocheFill,Binary[Secondary].RocheFill); Binary[Secondary].RocheFill = Binary[Primary].RocheFill; Flags.fill = ON; } else { Flags.fill = OFF; } /* Limits for Mass Ratio */ if (Flags.fill == ON) { if (Binary[Primary].Mq >= LIM_MQO_H) { Binary[Primary].Mq = LIM_MQO_H; Binary[Secondary].Mq = 1.0/Binary[Primary].Mq; WARNING(_("Overcontact --> Reduce Mass Ratio to limit")); } } else { if (Binary[Primary].Mq <= LIM_MQ_L ) { Binary[Primary].Mq = LIM_MQ_L; Binary[Secondary].Mq = 1.0/Binary[Primary].Mq; WARNING(_("Avoid Numerical Problems --> Increase Mass Ratio to limit")); } } /* no nonsynchroneous rotation if overcontact */ if (Flags.fill == ON) { if ( (fabs(Binary[Primary].Fratio - 1.0) >= FLT_EPSILON) || (fabs(Binary[Secondary].Fratio - 1.0) >= FLT_EPSILON) ) { WARNING(_("Overcontact --> Reset NonSync Rotation to Unity")); Binary[Primary].Fratio = 1.0; Binary[Secondary].Fratio = 1.0; Flags.asynchron1 = OFF; Flags.asynchron2 = OFF; } } /* ----------- set default band to plot -------------------- */ if (Flags.PlotBand == -1) Flags.PlotBand = Vmag; ComputeGravDark(); /* >>>>>>>>>>> consistency of mass,period, and distance <<<<<<<<<<<<<< */ if (fabs(Orbit.TruePeriod) >= FLT_EPSILON && fabs(Orbit.TrueDistance ) >= FLT_EPSILON ) { Orbit.TrueMass = (Orbit.TrueDistance * Orbit.TrueDistance) * Orbit.TrueDistance * 4.0 * PI * PI / (G_SI * (Orbit.TruePeriod * Orbit.TruePeriod)); } else if (fabs(Orbit.TruePeriod ) >= FLT_EPSILON && fabs(Orbit.TrueMass ) >= FLT_EPSILON ) { Orbit.TrueDistance = pow((Orbit.TruePeriod*Orbit.TruePeriod) * G_SI * Orbit.TrueMass / (4.0 * (PI*PI)), (1.0/3.0)); } else if (fabs(Orbit.TrueDistance) >= FLT_EPSILON && fabs(Orbit.TrueMass ) >= FLT_EPSILON ) { Orbit.TruePeriod = sqrt((Orbit.TrueDistance*Orbit.TrueDistance) * Orbit.TrueDistance * 4.0 * (PI*PI) / (G_SI * Orbit.TrueMass)); } else { /* -------------------- default -------------------------------- */ Orbit.TruePeriod = 2.2226262e7; /* 0.7 year */ Orbit.TrueDistance = 1.496e11; /* one astronomical unit */ Orbit.TrueMass = 3.978e30; /* two solar mass */ /* M*G/(4*PI*PI) = D*D*D/P*P */ /* (Keplers Third Law) */ } } /****************************************************************** @package nightfall @author Rainer Wichmann (rwichman@lsw.uni-heidelberg.de) @version 1.0 @short Print usage info @param (void) @return (void) @heading Main *******************************************************************/ void UsageInfo() { /* -------- PACKAGE, VERSION are declared in the makefile --- */ printf("\n"); /*@i@*/ printf(_(" Welcome to the %10s (version %5s)"), PACKAGE, VERSION); printf(_(" Light Curve Synthesis Program\n")); printf(_("\tCopyright (C) 1998 Rainer Wichmann\n")); printf("\n"); printf(_(" Usage: nightfall \n")); printf(_("\t (example: nightfall -G 0.8 85 0.8 1.0 5500 5800)\n")); printf(_("\t (interactive: nightfall -U 0.8 85 0.8 1.0 5500 5800)\n")); printf(_("\t (config_file: nightfall -U -C ty_boo.cfg)\n")); printf("\n"); printf(_(" Basic System Parameters:\n")); printf(_("\tQ = Mass(Secondary)/Mass(Primary),\n")); printf(_("\tI = Inclination angle of orbital plane (degree)\n")); printf(_("\tRF1, T1 = Primary: Roche_fill_factor, Temperature\n")); printf(_("\tRF2, T2 = Secondary: Roche_fill_factor, Temperature\n")); printf("\n"); printf(_(" Options:\n")); printf(_("\t<-U> Interactive (GUI)\n")); printf(_("\t<-A> Animated view\n")); #ifdef _WITH_OPENGL printf(_("\t<-g> Don't use OpenGL for animation\n")); #endif printf(_("\t<-C config_file> Input configuration file \n")); printf(_("\t<-I datafile> Input observed data \n")); printf(_("\t<-G[P,S,1,2]> Graph of lightcurve (default: 1)\n")); printf(_("\t\tP,S: close-up of Primary/Secondary eclipse\n")); printf(_("\t\t1,2: display 1/2 orbital cycles\n")); printf(_("\t<-B[U/B/V...]> Bandpass to display in graph \n")); printf("\t (default: V)\n"); printf(_("\t<-V[v,i,c,a]> Visualize geometry (default: v)\n")); printf(_("\t\tv: view of stars, \n")); printf(_("\t\ti: image of potential,\n")); printf(_("\t\tc: contourmap of potential\n")); printf(_("\t\ta: all of the above\n")); printf(_("\t<-H> Hardcopy (postscript plot)\n")); printf("\n"); printf(_(" Advanced System Parameters:\n")); printf(_("\t<-O[P/S] value> log g (surface gravity)\n")); printf(_("\t<-f[P/S] value> asynchroneous rotation (Omega/Omega_System)\n")); printf(_("\t<-s[P/S] longitude latitude radius dimfactor> Spot on Primary/Secondary\n")); printf(_("\t<-e eccentricity periastron_length> eccentric orbit parameter\n")); printf(_("\t<-t[P/M/D] value> period/total mass/distance in days/solar masses/radii\n")); printf("\n"); printf(_(" Debugging Options:\n")); printf(_("\t<-D[vwb]> Debug [verbose,warn,busy]\n")); /********************************************** printf(" <-D[vwbdgflmrsecxiqoa]> Debug [verbose,warn,busy,\n"); printf(" divide,geometry,fractional\n"); printf(" lagrange,minfind,rootfind,spot,\n"); printf(" eclipsetest,computeflux,\n"); printf(" xcentric,infile\n"); printf(" chisquare,optimize,atmosphere]\n"); *********************************************/ printf("\n"); printf(_(" Data Fitting Options:\n")); printf(_("\t<-X[fit_params] Tolerance>\n")); printf(_("\t\t where 'fit_params' is a string for selecting\n")); printf(_("\t\t the parameters to fit, as follows:\n")); printf("\t\t0123456789: Q, i, RF1, RF2, T1, T2, e, w, F1, F2,\n"); printf(_("\t\tA-HI-P: 2 Spots (Primary), 2 Spots (Secondary),\n")); printf(_("\t\tQR: Mass, Distance,\n")); printf(_("\t\ta-l: Third Light (UBVRIJHKuvby) \n")); printf(_("\t\t(set Tolerance = 0.001 to switch on Simulated Annealing)\n")); printf(_("\t<-Y[as above] Step1 Step2> Chi Square Map (2 Parameters)\n")); printf("\n"); printf(_(" Computation Options:\n")); printf(_("\t<-Plamda_zero> Profile of absorption line at\n")); printf(_("\t rest wavelength lamda_zero (nm)\n")); printf(_("\t<-Nnn> nn steps for lightcurve \n")); printf("\t (default 80)\n"); printf(_("\t<-M> use Model atmosphere \n")); printf(_("\t<-F> compute Fractional visibility \n")); printf(_("\t<-L[0-2]> Limb darkening method \n")); printf("\t (default: linear == 0)\n"); printf(_("\t\t0: linear, 1: quadratic, 2: square root, 3: detailed\n")); printf(_("\t<-R[1-9]> Reflection treatment \n")); printf(_("\t (default: sloppy == 0)\n")); printf(_("\t\t0: point source, 1-9: exact, up to 9 iterations\n")); printf(_("\t<-a[P/S] albedo> Albedo des Primary/Secondary\n")); printf(_("\t<-3CM> Third light, C: colour code, M: relative brightness \n")); printf("\n"); printf(_("Compiled with:\n")); printf(_("\t Support for interactive usage: ")); #ifdef _WITH_GTK printf(_("yes\n")); #else printf(_("no\n")); #endif printf(_("\t Support for graphics: \t\t")); #ifdef _WITH_PGPLOT #ifdef _WITH_GNUPLOT printf(_("gnuplot\n")); #else printf(_("pgplot\n")); #endif #else printf(_("none\n")); #endif printf(_("\t Support for OPENGL: \t\t")); #ifdef _WITH_OPENGL printf(_("yes\n")); #else printf(_("no\n")); #endif printf(_("\t Support for MPI: \t\t")); #ifdef _WITH_MPI printf(_("yes\n")); #else printf(_("no\n")); #endif } /****************************************************************** @package nightfall @author Rainer Wichmann (rwichman@lsw.uni-heidelberg.de) @version 1.0 @short Initialize @param (void) @return (void) @heading Main *******************************************************************/ void InitFlags() { int j; Flags.Step[0] = 0.01; Flags.Step[1] = 0.01; #ifdef _WITH_OPENGL Flags.use_opengl = ON; Flags.labels = OFF; #endif Flags.visualize = OFF; for (j = 0; j < 128; ++j) Flags.debug[j] = OFF; Flags.interactive = OFF; Flags.animate = OFF; /* added disk Flag MK 13.05.01 */ Flags.disk = OFF; Flags.warp = OFF; Flags.plot = OFF; Flags.PlotBand = -1; /* signals that nothing read yet */ for (j = 0; j < (NUM_MAG+2); ++j) Flags.Passbands[j] = OFF; for (j = 0; j < 64; ++j) Flags.simplex[j] = OFF; Flags.fractional = OFF; Flags.lineprofile = OFF; Flags.WantFit = OFF; Flags.SimplexTol = SIMPLEX_FIT_TOLERANCE; Flags.anneal = OFF; Flags.limb = 0; Flags.reflect = 0; Flags.Spots1 = 0; Flags.Spots2 = 0; Flags.fill = OFF; Flags.elliptic = OFF; Flags.asynchron1 = OFF; Flags.asynchron2 = OFF; Flags.blackbody = ON; Flags.GL = OFF; #ifdef _WITH_OPENGL Flags.axes = OFF; Flags.movie = OFF; Flags.frame = 0; Flags.wireframe = ON; Flags.points = OFF; Flags.texture = OFF; Flags.textype = IMAGE; Flags.GLdistancez = 1.5; Flags.GLdistancex = 0.0; #endif Flags.first_pass = ON; Flags.ProComputed = OFF; Flags.ConfFile[0] = '\0'; strcpy (Orbit.Name, "(void)"); Orbit.Dist = 1.0; Orbit.MinDist = 1.0; Orbit.Excentricity = 0.0; Orbit.Omega = 90.; Orbit.TruePeriod = 0.0; Orbit.TrueDistance = 0.0; Orbit.TrueMass = 0.0; Orbit.LambdaZero = LAMDAZERO_DEFAULT; /* should be nanometer */ for (j = 0; j < NUM_MAG; ++j) Orbit.Third[j] = 0.0; Binary[Primary].Fratio = 1.0; Binary[Secondary].Fratio = 1.0; /* see Light.h for more details */ #ifdef HAVE_DISK Binary[Disk].Tilt = 0.0; Binary[Disk].Warp = 0.0; Binary[Disk].Thick = 0.01; Binary[Disk].HR = 0.0; /* 0 for a flat disk */ Binary[Disk].Rout = 0.8; Binary[Disk].Rin = 0.2; Binary[Disk].Temperature = 8000.0; Binary[Disk].tempHS = Binary[Disk].Temperature; Binary[Disk].longitudeHS = 0.0; Binary[Disk].extentHS = 0.0; Binary[Disk].depthHS = 0.0; #endif Flags.plotOpened = OFF; strcpy(Out_Plot_File, "Nightfall.ps"); if (getenv("NIGHTFALL_PLOTFILE") != NULL) { strncpy( Out_Plot_File, getenv("NIGHTFALL_PLOTFILE"), (size_t) ((long) sizeof(Out_Plot_File) - 1) ); Out_Plot_File[sizeof(Out_Plot_File) - 1] = '\0'; } return; }