Disclaimer: This code is distributed freely to the community. The user accepts sole responsibility for the results produced by the code. Although every effort has been made to identify and eliminate errors, we accept no responsibility for erroneous model predictions.

The Starburst99 code and the models should be quoted as Leitherer et al. (1999, ApJS, 123, 3) and Vazquez & Leitherer (2005; ApJ, 621, 695). If reference to the Hillier or Pauldrach model atmospheres is made, please quote Smith, Norris, & Crowther (2002).

Download the source code v5.1 and the ancillary files. This link has all the files needed to run the UNIX/Solaris version of Starburst99, which is currently the default.

We also provide a link to a LINUX version of the Starburst99 source code. This version corresponds to v5.0 of the UNIX version and was kindly provided by Lihong Yao (Univ. of Toronto). If you are interested in this version, you should first download the UNIX version with all the ancillary files, and then replace the source code and all compiler specific files by the contents of the LINUX link. To compile use "g77 -O2 -o galaxy_linux galaxy_linux.f". Note that v5.1 of the source code contains several fixes that caused errors in v5.0. Read the Knowledge Base article of April 10, 2006 to determine if this may affect your results when using the available LINUX version.

File structure
The package comes with quite a few files. Some are essential whereas others are not, and you can do without them.

1. The README file -- you should read this first.
2. galaxy.f -- this is the code. It is standard Fortran 77. Nothing fancy, and some software wizards may consider the coding pretty pedestrian. Keep in mind the code was written for ourselves and not as shareware for the community. The main goal was to make it easily understandable and to make the structure be driven by astrophysical requirements. Elegance and speed were not our prime concern.
3. Makefile -- this file is used to compile the code. You may want to tailor this to your need.
4. go_galaxy -- this is a script to run the code. Its most important purpose is to assign the directories where the auxiliary files reside. You may also want to modify this depending on your directory structure.
5. save_output -- this is a script to save the output files and give them reasonable names. (The input file with the model parameters is also included here.) To get started, we suggest to keep the names as they are in this file.
6. mod****.dat -- 20 files containing the evolutionary models. There are 4 groups of 5 files. 1) "c" - Geneva tracks with standard mass-loss rates and chemical composition of 2x solar (040), solar (020), 0.4x solar (008), 0.2x solar (004), and 0.05x solar (001); 2) "e" - same as "c" but with high mass-loss rates; 3) "s" - Padova tracks as updated by Girardi (2000) for chemical composition 0.050, 0.020, 0.008, 0.004, and 0.0004; 4) "p" - same as "s" but with the inclusion of thermally pulsing AGB stars following the prescription of Vassiliadis & Wood (1993).
7. lcb97_***.flu -- 5 files containing the model atmospheres of Lejeune et al. They match the metallicities of the evolutionary tracks. Note that some models were interpolated from the original Lejeune set since the required metallicities were not available.
8. wr_beta*.fluxes -- 2 files with the WR model atmospheres of Schmutz et al. (1992).
9. CMFGEN*.dat -- 10 files containing the WR models computed by the UCL group (5 metallicities, 2 WR types).
10. WMbasic*.dat -- 5 files with the Pauldrach O-star atmospheres for 5 metallicities.
11. allstars*.txt -- 9 files in support of the optical high-resulution spectra. The files contain the line spectra ("flux"), theoretical continua ("cont"), and the wavelength grid ("wave") for the theoretical library of Martins et al. (2005). p03, p00, m05, and m10 indicate metallicities of twice solar, solar, 40% solar, and 10% solar.
12. sp.dat -- IUE spectral library of O, B, and WR stars used in the subroutine linesyn. This library is also available (in a more readable form) via the CD accompanying the article in PASP, 108, 996 (1996), or from our web page.
13. sp_low.dat -- FOS and STIS spectral library of LMC/SMC O stars used in the subroutine linesyn. The structure is identical to that of sp.dat. All stars other than O stars are the same as in sp.dat (Leitherer et al. 2001; ApJ, 550, 724).
14. fuse_high.dat -- FUSE spectral library of Galactic O and early B stars used in the subroutine fusesyn. The structure is identical to that of sp.dat. All stars other than O and early B have their continua set to 1.
15. fuse_low.dat -- same as fuse_high.dat but for LMC and SMC stars. Molecular hydrogen lines were removed. The reference for both FUSE libraries is Robert et al. (2003).
16. schkal.dat -- spectral-type calibration of Schmidt-Kaler (1982). The table simply contains a list of log L (in solar luminosities) and log Teff (in K). Each line corresponds to a certain spectral type.
17. irfeatures.dat -- contains data for the near-IR CO features at 1.62 and 2.29 microns, and the silicon feature at 1.59 microns from Origlia, Moorwood, and Oliva (1993; A&A, 280, 536).
18. standard.input1 -- this is the input file specifying the model parameters. It is explained in more detail below.
19. standard.* -- results of a standard model run. These files are obtained by using the parameters in the galaxy.input file. They can be used for test runs when implementing the code.

