Merge pull request #125 from andrew-saydjari/2023_03_01

2023 03 01

pipeline.jl

@@ -50,14 +50,17 @@ using LibGit2; git_branch, git_commit = initalize_git(src_dir); @passobj 1 worke
     runlist_range = 295:295 #1:600 #295, 245, 335, 101
     batchsize = 10 #40
 
-    cache_dir = "../local_cache_295_redo/"
-    inject_cache_dir = prior_dir*"2024_02_08/inject_local_cache"
+    RV_err_step = 4 # in subpix steps (used to be 8)
+
+    cache_dir = "../local_cache_inject295/"
+    inject_cache_dir = prior_dir*"2024_03_01/inject_local_cache_no_dibs"
 
     # Prior Dictionary
     prior_dict = Dict{String,String}()
 
     # Input List (not really a prior, but an input file we search for stars conditioned on)
-    prior_dict["runlists"] = prior_dir*"2024_01_19/outlists/dr17_dr17_star_input_lst_msked_"
+    prior_dict["runlists"] = prior_dir*"2024_03_01/inject_no_dibs_295/injection_input_lst_"
+    # prior_dict["runlists"] = prior_dir*"2024_01_19/outlists/dr17_dr17_star_input_lst_msked_"
 
     # Sky Priors
     prior_dict["skycont"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/sky_priors/APOGEE_skycont_svd_30_f"
@@ -72,10 +75,10 @@ using LibGit2; git_branch, git_commit = initalize_git(src_dir); @passobj 1 worke
     prior_dict["starLines_LSF"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/starLine_priors_norm94/APOGEE_stellar_kry_50_subpix_f"
 
     # DIB Priors
-    prior_dict["DIB_noLSF"] = prior_dir*"2023_07_22/dib_priors/precomp_dust_1_analyticDeriv_stiff.h5"
-    prior_dict["DIB_noLSF_soft"] = prior_dir*"2023_07_22/dib_priors/precomp_dust_3_analyticDeriv_soft.h5"
-    prior_dict["DIB_LSF"] = prior_dir*"2023_07_22/dib_priors/precomp_dust_1_analyticDerivLSF_stiff_"
-    prior_dict["DIB_LSF_soft"] = prior_dir*"2023_07_22/dib_priors/precomp_dust_3_analyticDerivLSF_soft_"
+    prior_dict["DIB_noLSF"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/dib_priors/precomp_dust_1_analyticDeriv_stiff.h5"
+    prior_dict["DIB_noLSF_soft"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/dib_priors/precomp_dust_3_analyticDeriv_soft.h5"
+    prior_dict["DIB_LSF"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/dib_priors/precomp_dust_1_analyticDerivLSF_stiff_"
+    prior_dict["DIB_LSF_soft"] = prior_dir*"2024_02_21/apMADGICS.jl/src/prior_build/dib_priors/precomp_dust_3_analyticDerivLSF_soft_"
 
     # Data for Detector Cals (not really a prior, but an input the results depend on in detail)
     prior_dict["chip_fluxdep_err_correction"] = src_dir*"data/chip_fluxdep_err_correction.jdat"
@@ -188,6 +191,8 @@ end
             fvec, fvarvec, cntvec, chipmidtimes, metaexport = deserialize(starcache)
             starscale,framecnts,a_relFlux,b_relFlux,c_relFlux,cartVisit = metaexport
         elseif tele[end]=='i'
+            println(starcache)
+            println(caching)
             @warn "Injections not found at injection cache dir!"
         else
             fvec, fvarvec, cntvec, chipmidtimes, metaexport = stack_out(release_dir,redux_ver,tele,field,plate,mjd,fiberindx,cache_dir=cache_dir)
@@ -273,7 +278,7 @@ end
         else
             Base.Fix2(chi2_wrapper,(rvmsk,Ctotinv_cur,Xd_obs,starCont_Mscale,V_subpix,pre_Vslice))
         end
-        lout = sampler_1d_hierarchy_var(chi2_wrapper_partial,slvl_tuple,minres=1//10,stepx=8)
+        lout = sampler_1d_hierarchy_var(chi2_wrapper_partial,slvl_tuple,minres=1//10,stepx=RV_err_step)
         svalc = lout[1][3]
         push!(out,lout) # 2
 
@@ -335,7 +340,7 @@ end
         skyscale1 = nanzeromedian(x_comp_out[4])
         dvec = (fvec .-(x_comp_out[2].+x_comp_out[3].+x_comp_out[4].+x_comp_out[5].*(1 .+ nanify(x_comp_lst[5],finalmsk))))./fvec;
         chi2res = x_comp_lst[1]'*(Ainv*x_comp_lst[1])
-        chi2r_fc = chi2red_fluxscale(chi2res./count(finalmsk), starscale1, fc=red_chi2_dict[tele])
+        chi2r_fc = chi2red_fluxscale(chi2res./count(finalmsk), starscale1, fc=red_chi2_dict[tele[1:6]])
         push!(out,(chi2res,chi2r_fc,nanzeroiqr(dvec),count(finalmsk),starscale1,skyscale1)) # 3
         push!(out,x_comp_out) # 4
         dflux_starlines = sqrt_nan.(get_diag_posterior_from_prior_asym(Ctotinv_fut, V_starlines_c, V_starlines_r))
@@ -393,7 +398,7 @@ end
                         x_comp_lst[5:end]...]
 
             chi2res = x_comp_lst[1]'*(Ainv*x_comp_lst[1])
-            chi2r_fc = chi2red_fluxscale(chi2res./count(finalmsk), starscale1, fc=red_chi2_dict[tele])
+            chi2r_fc = chi2red_fluxscale(chi2res./count(finalmsk), starscale1, fc=red_chi2_dict[tele[1:6]])
             push!(out,(chi2res,chi2r_fc)) # 8
 
             push!(out,x_comp_out) # 9

src/ingest.jl

@@ -158,7 +158,7 @@ function stack_out(release_dir,redux_ver,tele,field,plate,mjd,fiberindx; telluri
             Xd = read(f[2],:,fiberindx)
             Xd_stack[(1:2048).+(chipind-1)*2048] .= Xd[end:-1:1]; # ./thrptDict[chip] has already been done at the 2D -> 1D level
             Xd_std = read(f[3],:,fiberindx)
-            Xd_std_stack[(1:2048).+(chipind-1)*2048] .= Xd_std[end:-1:1].*err_factor.(Xd[end:-1:1],Ref(err_correct_Dict[join([tele,chip],"_")])); # ./thrptDict[chip] has already been done at the 2D -> 1D level
+            Xd_std_stack[(1:2048).+(chipind-1)*2048] .= Xd_std[end:-1:1].*err_factor.(Xd[end:-1:1],Ref(err_correct_Dict[join([tele[1:6],chip],"_")])); # ./thrptDict[chip] has already been done at the 2D -> 1D level
             pixmsk = read(f[4],:,fiberindx);
             pixmsk_stack[(1:2048).+(chipind-1)*2048] .= pixmsk[end:-1:1]
             waveobsa = read(f[5],:,fiberindx);

src/inject_test.jl

@@ -52,21 +52,20 @@ using LibGit2; git_branch, git_commit = initalize_git(src_dir); @passobj 1 worke
 
     RV_range_pix = (-68,68) # pixscale is ~4.14 km/s per pixel
 
-    dib_inject = true
-    dib_center_lambda_lst = [15273,15672]
-    dib_ew_range = (-1.5,0)
-    dib_sig_range = (0.7,3.7)
-    dib_vel_range = (-450, 450) # km/s
-
     skycont_only = false
     no_sky = false
     dibs_on=true
 
+    dib_center_lambda_lst = [15273] #,15672]
+    dib_ew_range = (-1.5,0)
+    dib_sig_range = (0.7,3.7)
+    dib_vel_range = (-450, 450) # km/s
+
     # Prior Dictionary
     prior_dict = Dict{String,String}()
 
-    prior_dict["out_dir"] = prior_dir*"2024_03_01/inject_dibs_295/"
-    prior_dict["inject_cache_dir"] = prior_dir*"2024_03_01/inject_cache_dir/"
+    prior_dict["out_dir"] = prior_dir*"2024_03_01/inject_no_dibs_295/"
+    prior_dict["inject_cache_dir"] = prior_dir*"2024_03_01/inject_local_cache_no_dibs/"
     prior_dict["local_cache"] = prior_dir*"2024_03_01/local_cache/"
 
     prior_dict["past_run"] = prior_dir*"2024_02_21/outdir_wu_295_th/apMADGICS_out.h5" # used for StarScale distribution only
@@ -198,7 +197,7 @@ end
             starcomp .*= (1 .+ starLines)
             if dibs_on
                 for (dib_ind, dib_center_lambda) in enumerate(dib_center_lambda_lst)
-                    dibcomp = gauss1d(ew[dib_ind],λ0[dib_ind],sigma[dib_ind],wavetarg)
+                    dibcomp = gauss1d_ew(ew[dib_ind],λ0[dib_ind],sigma[dib_ind],wavetarg)
                     starcomp .*= (1 .+ dibcomp)
                 end
             end
@@ -295,7 +294,7 @@ dib_ew = []
 for dib_center_lambda in dib_center_lambda_lst
     push!(dib_sig,rand(rng,Uniform(dib_sig_range...),nsamp));
     push!(dib_lam,rand(rng,Uniform(dib_center_lambda*(1+dib_vel_range[1]/c),dib_center_lambda*(1+dib_vel_range[2]/c)),nsamp));
-    if dib_inject
+    if dibs_on
         push!(dib_ew,rand(rng,Uniform(dib_ew_range...),nsamp))
     else
         push!(dib_ew,zeros(nsamp))

submit_np.sh

@@ -6,7 +6,7 @@
 
 #SBATCH --mem=0 #requesting all of the memory on the node
 
-#SBATCH --time=96:00:00
+#SBATCH --time=12:00:00
 #SBATCH --job-name=apMADGICS
 #SBATCH --output=%x_%j.out
 #SBATCH --err=%x_%j.err

test/gridSearch.jl

@@ -35,5 +35,22 @@
 
     out = sampler_2d_hierarchy_var(tstfun2d,tsttup_2d)
     @test out[1] == ((1.5, 2.5), 3.0, (3//2, 5//2), 3.0, CartesianIndex(3, 3), 0, 0.5, 0.5, -0.0, 0.5, 0.5)
-
+
+    alpha = 1
+    x = -10:1//10:10
+    y = alpha.*x.^2;
+    minind = findmin(y)[2]
+    @test err1d(x,y,minind;stepx=60) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=22) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=5) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=1) == 1//(2*alpha)
+
+    alpha = 4
+    x = -10:1//10:10
+    y = alpha.*x.^2;
+    minind = findmin(y)[2]
+    @test err1d(x,y,minind;stepx=60) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=22) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=5) == 1//(2*alpha)
+    @test err1d(x,y,minind;stepx=1) == 1//(2*alpha)
 end