diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000000000000000000000000000000000000..acae9ebd05d4c57c7135f61e4157a55f3b645ab3
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,6 @@
+*.hdf5
+M*/
+data/__pycache__/
+data/photometry/
+fiducial_files/run.sh
+fiducial_files/swift
diff --git a/data/create_ics_from_snapshot.py b/data/create_ics_from_snapshot.py
new file mode 100644
index 0000000000000000000000000000000000000000..1144d71283a2e57271a4f3b46a754b962e553dc5
--- /dev/null
+++ b/data/create_ics_from_snapshot.py
@@ -0,0 +1,43 @@
+
+import h5py
+import sys
+
+def ics_from_snapshot(inputfile, outputfile):
+    f_input  = h5py.File(inputfile, 'r')
+    f_output = h5py.File(outputfile, 'w')
+
+    f_input.copy('Header', f_output)
+    f_input.copy('Units', f_output)
+
+    g_PT0 = f_output.create_group('PartType0')
+    g_PT4 = f_output.create_group('PartType4')
+
+    f_input.copy('PartType0/Coordinates', g_PT0)
+    f_input.copy('PartType0/InternalEnergies', g_PT0)  # InternalEnergy
+    f_input.copy('PartType0/Masses', g_PT0)
+    f_input.copy('PartType0/ParticleIDs', g_PT0)
+    f_input.copy('PartType0/SmoothingLengths', g_PT0)   # SmoothingLength
+    f_input.copy('PartType0/Velocities', g_PT0)
+
+    f_input.copy('PartType4/Coordinates', g_PT4)
+    f_input.copy('PartType4/Masses', g_PT4)
+    f_input.copy('PartType4/ParticleIDs', g_PT4)
+    f_input.copy('PartType4/BirthTimes', g_PT4)  # StellarFormationTime 
+    f_input.copy('PartType4/Velocities', g_PT4)
+
+    f_output['PartType4/StellarFormationTime'] = f_output['PartType4/BirthTimes']
+    del f_output['PartType4/BirthTimes']
+
+    f_output['PartType0/SmoothingLength'] = f_output['PartType0/SmoothingLengths']
+    del f_output['PartType0/SmoothingLengths']
+
+    f_output['PartType0/InternalEnergy'] = f_output['PartType0/InternalEnergies']
+    del f_output['PartType0/InternalEnergies']
+
+    grp = f_output.create_group('Info')
+    grp.attrs["Input file"] = inputfile
+     
+    f_input.close()
+    f_output.close()
+
+    return
diff --git a/data/getEaglePhotometryTable.sh b/data/getEaglePhotometryTable.sh
new file mode 100755
index 0000000000000000000000000000000000000000..ee9c3b422f19518612416da0913b162fd4a120ff
--- /dev/null
+++ b/data/getEaglePhotometryTable.sh
@@ -0,0 +1,3 @@
+#!/bin/bash
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/YieldTables/EAGLE/photometry.tar.gz
+tar -xf photometry.tar.gz
diff --git a/data/getICs.sh b/data/getICs.sh
new file mode 100755
index 0000000000000000000000000000000000000000..389953f0fa576044507311fd8d23f65061365dff
--- /dev/null
+++ b/data/getICs.sh
@@ -0,0 +1,8 @@
+#!/bin/bash 
+
+# from Swift IsolatedGalaxy_feedback example
+
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/ICs/IsolatedGalaxies_COLIBRE/M4_disk.hdf5
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/ICs/IsolatedGalaxies_COLIBRE/M5_disk.hdf5
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/ICs/IsolatedGalaxies_COLIBRE/M6_disk.hdf5
+
diff --git a/data/getPS20CoolingTables.sh b/data/getPS20CoolingTables.sh
new file mode 100755
index 0000000000000000000000000000000000000000..0945dc565a4c8b1723b21e685b03369e782cd3be
--- /dev/null
+++ b/data/getPS20CoolingTables.sh
@@ -0,0 +1,2 @@
+#!/bin/bash
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/CoolingTables/COLIBRE/UV_dust1_CR1_G1_shield1.hdf5
diff --git a/data/getYieldTable.sh b/data/getYieldTable.sh
new file mode 100755
index 0000000000000000000000000000000000000000..40e24ef34ec47db2b5adde08b8c1a39cd44d6452
--- /dev/null
+++ b/data/getYieldTable.sh
@@ -0,0 +1,4 @@
+#!/bin/bash
+wget http://virgodb.cosma.dur.ac.uk/swift-webstorage/YieldTables/EAGLE/yieldtables.tar.gz
+tar -xf yieldtables.tar.gz
+
diff --git a/fiducial_files/isolated_galaxy.yml b/fiducial_files/isolated_galaxy.yml
new file mode 100644
index 0000000000000000000000000000000000000000..91710a6d55bf07a027857640f80895dfefd2cc19
--- /dev/null
+++ b/fiducial_files/isolated_galaxy.yml
@@ -0,0 +1,171 @@
+# Define some meta-data about the simulation
+MetaData:
+  run_name:   RUNNAME
+
+# Define the system of units to use internally.
+InternalUnitSystem:
+  UnitMass_in_cgs:     1.98848e43    # 10^10 M_sun in grams
+  UnitLength_in_cgs:   3.08567758e21 # 1 kpc in cm
+  UnitVelocity_in_cgs: 1e5           # 1 km/s in cm/s
+  UnitCurrent_in_cgs:  1             # Amperes
+  UnitTemp_in_cgs:     1             # Kelvin
+
+# Parameters for the self-gravity scheme
+Gravity:
+  eta:          0.025                 # Constant dimensionless multiplier for time integration.
+  MAC:          geometric
+  theta_cr:     0.7                   # Opening angle (Multipole acceptance criterion).
+  use_tree_below_softening:  1
+  max_physical_baryon_softening: SOFTENING # Physical softening length (in internal units).
+
+# Parameters governing the time integration (Set dt_min and dt_max to the same value for a fixed time-step run.)
+TimeIntegration:
+  time_begin:        0.    # The starting time of the simulation (in internal units).
+  time_end:          1.2276 # The end time of the simulation (in internal units).
+  dt_min:            1e-9  # The minimal time-step size of the simulation (in internal units).
+  dt_max:            1e-2  # The maximal time-step size of the simulation (in internal units).
+
+# Parameters governing the snapshots
+Snapshots:
+  basename:    output      # Common part of the name of output files
+  time_first:  0.          # Time of the first output if non-cosmological time-integration (in internal units)
+  delta_time:  0.01023     # Time difference between consecutive outputs (in internal units)
+  compression: 4           # Compress the snapshots
+
+# Parameters governing the conserved quantities statistics
+Statistics:
+  delta_time:           1e-2     # Time between statistics output
+  time_first:              0     # (Optional) Time of the first stats output if non-cosmological time-integration (in internal units)
+
+# Parameters related to the initial conditions
+InitialConditions:
+  file_name:               ../../data/INITIALCONDITIONS.hdf5 # The file to read
+  periodic:                0            # Are we running with periodic ICs?
+  stars_smoothing_length:  0.5
+  
+# Parameters for the hydrodynamics scheme
+SPH:
+  resolution_eta:        1.2348   # Target smoothing length in units of the mean inter-particle separation (1.2348 == 48Ngbs with the cubic spline kernel).
+  CFL_condition:         0.1      # Courant-Friedrich-Levy condition for time integration.
+  h_min_ratio:           HMINRATIO # Minimal smoothing in units of softening.
+  h_max:                 10.
+  minimal_temperature:   10.
+
+# Parameters for the stars neighbour search
+Stars:
+  resolution_eta:        1.1642  # Target smoothing length in units of the mean inter-particle separation
+  h_tolerance:           7e-3
+  overwrite_birth_time:    0     # Make sure the stars in the ICs do not do any feedback
+  birth_time:             -1.    # by setting all of their birth times to -1  
+  timestep_age_threshold_Myr: 10.   # Age at which stars switch from young to old (in Mega-years).
+  max_timestep_young_Myr:     0.1   # Maximal time-step length of young stars (in Mega-years).
+  max_timestep_old_Myr:       1.0   # Maximal time-step length of old stars (in Mega-years).
+  luminosity_filename:   ../../data/photometry
+
+# COLIBRE cooling parameters
+COLIBRECooling:
+  dir_name:                ../../data/UV_dust1_CR1_G1_shield1.hdf5 # Location of the Ploeckinger+20 cooling tables
+  H_reion_z:               7.5               # Redshift of Hydrogen re-ionization (Planck 2018)
+  H_reion_eV_p_H:          2.0
+  He_reion_z_centre:       3.5               # Redshift of the centre of the Helium re-ionization Gaussian
+  He_reion_z_sigma:        0.5               # Spread in redshift of the  Helium re-ionization Gaussian
+  He_reion_eV_p_H:         2.0               # Energy inject by Helium re-ionization in electron-volt per Hydrogen atom
+  delta_logTEOS_subgrid_properties: 0.3      # delta log T above the EOS below which the subgrid properties use Teq assumption
+  rapid_cooling_threshold:          0.333333 # Switch to rapid cooling regime for dt / t_cool above this threshold.
+
+EAGLEChemistry:
+  init_abundance_metal:     0.0133714      # Inital fraction of particle mass in *all* metals 
+  init_abundance_Hydrogen:  0.73738788833  # Inital fraction of particle mass in Hydrogen
+  init_abundance_Helium:    0.24924186942  # Inital fraction of particle mass in Helium
+  init_abundance_Carbon:    0.0023647215   # Inital fraction of particle mass in Carbon
+  init_abundance_Nitrogen:  0.0006928991   # Inital fraction of particle mass in Nitrogen
+  init_abundance_Oxygen:    0.00573271036  # Inital fraction of particle mass in Oxygen
+  init_abundance_Neon:      0.00125649278  # Inital fraction of particle mass in Neon
+  init_abundance_Magnesium: 0.00070797838  # Inital fraction of particle mass in Magnesium
+  init_abundance_Silicon:   0.00066495154  # Inital fraction of particle mass in Silicon
+  init_abundance_Iron:      0.00129199252  # Inital fraction of particle mass in Iron
+
+# Hernquist potential parameters
+HernquistPotential:
+  useabspos:       0       # 0 -> positions based on centre, 1 -> absolute positions 
+  position:        [0.,0.,0.]    # Location of centre of isothermal potential with respect to centre of the box (if 0) otherwise absolute (if 1) (internal units)
+  idealizeddisk:   1       # Run with an idealized galaxy disk
+  M200:            137.0   # M200 of the galaxy disk
+  h:               0.704   # reduced Hubble constant (value does not specify the used units!)
+  concentration:   9.0     # concentration of the Halo
+  diskfraction:    0.040   # Disk mass fraction
+  bulgefraction:   0.0     # Bulge mass fraction
+  timestep_mult:   0.01    # Dimensionless pre-factor for the time-step condition, basically determines the fraction of the orbital time we use to do the time integration
+  epsilon:         0.4      # Softening size (internal units)
+ 
+# EAGLE star formation parameters
+EAGLEStarFormation:
+  SF_threshold:                      Subgrid      # Zdep (Schaye 2004) or Subgrid
+  SF_model:                          SchmidtLaw  # PressureLaw (Schaye et al. 2008) or SchmidtLaw
+  star_formation_efficiency:         SFEFFICIENCY
+  min_over_density:                  100.0        # The over-density above which star-formation is allowed.
+  threshold_temperature1_K:          1000.        # When using subgrid-based SF threshold, subgrid temperature below which gas is star-forming.
+  threshold_temperature2_K:          1000.        # When using subgrid-based SF threshold, subgrid temperature below which gas is star-forming if also above the density limit.
+  threshold_number_density_H_p_cm3:  0.1          # When using subgrid-based SF threshold, subgrid number density above which gas is star-forming if also below the second temperature limit.
+  
+# Parameters for the EAGLE "equation of state"
+EAGLEEntropyFloor:
+  Jeans_density_threshold_H_p_cm3: EOSDENSITY # Physical density above which the EAGLE Jeans limiter entropy floor kicks in expressed in Hydrogen atoms per cm^3.
+  Jeans_over_density_threshold:    10.       # Overdensity above which the EAGLE Jeans limiter entropy floor can kick in.
+  Jeans_temperature_norm_K:        EOSTEMPERATURE # Temperature of the EAGLE Jeans limiter entropy floor at the density threshold expressed in Kelvin.
+  Jeans_gamma_effective:           EOSSLOPE # Slope the of the EAGLE Jeans limiter entropy floor
+  Cool_density_threshold_H_p_cm3: 1e-5       # Physical density above which the EAGLE Cool limiter entropy floor kicks in expressed in Hydrogen atoms per cm^3.
+  Cool_over_density_threshold:    10.        # Overdensity above which the EAGLE Cool limiter entropy floor can kick in.
+  Cool_temperature_norm_K:        10.        # Temperature of the EAGLE Cool limiter entropy floor at the density threshold expressed in Kelvin. (NOTE: This is below the min T and hence this floor does nothing)
+  Cool_gamma_effective:           1.         # Slope the of the EAGLE Cool limiter entropy floor
+
+# EAGLE feedback model with constant feedback energy fraction
+EAGLEFeedback:
+  use_SNII_feedback:                    1               # Global switch for SNII thermal (stochastic) feedback.
+  use_SNIa_feedback:                    1               # Global switch for SNIa thermal (continuous) feedback.
+  use_AGB_enrichment:                   1               # Global switch for enrichement from AGB stars.
+  use_SNII_enrichment:                  1               # Global switch for enrichement from SNII stars.
+  use_SNIa_enrichment:                  1               # Global switch for enrichement from SNIa stars.
+  filename:                             ../../data/yieldtables/  # Path to the directory containing the EAGLE yield tables.
+  IMF_min_mass_Msun:                    0.1             # Minimal stellar mass considered for the Chabrier IMF in solar masses.
+  IMF_max_mass_Msun:                  100.0             # Maximal stellar mass considered for the Chabrier IMF in solar masses.
+  SNII_min_mass_Msun:                   8.0             # Minimal mass considered for SNII stars in solar masses.
+  SNII_max_mass_Msun:                 100.0             # Maximal mass considered for SNII stars in solar masses.
+  SNII_feedback_model:                  MinimumDistance # Feedback modes: Random, Isotropic, MinimumDistance, MinimumDensity
+  SNII_sampled_delay:                   1               # Sample the SNII lifetimes to do feedback.
+  SNII_delta_T_K:                       3.16228e7       # Change in temperature to apply to the gas particle in a SNII thermal feedback event in Kelvin.
+  SNII_energy_erg:                      1.0e51          # Energy of one SNII explosion in ergs.
+  SNII_energy_fraction_function:        EAGLE           # Type of functional form to use for scaling the energy fraction with density and metallicity ('EAGLE', 'Separable', or 'Independent').
+  SNII_energy_fraction_min:             1.0             # Minimal fraction of energy applied in a SNII feedback event.
+  SNII_energy_fraction_max:             1.0             # Maximal fraction of energy applied in a SNII feedback event.
+  SNII_energy_fraction_Z_0:             0.0012663729    # Pivot point for the metallicity dependance of the SNII energy fraction (metal mass fraction).
+  SNII_energy_fraction_n_0_H_p_cm3:     1.4588          # Pivot point for the birth density dependance of the SNII energy fraction in cm^-3.
+  SNII_energy_fraction_n_Z:             0.8686          # Power-law for the metallicity dependance of the SNII energy fraction.
+  SNII_energy_fraction_n_n:             0.8686          # Power-law for the birth density dependance of the SNII energy fraction.
+  SNII_energy_fraction_use_birth_density: 0             # Are we using the density at birth to compute f_E or at feedback time?
+  SNII_energy_fraction_use_birth_metallicity: 0         # Are we using the metallicity at birth to compuote f_E or at feedback time?
+  SNIa_DTD:                             Exponential     # Functional form of the SNIa delay time distribution.
+  SNIa_DTD_delay_Gyr:                   0.04            # Stellar age after which SNIa start in Gyr (40 Myr corresponds to stars ~ 8 Msun).
+  SNIa_DTD_exp_timescale_Gyr:           2.0             # Time-scale of the exponential decay of the SNIa rates in Gyr.
+  SNIa_DTD_exp_norm_p_Msun:             0.002           # Normalisation of the SNIa rates in inverse solar masses.
+  SNIa_energy_erg:                     1.0e51           # Energy of one SNIa explosion in ergs.
+  AGB_ejecta_velocity_km_p_s:          10.0             # Velocity of the AGB ejectas in km/s.
+  stellar_evolution_age_cut_Gyr:        0.1             # Age in Gyr above which the enrichment is downsampled.
+  stellar_evolution_sampling_rate:       10             # Number of time-steps in-between two enrichment events for a star above the age threshold.
+  SNII_yield_factor_Hydrogen:           1.0             # (Optional) Correction factor to apply to the Hydrogen yield from the SNII channel.
+  SNII_yield_factor_Helium:             1.0             # (Optional) Correction factor to apply to the Helium yield from the SNII channel.
+  SNII_yield_factor_Carbon:             0.5             # (Optional) Correction factor to apply to the Carbon yield from the SNII channel.
+  SNII_yield_factor_Nitrogen:           1.0             # (Optional) Correction factor to apply to the Nitrogen yield from the SNII channel.
+  SNII_yield_factor_Oxygen:             1.0             # (Optional) Correction factor to apply to the Oxygen yield from the SNII channel.
+  SNII_yield_factor_Neon:               1.0             # (Optional) Correction factor to apply to the Neon yield from the SNII channel.
+  SNII_yield_factor_Magnesium:          2.0             # (Optional) Correction factor to apply to the Magnesium yield from the SNII channel.
+  SNII_yield_factor_Silicon:            1.0             # (Optional) Correction factor to apply to the Silicon yield from the SNII channel.
+  SNII_yield_factor_Iron:               0.5             # (Optional) Correction factor to apply to the Iron yield from the SNII channel.
+
+Restarts:
+  delta_hours:        2.0
+  max_run_time:      71.5
+  onexit:             1
+  resubmit_on_exit:   1
+  resubmit_command:   ./resubmit.sh
+
diff --git a/fiducial_files/resubmit.job b/fiducial_files/resubmit.job
new file mode 100644
index 0000000000000000000000000000000000000000..28eacbadc7a1c6d02d38644faf886ce927cb4dd2
--- /dev/null
+++ b/fiducial_files/resubmit.job
@@ -0,0 +1,19 @@
+#!/bin/bash
+#SBATCH -J RUNNAME
+#SBATCH -o output-%j.out
+#SBATCH -e error-%j.err
+#SBATCH -t 4320
+#SBATCH -N 1
+
+module purge
+module load hdf5/1.12.2-intel-19.1.3.304-eqsfbup
+module load libtool/2.4.7-intel-2021.5.0-uhp7fp3
+module load fftw/3.3.10-intel-19.1.3.304-kcawvtk
+module load metis/5.1.0-intel-2021.5.0-wijrhuu
+module load gsl/2.7.1-gcc-12.2.0-jgxkcb6
+module load python/3.10.6-intel-2021.5.0-mtvbpzb
+module load intel/19.1.3.304
+
+./swift --threads=24 --external-gravity --self-gravity --stars --star-formation --cooling --temperature --hydro --feedback --limiter --sync --restart isolated_galaxy.yml
+
+
diff --git a/fiducial_files/resubmit.sh b/fiducial_files/resubmit.sh
new file mode 100755
index 0000000000000000000000000000000000000000..122426a6e0dbb1179b236f19f60fcd58c35b0217
--- /dev/null
+++ b/fiducial_files/resubmit.sh
@@ -0,0 +1,3 @@
+#!/bin/sh
+
+sbatch resubmit.job
diff --git a/fiducial_files/submit.job b/fiducial_files/submit.job
new file mode 100644
index 0000000000000000000000000000000000000000..2bf295e739d42742c8f140efe400afcb7639ca2f
--- /dev/null
+++ b/fiducial_files/submit.job
@@ -0,0 +1,21 @@
+#!/bin/bash
+#SBATCH -J RUNNAME
+#SBATCH -o output-%j.out
+#SBATCH -e error-%j.err
+#SBATCH --time=TIMEINHOURS:00:00
+#SBATCH --partition=skylake_0096
+#SBATCH --qos=skylake_0096
+#SBATCH -N 1
+
+module purge
+module load hdf5/1.12.2-intel-19.1.3.304-eqsfbup
+module load libtool/2.4.7-intel-2021.5.0-uhp7fp3
+module load fftw/3.3.10-intel-19.1.3.304-kcawvtk
+module load metis/5.1.0-intel-2021.5.0-wijrhuu
+module load gsl/2.7.1-gcc-12.2.0-jgxkcb6
+module load python/3.10.6-intel-2021.5.0-mtvbpzb
+module load intel/19.1.3.304
+
+./swift --threads=24 --external-gravity --self-gravity --stars --star-formation --cooling --temperature --hydro --feedback --limiter --sync isolated_galaxy.yml
+
+
diff --git a/setup_reruns.py b/setup_reruns.py
new file mode 100644
index 0000000000000000000000000000000000000000..0ac53ba62d1765ef95335a196c7a11c08f337abc
--- /dev/null
+++ b/setup_reruns.py
@@ -0,0 +1,78 @@
+import os
+from pathlib import Path
+import shutil
+import glob
+import re
+
+from data.create_ics_from_snapshot import ics_from_snapshot
+
+##########################################
+# Select mass resolution and snapshotnumber
+##########################################
+mass_resolution_levels = ["M6"]
+snapshotnumber = 85
+##########################################
+
+
+
+localdir = os.getcwd()
+for mass in mass_resolution_levels:
+    os.chdir(os.path.join(localdir,mass))
+    for runfolder in glob.glob("Galaxy*3/"):
+        ##########################################
+        # setting runfolder name and creating rerunfolder
+        ##########################################
+        runfolder_absolute_path   = os.path.join(localdir, mass, runfolder)
+        rerunfolder_absolute_path = os.path.join(localdir, mass+"_reruns_snap%4.4i"%(snapshotnumber), "Rerun"+runfolder)
+        Path(rerunfolder_absolute_path).mkdir(parents=True, exist_ok=True)
+        print (runfolder_absolute_path)
+
+        ##########################################
+        # change working directory to rerunfolder
+        ##########################################
+        os.chdir(rerunfolder_absolute_path)
+
+        ##########################################
+        # copy files into rerunfolder
+        ##########################################
+        shutil.copy2(os.path.join(runfolder_absolute_path, "isolated_galaxy.yml"), "./")
+        shutil.copy2(os.path.join(runfolder_absolute_path, "swift"), "./")
+        shutil.copy2(os.path.join(runfolder_absolute_path, "output_%4.4i.hdf5"%(int(snapshotnumber))), "./")
+
+        ##########################################
+        # create ICs for rerun
+        ##########################################
+        ics_from_snapshot(os.path.join(rerunfolder_absolute_path, "output_%4.4i.hdf5"%(int(snapshotnumber))), \
+                          os.path.join(rerunfolder_absolute_path, "ICs.hdf5"))
+
+        ##########################################
+        # update yml file for rerun
+        ##########################################
+        f = open("isolated_galaxy.yml", "rt")
+        data = f.read()
+        data = re.sub(r'^.*disk\.hdf5.*$', '  file_name:               ICs.hdf5', data, flags=re.MULTILINE)
+        data = re.sub(r'^.*h_min_ratio.*$', '  h_min_ratio:           0.0001', data, flags=re.MULTILINE)
+        data = re.sub(r'^.*time_end.*$', '  time_end:         0.00001025', data, flags=re.MULTILINE)
+        data = re.sub(r'^.*delta_time\:  0.*$', '  delta_time:  0.00001023', data, flags=re.MULTILINE)
+        data = re.sub(r'^.*dt_max\:.*$', '  dt_max:            1e-8', data, flags=re.MULTILINE)
+
+        f.close()
+        f = open("isolated_galaxy.yml", "wt")
+        f.write(data)
+        f.close()        
+
+        ##########################################
+        # start swift for one timestep
+        ##########################################
+
+        os.system("./swift --threads=28 --external-gravity --self-gravity --stars --star-formation --cooling --temperature --hydro --feedback --limiter --sync -n 1 isolated_galaxy.yml")
+
+
+
+
+
+
+
+
+
+