How to run the code
Get organized first.

• Create a directory for the code and the unix scripts (files 1 -- 5).
• Create separate directories for the evolutionary tracks(6), for the model atmospheres (7 -- 11), and for the auxiliary files (12 -- 17).
• A separate directory for the input (18) and for the output (19) should also be created. The directory structure should correspond to what you have declared in 4 and 5. The output directory will be populated by files produced during a successful run. The source code should be pretty much machine independent, except for the names of the auxiliary files which are called in the code.
• Search for the string "*** THE FILE NAME IS LOCATION DEPENDENT ***" and modify accordingly. If you extract the code and all the files from the gzipped tar file, you will obtain the same filenames and folder structure as at STScI, and little renaming should be required.

The input
Once all the files are in place and the declarations are complete, the input parameters need to be specified. If this is your first attempt, we suggest to leave the parameters as they are. They produce reasonable results and should give you a first impression of what is in store. Once you have gained more experience and have become more adventurous, you can modify the parameters to suit your needs. These are the parameters to play with:

• MODEL DESIGNATION:                                           [NAME]
  standard  -- any identifier you want to assign to the model. You will find
               it in the header of each output file.                 
  

• CONTINUOUS STAR FORMATION (>0) OR FIXED MASS (<=0):          [ISF]
  -1  -- if this is a negative integer, star formation is instantaneous, 
         otherwise it is continuous.       
  

• TOTAL STELLAR MASS [106 SOLAR MASSES] IF 'FIXED MASS' IS CHOSEN: [TOMA]
  1.   -- this is the total stellar mass (spread out between the upper and lower
          cut-off mases). It is only used if an instantaneous burst is specified.
  

• SFR [SOLAR MASSES PER YEAR] IF 'CONT. SF' IS CHOSEN:         [SFR]
  1.   -- the star formation rate (only used for a continuous rate). The total
          accumulated mass is spread out between the upper and lower cut-off 
          masses. 
  

• NUMBER OF INTERVALS FOR THE IMF (KROUPA=2):                  [NINTERV]
  2    -- intervals of the multi-power-law IMF. If two intervals are specified, 
          the program expects two IMF exponents and three IMF boundaries in the
          next two input fields. Up to ten such intervals may be specified.
  

• IMF EXPONENT(S) (KROUPA=1.3,2.3):                            [XPONENT]
  1.3,2.3 -- one or more IMF exponents for a power-law can be specified. The 
             exponents refer to the individual power-law intervals, ordered 
             by increasing mass. 
  

• MASS BOUNDARIES FOR IMF (KROUPA=0.1,0.5,100) [SOLAR MASSES]: [XMASLIM]
  0.1,0.5,100. -- the boundaries of the IMF intervals corresponding to the
                  specified exponents.
  

• SUPERNOVA CUT-OFF MASS [SOLAR MASSES]:                       [SNCUT]
  8.0 -- stars with ZAMS masses of 8 M and higher form supernovae. This is the
         suggested standard value but can be modified if desired.
  

