@@ -43,6 +43,10 @@ function parse_commandline()
4343 arg_type = Int64
4444 #default = [ min_order(NEID2D()), max_order(NEID2D()) ]
4545 #default = [ first(orders_to_use_default(NEID2D())), last(orders_to_use_default(NEID2D())) ]
46+ "--mask_shape"
47+ help = "Specify CCF mask shape: gaussian, tophat, etc."
48+ arg_type = String
49+ default = "gaussian"
4650 "--mask_scale_factor"
4751 help = "Specify CCF mask width scale as multiple of NEID default v width " * string(default_ccf_mask_v_width(NEID2D()))
4852 arg_type = Float64
310314 end
311315 @assert isfile (line_list_filename) || islink (line_list_filename)
312316 if args[" sed_filename" ] != nothing
313- @warn (" DRP v1.1 now provides blaze in L2 file. This script has not been updated to use explicit an SED model." )
317+ @warn (" DRP >= v1.1 now provides blaze in L2 file. This script has not been updated to use explicit an SED model." )
314318 sed_filename = args[" sed_filename" ]
315319 # elseif !@isdefined sed_filename
316320 # sed_filename = joinpath("/home/eford/Code/RvSpectMLEcoSystem/NeidSolarScripts/data","neidMaster_HR_SmoothLampSED_20210101.fits")
342346 end
343347 @assert 50 <= v_step <= 2000 # m/s
344348
349+ mask_shape = symbol (args[" mask_shape" ])
350+ @assert in (mask_shape,[:gaussian ,:supergaussian ,:tophat ,:halfcos ])
345351 if args[" mask_scale_factor" ] != nothing
346352 mask_scale_factor = args[" mask_scale_factor" ]
347353 elseif ! @isdefined mask_scale_factor
@@ -446,17 +452,17 @@ if verbose println("# Reading manifest of files to process.") end
446452 df_files_use = df_files_use |>
447453 @filter ( args[" min_expmeter" ] == nothing || _. expmeter_mean >= args[" min_expmeter" ] ) |>
448454 DataFrame
449- if ! hasproperty (df_files_use,:mean_pyroflux ) || ! hasproperty (df_files_use,:rms_pyroflux ) || any (ismissing .(df_files_use. mean_pyroflux)) || any (ismissing .(df_files_use. rms_pyroflux))
455+ if ! hasproperty (df_files_use,:mean_pyroflux ) || ! hasproperty (df_files_use,:rms_pyroflux ) || all (ismissing .(df_files_use. mean_pyroflux)) || all (ismissing .(df_files_use. rms_pyroflux))
450456 @error " Manifest file doesn't include valid mean_pyroflux and/or rms_pyroflux." manifest_filename
451457 end
452458 df_files_use = df_files_use |>
453- @filter ( args[" min_pyrhelio" ] == nothing || _. mean_pyroflux >= args[" min_pyrhelio" ] ) |>
459+ @filter ( args[" min_pyrhelio" ] == nothing || (( _. mean_pyroflux >= args[" min_pyrhelio" ]) && ( ! ismissing (_ . mean_pyroflux))) ) |>
454460 DataFrame
455461 df_files_use = df_files_use |>
456462 @filter ( args[" max_expmeter_rms_frac" ] == nothing || _. expmeter_rms <= args[" max_expmeter_rms_frac" ]* _. expmeter_mean ) |>
457463 DataFrame
458464 df_files_use = df_files_use |>
459- @filter ( args[" max_pyrhelio_rms_frac" ] == nothing || _ . rms_pyroflux <= args[" max_pyrhelio_rms_frac" ]* _. mean_pyroflux ) |>
465+ @filter ( args[" max_pyrhelio_rms_frac" ] == nothing || ( ! ismissing (_ . mean_pyroflux) && ! ismissing (_ . rms_pyroflux) && (_ . rms_pyroflux <= args[" max_pyrhelio_rms_frac" ]* _. mean_pyroflux)) ) |>
460466 DataFrame
461467 df_files_use = df_files_use |>
462468 @filter ( args[" min_expmeter_to_pyrhelio" ] == nothing || _. expmeter_mean >= args[" min_expmeter_to_pyrhelio" ]* _. mean_pyroflux ) |>
@@ -763,12 +769,40 @@ end
763769if verbose println (now ()) end
764770line_width_50 = line_width_50_default
765771 # order_list_timeseries = extract_orders(all_spectra,pipeline_plan, orders_to_use=orders_to_use, remove_bad_chunks=false, recalc=true )
772+ @time (order_ccfs, order_ccf_vars, v_grid_order_ccfs) = ccf_orders (order_list_timeseries, line_list_espresso, pipeline_plan,
773+ mask_type= mask_shape, Δfwhm= Δfwhm,
774+ mask_scale_factor= mask_scale_factor, range_no_mask_change= 5 * line_width_50, ccf_mid_velocity= ccf_mid_velocity, v_step= v_step,
775+ v_max= max (range_no_mask_change* line_width_50,2 * max_bc), orders_to_use= orders_to_use, allow_nans= true , calc_ccf_var= true ,
776+ recalc= true )
777+ #=
778+ if mask_shape == "gaussian"
766779 @time (order_ccfs, order_ccf_vars, v_grid_order_ccfs) = ccf_orders(order_list_timeseries, line_list_espresso, pipeline_plan,
767780 mask_type=:gaussian, Δfwhm=Δfwhm,
768781 mask_scale_factor=mask_scale_factor, range_no_mask_change=5*line_width_50, ccf_mid_velocity=ccf_mid_velocity, v_step=v_step,
769782 v_max=max(range_no_mask_change*line_width_50,2*max_bc), orders_to_use=orders_to_use, allow_nans=true, calc_ccf_var=true,
770783 recalc=true)
771-
784+ elseif mask_shape == "tophat"
785+ @time (order_ccfs, order_ccf_vars, v_grid_order_ccfs) = ccf_orders(order_list_timeseries, line_list_espresso, pipeline_plan,
786+ mask_type=:tophat, Δfwhm=Δfwhm,
787+ mask_scale_factor=mask_scale_factor, range_no_mask_change=5*line_width_50, ccf_mid_velocity=ccf_mid_velocity, v_step=v_step,
788+ v_max=max(range_no_mask_change*line_width_50,2*max_bc), orders_to_use=orders_to_use, allow_nans=true, calc_ccf_var=true,
789+ recalc=true)
790+ elseif mask_shape == "supergaussian"
791+ @time (order_ccfs, order_ccf_vars, v_grid_order_ccfs) = ccf_orders(order_list_timeseries, line_list_espresso, pipeline_plan,
792+ mask_type=:supergaussian, Δfwhm=Δfwhm,
793+ mask_scale_factor=mask_scale_factor, range_no_mask_change=5*line_width_50, ccf_mid_velocity=ccf_mid_velocity, v_step=v_step,
794+ v_max=max(range_no_mask_change*line_width_50,2*max_bc), orders_to_use=orders_to_use, allow_nans=true, calc_ccf_var=true,
795+ recalc=true)
796+ elseif mask_shape == "halfcos"
797+ @time (order_ccfs, order_ccf_vars, v_grid_order_ccfs) = ccf_orders(order_list_timeseries, line_list_espresso, pipeline_plan,
798+ mask_type=:halfcos, Δfwhm=Δfwhm,
799+ mask_scale_factor=mask_scale_factor, range_no_mask_change=5*line_width_50, ccf_mid_velocity=ccf_mid_velocity, v_step=v_step,
800+ v_max=max(range_no_mask_change*line_width_50,2*max_bc), orders_to_use=orders_to_use, allow_nans=true, calc_ccf_var=true,
801+ recalc=true)
802+ else
803+ @warn "Invalid mask_shape, using gaussian."
804+ end
805+ =#
772806# orders_to_use2 = orders_to_use[map(i->!iszero(order_ccfs[:,i,:]),1:length(orders_to_use))]
773807# order_list_timeseries2 = extract_orders(all_spectra,pipeline_plan, orders_to_use=orders_to_use2, remove_bad_chunks=false, recalc=true )
774808
0 commit comments