diff --git a/Gtools/gscripts/ghelp b/Gtools/gscripts/ghelp index 5d3bde8..648119f 100755 --- a/Gtools/gscripts/ghelp +++ b/Gtools/gscripts/ghelp @@ -1,268 +1,272 @@ #!/usr/bin/env python print """ ################################################ gscripts/ ################################################ ################################################ # utilities gatime2Myr --param params --a1 0.01 --a2 0.02 gatime2Myr --param params --a 0.1914 --da 7.56271e-09 gatime2Myr --param params --a 0.1914 --Dloga 3.95127e-08 ################################################ # run preparation galloc : determine the memory allocation for gadget gmktimesteps : set output timesteps for cosmological simulations gmktimesteps -p 10 > outputtimesteps.txt gread_header : return some info on the header (incomplete, actually return the atime) gread_header snap.dat gheader : return some info on the header gheader snap.dat addgas : add gas to a cosmological box ################################################ # global statitics on the run gcpu : plot cpu info gcpu cpu.txt gcpu -o budget --legend cpu.txt gcpu -o gravity --legend cpu.txt gcpu -o hydro --legend cpu.txt genergy : plot info relative to the energy genergy --relative energy.txt genergy -o EnergyKineticFeedback energy.txt gplot_energy : plot info relative to the energy budget gplot_energy --legend energy.txt gsteps : plot info relative to the steps gsteps info.txt gsteps -rf info.txt gsteps --param params --log y --rf 100 info.txt ################################################ pscripts/ ################################################ ################################################ # global output analysis # star formation rate and stellar mass (sfr.txt) pSfrvsTime --param params -o SFR sfr.txt pSfrvsTime --param params -o MSTARS sfr.txt pSfrvsTime --param params -o MNEWSTARS sfr.txt pSfrvsTime --param params -o sfr sfr.txt pSfrvsTime --param params -o mnewstars --integrate --derive --rf 100 sfr.txt pSfrvsTime --param params -o mstars sfr.txt # supernovae (chimie.txt) gSNvsTime -o NSNII,NSNIa --rf 1000 chimie.txt gSNvsTime -o EgySNThermal,EgySNKinetic --rf 1000 chimie.txt gSNvsTime -o EgySNThermal,EgySNKinetic --integrate chimie.txt gSNvsTime -o Nflag --rf 1000 chimie.txt gSNvsTime -o Nflag --integrate chimie.txt # accretion (accretion.txt) gAccretionvsTime -o MT accretion.txt gAccretionvsTime -o MT_th accretion.txt gAccretionvsTime -o DMDT accretion.txt gAccretionvsTime -o DMDT_th accretion.txt ################################################ # unit conversion # change units, including the little h (Hubble parameter) if it is different than 1 gchangeunits --param1 params --param2 params.dSph snap1.dat -o snap2.dat # correct from scaling factor (and add the correct time in code unit) gcosmo2giso --param1 params --param2 params.dSph snap2.dat -o snap3.dat +################################################ +# Jeans Informations + +gJeans --Density 1e3 --Temperature 10 ################################################ # snapshot analysis # luminosity pget_Luminosity --param params snap.dat pget_Luminosity --param params --tnow 3000 snap.dat # virial radius pget_VirialRadius --param params -dX 200 snap.dat # star formation rate pSFRvsTime --param params -o SFR --xmin=0 --xmax=14 snap.dat # mass vs time pSFRvsTime --param params -o MSTARS --xmin=0 --xmax=14 snap.dat pSFRvsTime --param params -o MSTARSz --xmin=0 --xmax=10 snap.dat pSFRvsTime --param params -o MSTARSa --xmin=0 --xmax=1 snap.dat # logRho - logT pplot --param params --select=gas --x logrho --y logT --xmin -3 --xmax 0 --ymin 1 --ymax 5 snap.dat pplot --param params --select=gas --x logrho --y logT --xmin -3 --xmax 0 --ymin 1 --ymax 5 --map snap.dat pplot --param params --select=gas --colorbar --log z --x logrho --y logT --z Tcool --xmin -3 --xmax 0 --ymin 1 --ymax 5 --zmin 1 --zmax 1e16 snap.dat # rho - T pplot --param params --select=gas --x rho --y T --log xy --xmin 1e-3 --xmax 1e0 --ymin 1e1 --ymax 1e5 snap.dat # R - T pplot --param params --select=gas --x R --y logT snap.dat # R - Fe pplot --param params --select=gas --x R --y Fe snap.dat # MgFe - Fe pplot --param params --select=stars1 --x Fe --y MgFe --xmin=-4 --xmax=0.5 --ymin=-1 --ymax=1.5 snap.dat pplot --param params --select=stars1 --x Fe --y MgFe --xmin=-4 --xmax=0.5 --ymin=-1 --ymax=1.5 --map snap.dat # Fe - Age pplot --param params --select=stars1 --x Age --y Fe --xmin=0 --xmax=14 --ymin=-4 --ymax=0.5 snap.dat # Tcool pplot --param params --select=gas --x logrho --y Tcool --xmin=-5 --xmax=2 --ymin 0 --ymax 1e16 --log y snap.dat ################################ # metallicite ################################ # MgFevsFe pMgFevsFe --select=stars1 --xmin=-4 --xmax=0.5 --ymin=-1 --ymax=1.5 --map snap.dat # BaFevsFe pBaFevsFe --select=stars1 --xmin=-4 --xmax=0.5 --ymin=-1 --ymax=1.5 --map snap.dat # [E1/E2] vs [E3/E4] pplot --select=stars1 --xmin=-4 --xmax=0.5 --ymin=-1 --ymax=1.5 --map --x Fe/H --y Ba/Fe snap.dat # NvsAge pNvsX --param params --select stars1 --x Age --nx 14 --xmin=0 --xmax=14 snap.dat # NvsFe pNvsX --param params --select stars1 --x Fe --nx 40 --xmin=-4 --xmax=0.5 --ymin=0 --ymax=0.15 snap.dat ################################ # profiles ################################ # grille spherique # density pSphericalProfile -o density --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # masse pSphericalProfile -o mass --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # integrated mass pSphericalProfile -o imass --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # velocity disp. in z pSphericalProfile -o sigmaz --xmax 10 --rmax 10 --nr 64 --center histocenter --param params snap.dat # circular velocity pSphericalProfile -o vcirc --xmax 10 --rmax 10 --nr 64 --center histocenter --param params --eps 0.1 snap.dat # grille cylindrique # density pCylindricalProfile -o sdens --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # masse pCylindricalProfile -o mass --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # integrated mass pCylindricalProfile -o imass --xmax 10 --rmax 10 --nr 64 --center histocenter --log xy --param params snap.dat # velocity disp. in z pCylindricalProfile -o sigmaz --xmax 10 --rmax 10 --nr 64 --center histocenter --param params snap.dat # circular velocity pCylindricalProfile -o vcirc --xmax 10 --rmax 10 --nr 64 --center histocenter --param params --eps 0.1 snap.dat # luminosity profile pCylindricalProfile -o Lum --xmax 10 --rmax 10 --nr 64 --center histocenter --select stars1 --param params snap.dat ############################################################ # velocities and velocity dispertions in cylindrical coords ############################################################ # extract all gExtractCylindrical --Rmax 10 --components="('gas','stars1','halo1')" --disk="('gas','stars1')" --nR 32 --eps=0.01 --nmin 3 --param params snap.dat # number of points per bins pplotdmp --mode=Numbers gas.dmp # circular velocity per component pplotdmp --mode=Velocities --legend gas.dmp stars1.dmp halo1.dmp # velocities and velocity disp pplotdmp --mode=Dispersions --legend gas.dmp # stability Toomre paramter Q pplotdmp --mode=Q gas.dmp """ diff --git a/Gtools/pscripts/pSFRvsTime b/Gtools/pscripts/pSFRvsTime index a33026c..f8ae007 100755 --- a/Gtools/pscripts/pSFRvsTime +++ b/Gtools/pscripts/pSFRvsTime @@ -1,353 +1,362 @@ #!/usr/bin/python from optparse import OptionParser import Ptools as pt from pNbody import * from pNbody import units from pNbody import ctes from pNbody import thermodyn import string from scipy import optimize def parse_options(): usage = "usage: %prog [options] file" parser = OptionParser(usage=usage) parser = pt.add_postscript_options(parser) parser = pt.add_ftype_options(parser) parser = pt.add_reduc_options(parser) parser = pt.add_center_options(parser) parser = pt.add_select_options(parser) parser = pt.add_cmd_options(parser) parser = pt.add_display_options(parser) parser = pt.add_info_options(parser) parser = pt.add_limits_options(parser) parser = pt.add_log_options(parser) parse = pt.add_units_options(parser) parser.add_option("-o", action="store", dest="obs", type="string", default = 'SFR', help="observable name", metavar=" NAME") parser.add_option("--x", action="store", dest="x", type="string", default = 'r', help="x value to plot", metavar=" STRING") parser.add_option("--nx", action="store", type="int", dest="nx", default = 500, help="number of bins") parser.add_option("--legend", action="store_true", dest="legend", default = False, help="add a legend") parser.add_option("--forceComovingIntegrationOn", action="store_true", dest="forceComovingIntegrationOn", default = False, help="force the model to be in in comoving integration") + parser.add_option("--forceComovingIntegrationOff", + action="store_true", + dest="forceComovingIntegrationOff", + default = False, + help="force the model not to be in in comoving integration") + parser.add_option("--curve", action="store_true", dest="curve", default = False, help="draw curve instead of histogram") (options, args) = parser.parse_args() pt.check_files_number(args) files = args return files,options def get_value_and_label(nb,val): if val == 'Age': label = r'$\rm{Age}\,\left[ \rm{Gyr} \right]$' out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Gyr,units.Unit_K]) val = nb.StellarAge(units=out_units.UnitTime) return val,label elif val == 'CosmicTime': label = r'$\rm{Age}\,\left[ \rm{Gyr} \right]$' out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Gyr,units.Unit_K]) val = nb.CosmicTime(units=out_units.UnitTime) return val,label elif val == 'Redshift': label = r'$\rm{Redshift}$' out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Gyr,units.Unit_K]) val = nb.Redshift() return val,label elif val == 'ScalingFactor': label = r'$\rm{a}$' val = nb.tstar return val,label else: label = r'%s'%val print "val = %s"%val exec("val = %s"%val) return val, label ####################################### # MakePlot ####################################### def MakePlot(files,opt): # some inits colors = pt.Colors(n=len(files)) datas = [] # read files for file in files: nb = Nbody(file,ftype=opt.ftype) ################ # units ################ # define local units unit_params = pt.do_units_options(opt) nb.set_local_system_of_units(params=unit_params) # define output units # nb.ToPhysicalUnits() if opt.forceComovingIntegrationOn: nb.setComovingIntegrationOn() + if opt.forceComovingIntegrationOff: + nb.setComovingIntegrationOff() + ################ # apply options ################ nb = nb.select(1) nb = pt.do_reduc_options(nb,opt) nb = pt.do_center_options(nb,opt) nb = pt.do_cmd_options(nb,opt) nb = pt.do_info_options(nb,opt) nb = pt.do_display_options(nb,opt) ################ # some info ################ print "---------------------------------------------------------" nb.localsystem_of_units.info() nb.ComovingIntegrationInfo() print "---------------------------------------------------------" if opt.obs=='SFR': opt.x = "CosmicTime" opt.ylabel = r'$\rm{Sfr}\,\left[ M_\odot/yr \right]$' opt.xlabel = r'$\rm{Time}\,\left[ Gyr \right]$' if opt.obs=='MSTARS': opt.x = "CosmicTime" opt.ylabel = r'$\rm{Stellar\,\,Mass}\,\left[ M_\odot \right]$' opt.xlabel = r'$\rm{Time}\,\left[ Gyr \right]$' if opt.obs=='MSTARSz': opt.x = "Redshift" opt.ylabel = r'$\rm{Stellar\,\,Mass}\,\left[ M_\odot \right]$' opt.xlabel = r'$\rm{Redshift}$' if opt.obs=='MSTARSa': opt.x = "ScalingFactor" opt.ylabel = r'$\rm{Stellar\,\,Mass}\,\left[ M_\odot \right]$' opt.xlabel = r'$\rm{a}$' ################ # get values ################ x,xlabel = get_value_and_label(nb,opt.x) # if need the current time if opt.x == "CosmicTime": out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Gyr,units.Unit_K]) CurrentTime = nb.CosmicTime(age=nb.atime,units=out_units.UnitTime) # do the 1d histogram G = libgrid.Generic_1d_Grid(opt.xmin,opt.xmax,opt.nx) y = G.get_MassMap(x,nb.mass) # mass weighted histogram x = G.get_r() if opt.obs=='SFR': y_out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_yr,units.Unit_K]) y = y * nb.localsystem_of_units.convertionFactorTo(units.Unit_Ms) dt = (x[1:]-x[:-1])*1e9 # Gyrs to yrs dMdt = y[:-1] / dt y = dMdt x = x[1:] # comoving conversion if nb.isComovingIntegrationOn(): print " converting to physical units (a=%5.3f h=%5.3f)"%(nb.atime,nb.hubbleparam) - y = y/nb.hubbleparam # mass conversion + y = y/nb.hubbleparam # mass conversion if opt.obs=='MSTARS': y_out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_yr,units.Unit_K]) y = y * nb.localsystem_of_units.convertionFactorTo(units.Unit_Ms) y = add.accumulate(y) # comoving conversion if nb.isComovingIntegrationOn(): print " converting to physical units (a=%5.3f h=%5.3f)"%(nb.atime,nb.hubbleparam) - y = y/nb.hubbleparam # mass conversion + y = y/nb.hubbleparam # mass conversion if opt.obs=='MSTARSz': y_out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_yr,units.Unit_K]) y = y * nb.localsystem_of_units.convertionFactorTo(units.Unit_Ms) y = y[::-1] y = add.accumulate(y) y = y[::-1] # comoving conversion if nb.isComovingIntegrationOn(): print " converting to physical units (a=%5.3f h=%5.3f)"%(nb.atime,nb.hubbleparam) - y = y/nb.hubbleparam # mass conversion + y = y/nb.hubbleparam # mass conversion if opt.obs=='MSTARSa': y_out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_yr,units.Unit_K]) y = y * nb.localsystem_of_units.convertionFactorTo(units.Unit_Ms) y = add.accumulate(y) # comoving conversion if nb.isComovingIntegrationOn(): print " converting to physical units (a=%5.3f h=%5.3f)"%(nb.atime,nb.hubbleparam) - y = y/nb.hubbleparam # mass conversion + y = y/nb.hubbleparam # mass conversion data = pt.DataPoints(x,y,color=colors.get(),label=file) datas.append(data) #if file == files[0]: # xmin,xmax,ymin,ymax = pt.SetLimits(opt.xmin,opt.xmax,opt.ymin,opt.ymax,x,y,opt.log) for d in datas: pt.plot(d.x,d.y,color=d.color) # set limits xmin,xmax,ymin,ymax = pt.SetLimitsFromDataPoints(opt.xmin,opt.xmax,opt.ymin,opt.ymax,datas,opt.log) # plot axis pt.SetAxis(xmin,xmax,ymin,ymax,opt.log) pt.xlabel(opt.xlabel,fontsize=pt.labelfont) pt.ylabel(opt.ylabel,fontsize=pt.labelfont) pt.grid(False) if opt.legend: pt.LegendFromDataPoints(datas) # add current Age if opt.x == "CosmicTime": x = [CurrentTime,CurrentTime] y = [0,0] pt.scatter(x,y,c='r',s=25) ######################################################################## # MAIN ######################################################################## if __name__ == '__main__': files,opt = parse_options() pt.InitPlot(files,opt) pt.pcolors MakePlot(files,opt) pt.EndPlot(files,opt)