• BLACK HOLE CUT-OFF MASS [SOLAR MASSES]:                      [BHCUT]
  120. -- stars with ZAMS masses of 120 M and lower form supernovae. An 
          alternative scenario would be to let stars above a certain threshold
          form a black hole. For instance, BHCUT=40. results in SNe only from the
          mass range 40 to 8 M.
  

• METALLICITY + TRACKS:                                        [IZ]
  GENEVA STD: 11=0.001;  12=0.004; 13=0.008; 14=0.020; 15=0.040  
  GENEVA HIGH:21=0.001;  22=0.004; 23=0.008; 24=0.020; 25=0.040
  PADOVA STD: 31=0.0004; 32=0.004; 33=0.008; 34=0.020; 35=0.050
  PADOVA AGB: 41=0.0004; 42=0.004; 43=0.008; 44=0.020; 45=0.050
  24 -- this integer indicates the evolutionary tracks to be used. The
        choices are 11-15, 21-25, 31-35, and 41-45, where the numbers indicate
        the metallicity of the four sets of tracks that are available. Example: 
        "23" selects 40% solar metallicity from the Geneva high mass-loss 
        tracks.
  

• WIND MODEL (0: EVOLUTION; 1: EMP.; 2: THEOR.; 3: ELSON):        [IWIND]
  2 -- this selects the wind model to be used for the calculation of the
       wind power. The four models are discussed in ApJ, 401, 596 (1992). "0"
       is the suggested default parameter.
  

• INITIAL TIME [1.E6 YEARS]:                                   [TIME1]
  0.01 -- the epoch of the onset of the star formation. In almost all cases you
          want this to be close to 0. It should not be exactly 0 for numerical
          reasons. 0.01 (i.e. 10e4 yr) is a good number.
  

• TIME SCALE: LINEAR (=0) OR LOGARITHMIC (=1)                  [JTIME]
  0 --  a switch to select linear or logarithmic time intervals.
  

• TIME STEP [1.e6 YEARS] (ONLY USED IF JTIME=0):               [TBIV]  
  0.1 --  this is the timestep used for the calculations. It is a very important
          parameter. On the one hand, the computing time scales with STEP, so
          you want to avoid too high resolution, but on the other, short
          evolutionary phases can be missed. 0.1 (i.e. 10e5 yr) is a good 
          value if you use full isochrone synthesis. If full isochrone
          synthesis is not used, 0.1 or large is suggested only for tests
          --- be aware that WR or RSG numbers are no longer properly calculated
          for a STEP of 0.1 unless full ischrone synthesis is selected!       
  

• NUMBER OF STEPS        (ONLY USED IF JTIME=1):               [ITBIV]
  1000 -- if a logarithmic scaling is selected, the time step size varies with 
          time and is no longer specified via the TIME STEP field. In this case, 
          we enter the total number of time steps, which will then be distributed
          logarithmically between the first and the last time point. As before, 
          users should beware of too small steps during short-lived evolutionary 
          phases. 
  

• LAST GRID POINT [1.e6 YEARS]:                                [TMAX]
  100. -- the oldest age of the model. 
  

• SMALL (=0) OR LARGE (=1) MASS GRID;
  ISOCHRONE ON  LARGE GRID (=2) OR FULL ISOCHRONE (=3):        [JMG]
  3 -- these are four options for the interpolation in mass. They are explained
       in the code. Shortly: 0 -- evolutionary synthesis with a mass 
       resolution of 5 M (only recommended for tests); 1 -- same as 0, but with
       a resolution of 1 M. This method was used in Leitherer & Heckman (1995); 
       2 -- isochrone synthesis with a fixed mass resolution of 1 M; 
       3 -- isochrone synthesis with a variable mass grid. This is the recommended
       method. In particular, FULL ISOCHRONE must be used if masses below 1 solar 
       masses from the Padova tracks are to be included in the modeling. 
  

