{ Cartesian Stardata Generator                                       }
{ Turbo Pascal 5.0                                                   }
{ This program accepts as input the Gliese Near Star Catalog         }
{ then allows the user to input filter criteria to exclude certain   }
{ types of stars.                                                    }
{ Please note that it works with the 2.0 catalog, not the 3.0.}
{ The program then produces a text file containing the stars that    }
{ pass the filter, and the following information:                    }
{ Gliese catalog number, Star name, Spectral class, distance from    }
{ the sun in parsecs, cartesian x,y,z co-ordinates with the x-y      }
{ plane being the celestial equator, cartesian xg,yg,zg co-ordinates }
{ with the x-y plane being the galactic plane, and the center of the }
{ galaxy being about 10 kiloparsecs down the positive x-axis,        }
{ and star's luminosity.                                             }
{ The Gliese near star catalog should be available at the following  }
{ internet sites:                                                    }
{ http://cdsweb.u-strasbg.fr/cats.html (query the database for Gliese) }
{ by ftp at cdsarc.u-strasbg.fr in /pub/cats                           }
{ http://nssdc.gsfc.nasa.gov/adc/adc.html                              }
{ Program written by Winchell Chung, Feb 20, 1995                      }
{ This program is in the public domain                                 }
{+=======================(:)=== ^ ===(:)======================================+}
{|    WINCHELL CHUNG     |=|   /_\   |=| I'm nobody. Nobody at all. But the   |}
{|Nyrath the nearly wise |\|  <(*)>  |\| secrets of the universe don't mind.  |}
{|  nyrath@clark.net     |=| /_/|\_\ |=| They reveal themselves to nobodies   |}
{| winchchung@delphi.com |\|  //|\\  |\| that care. OUTER LIMITS: Galaxy Being|}
{+=======================(:)=///|\\\=(:)======================================+}

program StarMap;
uses
   crt,dos;
Const
   degreesToRadians = 0.0174532925;
Type
   KeyListType = string[10];
Var
   starfilespec: string[80];
   starDataFile: text;
   resultfilespec: string[80];
   resultFile: text;
   iocode: integer;
      gnumS:          string[3];   { Gliese catalog number }
      catCodeS:       string[1];
      componentCodeS: string[1];   { multiple star designation }
      closeCodeS:     string[1];
      starNameS:      string[12];  { star name }
      hdNumS:         string[6];
      raHoursS:       string[2];   { Right Ascension hours }
      raMinutesS:     string[2];   { Right Ascension minutes }
      raSecondsS:     string[2];   { Right Ascension seconds }
      decSignS:       string[1];
      decDegreesS:    string[2];   { format F3.1 }
      blank1S:        string[1];
      decMinutesS:    string[3];
      blank2S:        string[1];
      pmRAS:          string[5];
      pmDecS:         string[5];
      rvS:            string[4];
      wdRVCodeS:      string[1];
      rvQS:           string[1];
      luminClassS:    string[1];
      spectralClassS: string[3];
      mkLuminClassS:  string[1];
      spectralPeculS: string[1];
      spectralCompS:  string[1];
      spectralSrcS:   string[1];
      blank3S:        string[1];
      appMagS:        string[4];
      ipvCodeS:       string[1];
      bvColorS:       string[4];
      bvJointS:       string[1];
      ubColorS:       string[4];
      ubJointS:       string[1];
      riColorS:       string[4];
      riJointS:       string[1];
      trigParallaxS:  string[4];
      paraErrorS:     string[3];
      absMag1S:       string[4];
      absMagErrorS:   string[1];
      spectParallaxS: string[4];
      spectParErrorS: string[1];
      photoParallaxS: string[4];
      photoParErrorS: string[1];
      bestParallaxS:  string[4];
      bestParaErrorS: string[3];
      blank4S:        string[1];
      absMag2S:       string[4];
      absMagPhotoS:   string[1];
      absMagError2S:  string[1];
      uSpaceVelS:     string[4];
      vSpaceVelS:     string[4];
      wSpaceVelS:     string[4];
   RA:   real;
   RAHours, RAminutes, RAseconds: real;
   DECDegrees, DECminutes: real;
   DEC:  real;
   parallax: real;
   dist:     real;
   x,y,z:    real;
   xg,yg,zg: real;
   phi,theta,rho,rvect: real;
   theError: integer;
   a: integer;
   luminclass: integer;
   mkLuminClass: integer;
   spectralPecul:integer;
   tabbing:integer;
   absMag2:  real;
   lum: real;
   minDist, maxDist, minLum, maxLum: real;
   noBinaries: boolean;
   theIndex: integer;
   thisBinaryFlag: boolean;
{ ***************************** }
{ raises NUMBER to the power of EXPONENT }
function Power(Number, Exponent: real): real;
begin
   If Number > 0.0 then Power := exp(Exponent * ln(Number))
   else                 Power := 0.0
end;

{ ***************************** }
{ Code = -2  Only ENTER key was input }
{      = -1  Nonnumeric key input, value is in CharFlag }
{      =  0  Real Number input, value is in Number }
{      >  1  Illegal entry. Code contains position of first illegal char }
procedure GetNumR(var Number: real; var CharFlag: char;
                  var Code: integer);
const
   PromptString = 'Entry?';
var
   Entry: string[30];
begin
   write(PromptString);
   readln(Entry);
   val(Entry, Number, Code);
   case length(Entry) of
      0: Code := -2; { only ENTER key was input }
      1: if Code > 0 then begin
         Code := -1;
         CharFlag := Entry[1];
         end;
      end;
end;

{ ***************************** }
{ Gets single character user reply from the keyboard, }
{ limiting it to characters in the list KeyList       }
{ KeyListType must be defined in Types above          }
procedure GetKey(KeyList: KeyListType; var Reply, Reply2: char);
begin
   Reply2 := chr(0);
   repeat
      Reply := readkey;
      if (Reply = #27) and (keypressed) then Reply2 := readkey;
      if length(KeyList) = 0 then exit;
      until pos(Reply, KeyList) > 0
end;


{ ***************************** }
{ Gets single character user reply from the keyboard, }
{ limiting it to range FIRST to LAST                  }
procedure GetReply(First, last:char; var reply:char);
var
   temp: char;
begin
   if first > last then begin
      temp := first;
      first := last;
      last := temp;
      end;
   writeln('Enter reply from ',first,' to ',last,'.');
   repeat
      reply := readkey;
      if (reply = #27) and (keypressed) then begin
         reply := readkey;
         reply := #0
         end
      until (reply >= first) and (reply <= last)
end;

{ ***************************** }
procedure PrintMenu;
begin
   ClrScr;
   writeln('[][][][] FILTERING PARAMETERS [][][][]');
   writeln('only include stars that fit the following:');
   writeln;
   writeln('[1] Minimum distance from the sun = ',minDist:5:2,' parsecs');
   writeln('[2] Maximum distance from the sun = ',maxDist:5:2,' parsecs');
   writeln('[3] Minimum luminosity = ',minLum:3:1,' (sun=1.0)');
   writeln('[4] Maximum luminosity = ',maxLum:3:1);
   if (nobinaries = true) then
      writeln('[5] *NO* binaries and trinaries allowed')
   else
      writeln('[5] Binaries and trinaries allowed');
   writeln('[6] All Done');
   writeln;
end;

{ ***************************** }
function GetUserInput: boolean;
var
   reply, reply2: char;
   alldone: boolean;
   Number: real;
   CharFlag: char;
   Code: integer;
begin
   alldone := false;
   PrintMenu;
   GetReply('1','6',reply);
   case reply of
      '1': begin { minimum sun distance }
         writeln('input new minimum distance from the Sun in parsecs');
         GetNumR(Number, CharFlag, Code);
         if Code = 0 then begin
            if Number > maxDist then maxDist := Number;
            minDist := Number;
            end;
         end;
      '2': begin { maximum sun distance }
         writeln('input new maximum distance from the Sun in parsecs');
         GetNumR(Number, CharFlag, Code);
         if Code = 0 then begin
            if Number < minDist then minDist := Number;
            maxDist := Number;
            end;
         end;
      '3': begin { minimum luminosity }
         writeln('input new minimum luminosity (sun=1.0)');
         GetNumR(Number, CharFlag, Code);
         if Code = 0 then begin
            if Number > maxLum then maxLum := Number;
            minLum := Number;
            end;
         end;
      '4': begin { maximum luminosity }
         writeln('input new maximum luminosity (sun=1.0)');
         GetNumR(Number, CharFlag, Code);
         if Code = 0 then begin
            if Number < minLum then minLum := Number;
            maxLum := Number;
            end;
         end;
      '5': begin { no binaries allowed? }
         writeln('Are *no* binaries and trinaries to be allowed?');
         writeln('Press Y for yes or N for no');
         GetKey('YyNn', Reply, Reply2);
         writeln(Reply);
         if (Reply = 'y') or (Reply = 'Y') then noBinaries := TRUE;
         if (Reply = 'n') or (Reply = 'N') then noBinaries := FALSE;
         end;
      '6': begin { all done }
         alldone := true;
         end;
      end;
   GetUserInput := alldone;
end;

{ ***************************** }
procedure InitSettings;
begin
   minDist := 0.0;
   maxDist := 16.5;
   minLum := 0.4;
   maxLum := 2.0;
   noBinaries := TRUE;
end;

{ ***************************** }
procedure GetSettings;
begin
   InitSettings;
   While (GetUserInput = false) Do
   begin
   end;
end;


{ ***************************** }
   { MAIN PROGRAM }

BEGIN

   GetSettings;

   clrscr;

   Writeln('Enter filespec of Gliese input data file (GLIESE2.DAT)');
   readln(starfilespec);
   {$I-}
   Assign(starDataFile, starfilespec);
   Reset(starDataFile);
   iocode := ioresult;
   if iocode <> 0 then begin
      writeln(chr(7));
      writeln('file not found');
      exit;
      end;
   {$I+}

   Writeln('Enter filespec of output data file');
   readln(resultfilespec);
   {$I-}
   Assign(resultFile, resultfilespec);
   rewrite(resultfile);
   iocode := ioresult;
   if iocode <> 0 then begin
      writeln(chr(7));
      writeln('file not found');
      exit;
      end;
   {$I+}

   write(resultfile,'Stars ');
   if (noBinaries = TRUE) then
   write(resultfile,'with no binaries, ');
   writeln(resultfile,'farther than ',minDist:5:2,'pc and closer than ',maxDist:5:2,'pc,');
   writeln(resultfile,'Lum between ',minLum:3:1,' and ',maxLum:3:1);
   writeln(resultfile,'Gliese StarName     Spec  DistPC   X      Y      Z      Xg     Yg     Zg  Lum');

   theIndex := 1;
   while not eof(starDataFile) do begin

      Readln(starDataFile,gnumS,catCodeS,componentCodeS,closeCodeS,starNameS,
         hdNumS,raHoursS,raMinutesS,raSecondsS,decSignS,decDegreesS,blank1S,
         decMinutesS,blank2S,pmRAS,pmDecS,rvS,wdRVCodeS,rvQS,luminClassS,
         spectralClassS,mkLuminClassS,spectralPeculS,spectralCompS,spectralSrcS,
         blank3S,appMagS,ipvCodeS,bvColorS,bvJointS,ubColorS,ubJointS,riColorS,
         riJointS,trigParallaxS,paraErrorS,absMag1S,absMagErrorS,spectParallaxS,
         spectParErrorS,photoParallaxS,photoParErrorS,bestParallaxS,bestParaErrorS,
         blank4S,absMag2S,absMagPhotoS,absMagError2S,uSpaceVelS,vSpaceVelS,wSpaceVelS);

      Writeln('[',theIndex,'] ',starNames);
      theIndex := theIndex + 1;

      Val(bestparallaxs, parallax, theError);
      parallax := parallax / 1000;
      if parallax = 0 then parallax := 1;
      dist := 1/parallax;

      Val(absMag2S, absMag2, theError);
      absMag2 := absMag2/100;
      lum := Power(10,(0.4 * (4.9 - absMag2)));

      thisBinaryFlag := TRUE;
      if (noBinaries = TRUE) and (componentCodes <> ' ') then
         thisBinaryFlag := FALSE;
      if  ((dist < maxDist) and (dist > minDist)) and
          (thisBinaryFlag = TRUE) and
          ((lum > minLum) and (lum < maxLum)) then begin

         if catCodes ='0' then
            Write(resultfile, gnums,componentCodes,'  ')
         else
            Write(resultfile, gnums,'.',catCodes,componentCodes);
         Write(resultfile, ' ',starNames);
         Write(resultfile, ' ',Copy(spectralclasss,1,2));
         tabbing := 5;
         Val(mkluminClasss, mkluminclass, theError);
         case mkluminClass of
            1: begin Write(resultfile, 'I+'); tabbing := tabbing - 2; end;
            2: begin Write(resultfile, 'III'); tabbing := tabbing - 3; end;
            3: begin Write(resultfile, '-IV'); tabbing := tabbing - 3; end;
            4: begin write(resultfile, 'IV'); tabbing := tabbing - 2; end;
            5: begin write(resultfile, 'IV+'); tabbing := tabbing - 3; end;
            6: begin write(resultfile, 'IV'); tabbing := tabbing - 2; end;
            7: begin Write(resultfile, 'V'); tabbing := tabbing - 1; end;
            8: begin Write(resultfile, 'V+'); tabbing := tabbing - 2; end;
            9: begin write(resultfile, 'VI'); tabbing := tabbing - 2; end;
            else begin
               Val(luminClasss, luminclass, theError);
               case luminClass of
                  3: begin Write(resultfile, 'g'); tabbing := tabbing - 1; end;
                  4: begin write(resultfile, 'sg'); tabbing := tabbing - 2; end;
                  5: begin write(resultfile, 'd'); tabbing := tabbing - 1; end;
                  6: begin write(resultfile, 'sd'); tabbing := tabbing - 2; end;
                  7: begin write(resultfile, 'wD'); tabbing := tabbing - 2; end;
                  end;
               end;
            end;
         Val(spectralPeculs, spectralPecul, theError);
         case spectralPecul of
            1: begin Write(resultfile, 'p'); tabbing := tabbing - 1; end;
            2: begin Write(resultfile, 'e'); tabbing := tabbing - 1; end;
            3: begin Write(resultfile, 'm'); tabbing := tabbing - 1; end;
            4: begin write(resultfile, 'n'); tabbing := tabbing - 1; end;
            5: begin write(resultfile, 's'); tabbing := tabbing - 1; end;
            6: begin write(resultfile, 'ss'); tabbing := tabbing - 2; end;
            7: begin write(resultfile, 'ep'); tabbing := tabbing - 2; end;
            9: begin Write(resultfile, 'wk'); tabbing := tabbing - 2; end;
            end;
         case tabbing of
            1: Write(resultfile,' ');
            2: Write(resultfile,'  ');
            3: Write(resultfile,'   ');
            4: Write(resultfile,'    ');
            5: Write(resultfile,'     ');
            6: Write(resultfile,'      ');
            end;

         Val(rahourss, rahours, theError);
         Val(raminutess, raminutes, theError);
         Val(rasecondss, raseconds, theError);
         ra := (rahours*15) + (raminutes*0.25) + (raseconds*0.004166666);

         Val(decdegreess, decdegrees, theError);
         Val(decminutess, decminutes, theError);
         decminutes := decminutes/10;
         dec := (decdegrees) + (decminutes / 60);
         if (decsigns = '-') then dec := dec * -1;


         if (dist = 1.0) then write(resultfile,'*****')
         else                 write(resultfile, dist:5:2);

         phi := ra * degreesToRadians;
         theta := dec * degreesToRadians;
         rho := dist;

         rvect := rho * Cos(theta);
         x := rvect * Cos(phi);
         y := rvect * Sin(phi);
         z := rho * Sin(theta);

         xg := -(0.0672 * x) - (0.8727 * y) - (0.4835 * z);
         yg :=  (0.4927 * x) - (0.4504 * y) + (0.7445 * z);
		 zg := -(0.8676 * x) - (0.1884 * y) + (0.4602 * z);

         if (dist = 1.0) then begin
            write(resultfile,'  0.00,  0.00,  0.00|  0.00,  0.00,  0.00');
            end
         else begin
            Write(resultfile, x:6:2,',',y:6:2,',',z:6:2);
            Write(resultfile, '|',xg:6:2,',',yg:6:2,',',zg:6:2);
            end;
         writeln(resultfile,' ',lum:4:2);

         end;
      end;

   Close(starDataFile);
   Close(resultfile);
END.