• LMIN, LMAX (ALL=0):                                          [LMIN,LMAX]
  0 -- LMIN and LMAX are the indices of the evolutionary tracks, sorted by mass.
       Normally you do not want to mess with the variable and leave it at 0. 
       However, if you want to track down some peculiarity of the output, you
       may want to compute the parameters for only one track. For instance, 
       specifying 21,21 indicates that only a 100 M star should be used, and
       everything else is suppressed. The cross-ID's between index and mass
       are at the bottom of the input file. The example here refers to JMG=1 or
       2. For JMG=0, you would have chosen 5,5. This does not apply to JMG=3
       since the mass grid is variable. If JMG=3, LMIN and LMAX are not used.
  

• TIME STEP FOR PRINTING OUT THE SYNTHETIC SPECTRA [1.e6YR]:   [TDEL]
  1.0 -- the file containing the output spectrum can be pretty big. This
         parameter controls the time step to print out the spectrum. This is 
         independent of the time resolution -- only the print out is affected!
         1 Myr is usually a good value but if you compute the starburst up to
         100 Myr, you may prefer TDEL=5 Myr unless you have many Mb of disk
         space.
  

• ATMOSPHERE FOR THE LOW-RES SPECTRUM: 1=PLA, 2=LEJ, 3=LEJ+SCH, 4=LEJ+HIL, 5=PAU+HIL [IATMOS]
  5 -- this is the choice of the model atmosphere. 1 is a bare-bone version with
       black bodies, good only for tests. 2 uses the Kurucz models as compiled
       by Lejeune for all stars. 3 uses Lejeune for stars with plane-parallel
       atmospheres and Schmutz for stars with strong winds. 4 uses Lejeune, 
       but replaces the Schmutz by the Hillier atmospheres. 5 is like 4, except
       for the O atmospheres, for which we use the Pauldrach models. 5 is the 
       recommended value. 
  

• METALLICITY OF THE HIGH RESOLUTION MODELS                    [ILIB]
  (1=0.001, 2= 0.008, 3=0.020, 4=0.040):
  3 -- a switch to choose the metallicity of the optical high-res spectra. 4
       choices are offered, and the user can decide how to match them to the    
       evolution models.
  

• LIBRARY FOR THE UV LINE SPECTRUM: (1=SOLAR, 2=LMC/SMC)       [ILINE]
  1 -- a switch for the choice of the UV spectral library. This switch applies
       to both the FUSE and the HST/IUE libraries. It is independent of the 
       metallicity of the tracks/atmospheres. Normally one would use
       ILINE=1 with IZ=24 and ILINE=2 with IZ=22.
  

• RSG FEATURE: MICROTURB. VEL (1-6), SOL/NON-SOL ABUND (0,1)   [IVT,IRSG]
  3,0 -- atmospheric parameters used for the spectral features in the near-IR.
         Detailed explanations are in the sp-feature subroutine. Defaults are
         3,0, i.e. microturbulent velocities of 3 km/sec and solar abundance
         ratios for alpha-element/Fe.
  

• OUTPUT FILES (NO<0, YES>=0)                                  [IO1,...]
  +1,+1,-1,+1,+1,+1,+1,+1,+1,+1,+1,+1,+1,+1,-1
  	These are options to generate various outputs. We recommend to use 
          the default setting for the flags, at least until you become more 
          familiar with the code. Some of the subroutines are interrelated. If 
          you choose such a subroutine but not the other, required one, a warning 
          will be issued. The 15 output flags are explained in the next section.
          They are discussed in the order as they appear above. For instance,
          "(7)" refers to the 7th of the 15 flags.
  

The output

1. Computation of the number of ionizing photons. (7) must be set to "yes" since the spectrum below 912 A is needed. Output is the number of ionizing photons in the HI, HeI, and HeII continuum, their fractions relative to the total luminosity, and the total luminosity.
2. Calculation of the supernova rate and the mechanical luminosities (winds and supernovae). It requires (4) to obtain the stellar wind luminosities. Otherwise it is independent of other subroutines.
3. HRD with a few evolutionary tracks. This is mostly useful for test purposes. This part is independent of all other subroutines and can be turned on/off without doing any harm. The HRD cannot be produced in Isochrone Synthesis mode, therefore it is turned off by default.
4. Mechanical luminosity and related quantities due to winds. (No supernovae!) It does not depend on any other subroutine since no information on the energy distribution is needed.
5. Two output files containing the stellar spectral types during each time step and the relative numbers of WR stars. The spectral types follow the scheme by Schmidt-Kaler, oversampled by a factor of 2. For instance, there are 18 entries for spectral type B. They are the number of stars for types B0, B0.5, B1,...B9.5 (total of 18). Schmidt-Kaler's table has B0, B1,....B9 (total of 9). The spectral types are printed out only every TDEL. Otherwise it is too bulky.
6. The mass in individual elements released via stellar winds and supernovae. No other subroutines are needed. The supernova yields for type II supernovae are taken into account.
7. The continuous spectrum of the stellar population for each time step. The columns are time, wavelength, stellar+nebular, stellar only, and nebular only fluxes. (1) is needed in order to calculate the nebular continuum.
8. The ultraviolet line spectrum at 0.75 A resolution from 1200 to 1600 A (LMC/SMC library) or to 1800 A (Milky Way library). The subroutine needs (7) to compute the stellar continuum and (1) for the nebular continuum. If (1) is turned off, the nebular contribution can not be added (it is often small, though). The columns have time, wavelength, absolute luminosity, and rectified (continuum=1) luminosity.
9. Calculation of colors and magnitudes. The subroutine needs (7) to compute the stellar continuum and (1) for the nebular continuum. If (1) is turned off, the nebular contribution can not be added and the computed colors are for stars only (this may sometimes be desirable). The filter system is defined in the code.
10. Calculation of the strengths of H_alpha, H_beta, Pa_beta, and Br_gamma. For each line we give the continuum luminosity, the line luminosity, and the equivalent width (everything logarithmic). The subroutine needs (7) to compute the stellar continuum and (1) for the number of ionizing photons.
11. Calculation of the strengths of various IR spectral features. First is the CO index as computed by Doyon et al. (1994, ApJ, 421, 101). (Please note that this calculation has no metallicity dependence. A later version of this routine will compute the CO index using the model atmospheres themselves and give metallicity-dependent results.) Next are two computations of the CaII IR triplet using the relations of Diaz et al. (1989, MNRAS, 239, 325). The relations from Diaz et al. have no temperature dependence; the first calculation has the feature present in stars of all temperatures; the second has the index set to zero strength for stars with T>7200K (spectral type A or earlier). Next come the 1.62 and 2.29 micron CO features, and the 1.59 micron Si feature, which were modelled for individual stars by Origlia et al. (1993, A&A, 280, 536.) The indices can be computed for solar [Si/Fe] and [C/Fe], or a model with enhanced [Si/Fe] and depleted [C/Fe] (as for young systems enriched primarily by Type II SNe), and for stellar atmospheric microturbulent velocities (MTVs) of 1-6 km/s. (Note that the changes to the abundance ratios and MTVs are self-contained in this routine and have no effect upon the other outputs, e.g., colors, of the code.)
12. This subroutine is equivalent to (8), but it computes the spectral region between 1000 and 1180 A.
13. Calculates fully theoretical spectra between 3000 and 7000 A at 0.3 A resolution. These spectra are independent of the calculation of the low-resolution spectra in output 7. The file structure is the same as in (8).
14. Calculation of the most important WR emission lines using the line luminosities of Schaerer & Vacca (1998, ApJ, 497, 658). These are only those lines originating in WR winds --- not the nebular lines in the HII region. Quantities given are the line fluxes and the equivalent widths. The subroutine needs (7) to compute the stellar continuum and (1) for the nebular continuum. If (1) is turned off, the nebular contribution can not be added.
15. This subroutine will be supported in a future realease. Eventually it will compute a high-resolution UV line spectrum from model atmosheres, as opposed to using an empirical library. We are currently testing these models and find them not quite ready for general release. A pilot project is discussed in Rix et al. (2004; ApJ, 615, 98). Stay tuned!