jfjoch_process: Add option to use azimuthal integration as background for Bragg integration

tools/jfjoch_process.cpp

@@ -45,39 +45,75 @@ void print_usage() {
 
     std::cout << "  Spot finding" << std::endl;
     std::cout << "   --spot-sigma <num>               Noise sigma level for spot finding (default: 3.0)" << std::endl;
-    std::cout << "   --spot-threshold <num>           Photon count threshold for spot finding (default: 10)" << std::endl;
-    std::cout << "   --spot-high-resolution <num>     High resolution limit for spot finding (default: 1.5)" << std::endl;
+    std::cout << "   --spot-threshold <num>           Photon count threshold for spot finding (default: 10)" <<
+            std::endl;
+    std::cout << "   --spot-high-resolution <num>     High resolution limit for spot finding (default: 1.5)" <<
+            std::endl;
     std::cout << "   --max-spots <num>                Max spot count (default: 250)" << std::endl;
     std::cout << std::endl;
 
+    std::cout << "  Azimuthal integration" << std::endl;
+    std::cout << "    -Q, --azim-q-spacing <num>            Q spacing for azimuthal integration (default: 0.01)" <<
+            std::endl;
+    std::cout << "    --azim-q-min <num>                    Minimum Q value for azimuthal integration (default: 0.0)" <<
+            std::endl;
+    std::cout << "    --azim-q-max <num>                    Maximum Q value for azimuthal integration (default: 5.0)" <<
+            std::endl;
+    std::cout << std::endl;
+
     std::cout << "  Indexing" << std::endl;
-    std::cout << "   -R, --two-pass-rotation[=num]    Two-pass offline rotation indexing (optional: number of images, default: 30)" << std::endl;
-    std::cout << "   --single-pass-rotation[=num]     Use online-like single-pass rotation indexing (optional: min angular range deg)" << std::endl;
+    std::cout <<
+            "   -R, --two-pass-rotation[=num]    Two-pass offline rotation indexing (optional: number of images, default: 30)"
+            << std::endl;
+    std::cout <<
+            "   --single-pass-rotation[=num]     Use online-like single-pass rotation indexing (optional: min angular range deg)"
+            << std::endl;
     std::cout << "   --redo-rotation-spots            Redo spot finding for two-pass rotation indexing" << std::endl;
-    std::cout << "   --force-rotation-lattice <vec>   Force rotation indexer with external lattice (in Angstrom) : \"a0x,a0y,a0z,a1x,a1y,a1z,a2x,a2y,a2z\" (9 floats, skips first pass)" << std::endl;
+    std::cout <<
+            "   --force-rotation-lattice <vec>   Force rotation indexer with external lattice (in Angstrom) : \"a0x,a0y,a0z,a1x,a1y,a1z,a2x,a2y,a2z\" (9 floats, skips first pass)"
+            << std::endl;
     std::cout << "   -X, --indexing-algorithm <txt>   Indexing algorithm (FFBIDX|FFT|FFTW|Auto|None)" << std::endl;
-    std::cout << "   -S, --space-group <num>          Space group number - used for both indexing and scaling" << std::endl;
+    std::cout << "   -S, --space-group <num>          Space group number - used for both indexing and scaling" <<
+            std::endl;
     std::cout << "   -C, --unit-cell <cell>           Fix reference unit cell: \"a,b,c,alpha,beta,gamma\"" << std::endl;
-    std::cout << "   -r, --refine <txt>               Geometry refinement algorithm (none|orientation|beam_and_lattice|pixelrefine)" << std::endl;
+    std::cout <<
+            "   -r, --refine <txt>               Geometry refinement algorithm (none|orientation|beam_and_lattice|pixelrefine)"
+            << std::endl;
+    std::cout << std::endl;
+
+    std::cout << "  Integration" << std::endl;
+    std::cout <<
+            "    --integration-use-azim         Background from Bragg peak integration is based on azimuthal integration results"
+            << std::endl;
     std::cout << std::endl;
 
     std::cout << "  Scaling and merging" << std::endl;
-    std::cout << "   -M, --scale-merge               Scale and merge (refine mosaicity) and write scaled.hkl + image.dat" << std::endl;
-    std::cout << "   -P, --partiality <txt>          Partiality refinement fixed|rot|unity (default: fixed)" << std::endl;
+    std::cout <<
+            "   -M, --scale-merge               Scale and merge (refine mosaicity) and write scaled.hkl + image.dat" <<
+            std::endl;
+    std::cout << "   -P, --partiality <txt>          Partiality refinement fixed|rot|unity (default: fixed)" <<
+            std::endl;
     std::cout << "   -A, --anomalous                 Anomalous mode (don't merge Friedel pairs)" << std::endl;
     std::cout << "   -B, --refine-bfactor            Refine per image B-factor" << std::endl;
-    std::cout << "   -w, --wedge[=num]               Refine image wedge during scaling with starting wedge value" << std::endl;
-    std::cout << "   --scaling-high-resolution <num> High resolution limit for spot finding (default: no limit)" << std::endl;
-    std::cout << "   --min-partiality <num>          Minimum partiality to accept reflection (default: 0.02)" << std::endl;
+    std::cout << "   -w, --wedge[=num]               Refine image wedge during scaling with starting wedge value" <<
+            std::endl;
+    std::cout << "   --scaling-high-resolution <num> High resolution limit for spot finding (default: no limit)" <<
+            std::endl;
+    std::cout << "   --min-partiality <num>          Minimum partiality to accept reflection (default: 0.02)" <<
+            std::endl;
     std::cout << "   --min-image-cc <num>            Per-image CC limit in percent (default: no limit)" << std::endl;
-    std::cout << "   --scaling-iterations <num>      Number of scaling iterations with no reference data (default: 3)" << std::endl;
-    std::cout << "   --scaling-output <txt>              Output format for scaling results mtz|cif|txt (default: mtz)" << std::endl;
+    std::cout << "   --scaling-iterations <num>      Number of scaling iterations with no reference data (default: 3)"
+            << std::endl;
+    std::cout << "   --scaling-output <txt>              Output format for scaling results mtz|cif|txt (default: mtz)"
+            << std::endl;
     std::cout << "   -z, --reference-mtz <file>      Reference MTZ file" << std::endl;
     std::cout << std::endl;
 
     std::cout << "  Pixel refinement (experimental, select via -r pixelrefine, needs --reference-mtz)" << std::endl;
-    std::cout << "   --bandwidth <num>               Relative X-ray bandwidth FWHM (e.g. 0.01 for 1% DMM); default from file or 0" << std::endl;
-    }
+    std::cout <<
+            "   --bandwidth <num>               Relative X-ray bandwidth FWHM (e.g. 0.01 for 1% DMM); default from file or 0"
+            << std::endl;
+}
 
 enum {
     OPT_SPOT_SIGMA = 1000,
@@ -92,7 +128,10 @@ enum {
     OPT_SINGLE_PASS_ROTATION,
     OPT_REDO_ROTATION_SPOTS,
     OPT_FORCE_ROTATION_LATTICE,
-    OPT_BANDWIDTH
+    OPT_BANDWIDTH,
+    OPT_INTEGRATION_USE_AZIM,
+    OPT_AZIM_Q_MIN,
+    OPT_AZIM_Q_MAX
 };
 
 static option long_options[] = {
@@ -112,13 +151,16 @@ static option long_options[] = {
     {"wedge", optional_argument, nullptr, 'w'},
     {"scale-merge", no_argument, nullptr, 'M'},
     {"refine", required_argument, nullptr, 'r'},
+    {"azim-q-spacing", required_argument, nullptr, 'Q'},
+    {"azim-q-min", required_argument, nullptr, OPT_AZIM_Q_MIN},
+    {"azim-q-max", required_argument, nullptr, OPT_AZIM_Q_MAX},
 
     {"two-pass-rotation", optional_argument, nullptr, 'R'},
     {"single-pass-rotation", optional_argument, nullptr, OPT_SINGLE_PASS_ROTATION},
     {"redo-rotation-spots", no_argument, nullptr, OPT_REDO_ROTATION_SPOTS},
     {"force-rotation-lattice", required_argument, nullptr, OPT_FORCE_ROTATION_LATTICE},
 
-
+    {"integration-use-azim", no_argument, nullptr, OPT_INTEGRATION_USE_AZIM},
     {"spot-sigma", required_argument, nullptr, OPT_SPOT_SIGMA},
     {"spot-threshold", required_argument, nullptr, OPT_SPOT_THRESHOLD},
     {"spot-high-resolution", required_argument, nullptr, OPT_SPOT_RESOLUTION},
@@ -181,7 +223,6 @@ std::optional<UnitCell> parse_unit_cell_arg(const char *arg) {
         return std::nullopt;
 
 
-
     UnitCell uc{};
     if (!parse_float_strict(parts[0], uc.a)) return std::nullopt;
     if (!parse_float_strict(parts[1], uc.b)) return std::nullopt;
@@ -297,6 +338,11 @@ int main(int argc, char **argv) {
     double min_partiality = 0.02;
     double min_image_cc = 0.0;
     int64_t scaling_iter = 3;
+    std::optional<double> azim_q_spacing;
+    std::optional<double> azim_q_min;
+    std::optional<double> azim_q_max;
+    bool use_azim_for_integration = false;
+
     std::optional<CrystalLattice> forced_rotation_lattice;
 
     std::optional<float> bandwidth_fwhm; // relative FWHM of dlambda/lambda
@@ -317,225 +363,237 @@ int main(int argc, char **argv) {
 
     int opt;
     int option_index = 0;
-    const char *short_opts = "vo:N:s:e:t:R::X:C:z:FABw::S:MP:r:";
+    const char *short_opts = "vo:N:s:e:t:R::X:C:z:FABw::S:MP:r:Q:";
 
-  while ((opt = getopt_long(argc, argv, short_opts, long_options, &option_index)) != -1) {
-    switch (opt) {
-        case 'o':
-            output_prefix = optarg;
-            break;
-        case 'v':
-            verbose = true;
-            break;
-        case 'N':
-            nthreads = atoi(optarg);
-            break;
-        case 's':
-            start_image = atoi(optarg);
-            break;
-        case 'e':
-            end_image = atoi(optarg);
-            break;
-        case 't':
-            image_stride = atoi(optarg);
-            break;
-        case 'R':
-            if (rotation_indexing) {
-                logger.Error("Rotation indexing already enabled");
-                exit(EXIT_FAILURE);
+    while ((opt = getopt_long(argc, argv, short_opts, long_options, &option_index)) != -1) {
+        switch (opt) {
+            case 'o':
+                output_prefix = optarg;
+                break;
+            case 'v':
+                verbose = true;
+                break;
+            case 'N':
+                nthreads = atoi(optarg);
+                break;
+            case 's':
+                start_image = atoi(optarg);
+                break;
+            case 'e':
+                end_image = atoi(optarg);
+                break;
+            case 't':
+                image_stride = atoi(optarg);
+                break;
+            case 'R':
+                if (rotation_indexing) {
+                    logger.Error("Rotation indexing already enabled");
+                    exit(EXIT_FAILURE);
+                }
+                rotation_indexing = true;
+                two_pass_rotation = true;
+                if (optarg)
+                    rotation_indexing_image_count = atoi(optarg);
+
+                break;
+            case 'Q':
+                azim_q_spacing = atof(optarg);
+                break;
+            case OPT_AZIM_Q_MIN:
+                azim_q_min = atof(optarg);
+                break;
+            case OPT_AZIM_Q_MAX:
+                azim_q_max = atof(optarg);
+                break;
+            case OPT_INTEGRATION_USE_AZIM:
+                use_azim_for_integration = true;
+                break;
+            case OPT_SINGLE_PASS_ROTATION:
+                if (rotation_indexing) {
+                    logger.Error("Rotation indexing already enabled");
+                    exit(EXIT_FAILURE);
+                }
+                rotation_indexing = true;
+                two_pass_rotation = false;
+
+                if (optarg)
+                    rotation_indexing_range = atof(optarg);
+                break;
+            case OPT_REDO_ROTATION_SPOTS:
+                reuse_rotation_spots = false;
+                break;
+            case OPT_FORCE_ROTATION_LATTICE: {
+                if (rotation_indexing) {
+                    logger.Error("Rotation indexing already enabled");
+                    exit(EXIT_FAILURE);
+                }
+                rotation_indexing = true;
+
+                auto latt = parse_lattice_arg(optarg);
+                if (!latt.has_value()) {
+                    logger.Error(
+                        "Invalid rotation lattice. Expected: \"a0x,a0y,a0z,a1x,a1y,a1z,a2x,a2y,a2z\" (9 floats, comma-separated). Got: {}",
+                        optarg ? optarg : "<null>");
+                    print_usage();
+                    exit(EXIT_FAILURE);
+                }
+                forced_rotation_lattice = latt;
+                auto uc = latt->GetUnitCell();
+                logger.Info(
+                    "Forced rotation lattice set: a={:.3f} b={:.3f} c={:.3f} alpha={:.3f} beta={:.3f} gamma={:.3f}",
+                    uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma);
+                break;
             }
-            rotation_indexing = true;
-            two_pass_rotation = true;
-            if (optarg)
-                rotation_indexing_image_count = atoi(optarg);
+            case 'X': {
+                std::string alg = optarg ? optarg : "";
+                std::transform(alg.begin(), alg.end(), alg.begin(),
+                               [](unsigned char c) { return static_cast<char>(std::tolower(c)); });
 
-            break;
-        case OPT_SINGLE_PASS_ROTATION:
-            if (rotation_indexing) {
-                logger.Error("Rotation indexing already enabled");
-                exit(EXIT_FAILURE);
+                if (alg == "ffbidx")
+                    indexing_algorithm = IndexingAlgorithmEnum::FFBIDX;
+                else if (alg == "fft")
+                    indexing_algorithm = IndexingAlgorithmEnum::FFT;
+                else if (alg == "fftw")
+                    indexing_algorithm = IndexingAlgorithmEnum::FFTW;
+                else if (alg == "auto")
+                    indexing_algorithm = IndexingAlgorithmEnum::Auto;
+                else if (alg == "none")
+                    indexing_algorithm = IndexingAlgorithmEnum::None;
+                else {
+                    logger.Error("Invalid indexing algorithm: {}", alg);
+                    print_usage();
+                    exit(EXIT_FAILURE);
+                }
+                break;
             }
-            rotation_indexing = true;
-            two_pass_rotation = false;
-
-            if (optarg)
-                rotation_indexing_range = atof(optarg);
-            break;
-        case OPT_REDO_ROTATION_SPOTS:
-            reuse_rotation_spots = false;
-            break;
-        case OPT_FORCE_ROTATION_LATTICE: {
-            if (rotation_indexing) {
-                logger.Error("Rotation indexing already enabled");
-                exit(EXIT_FAILURE);
+            case 'r': {
+                std::string alg = optarg ? optarg : "";
+                std::transform(alg.begin(), alg.end(), alg.begin(),
+                               [](unsigned char c) { return static_cast<char>(std::tolower(c)); });
+                if (alg == "none")
+                    refinement_algorithm = GeomRefinementAlgorithmEnum::None;
+                else if (alg == "beam_and_lattice")
+                    refinement_algorithm = GeomRefinementAlgorithmEnum::BeamCenter;
+                else if (alg == "orientation")
+                    refinement_algorithm = GeomRefinementAlgorithmEnum::OrientationOnly;
+                else if (alg == "pixelrefine")
+                    refinement_algorithm = GeomRefinementAlgorithmEnum::PixelRefine;
+                else {
+                    logger.Error("Invalid geom refinement algorithm: {}", alg);
+                    print_usage();
+                    exit(EXIT_FAILURE);
+                }
+                break;
             }
-            rotation_indexing = true;
-
-            auto latt = parse_lattice_arg(optarg);
-            if (!latt.has_value()) {
-                logger.Error(
-                    "Invalid rotation lattice. Expected: \"a0x,a0y,a0z,a1x,a1y,a1z,a2x,a2y,a2z\" (9 floats, comma-separated). Got: {}",
-                    optarg ? optarg : "<null>");
-                print_usage();
-                exit(EXIT_FAILURE);
+            case 'C': {
+                auto uc = parse_unit_cell_arg(optarg);
+                if (!uc.has_value()) {
+                    logger.Error(
+                        "Invalid unit cell. Expected: \"a,b,c,alpha,beta,gamma\" (6 floats, comma-separated, no spaces). Got: {}",
+                        optarg ? optarg : "<null>");
+                    print_usage();
+                    exit(EXIT_FAILURE);
+                }
+                fixed_reference_unit_cell = uc;
+                logger.Info(
+                    "Fixed reference unit cell set: a={:.3f} b={:.3f} c={:.3f} alpha={:.3f} beta={:.3f} gamma={:.3f}",
+                    uc->a, uc->b, uc->c, uc->alpha, uc->beta, uc->gamma);
+                break;
             }
-            forced_rotation_lattice = latt;
-            auto uc = latt->GetUnitCell();
-            logger.Info(
-                "Forced rotation lattice set: a={:.3f} b={:.3f} c={:.3f} alpha={:.3f} beta={:.3f} gamma={:.3f}",
-                uc.a, uc.b, uc.c, uc.alpha, uc.beta, uc.gamma);
-            break;
-        }
-        case 'X': {
-            std::string alg = optarg ? optarg : "";
-            std::transform(alg.begin(), alg.end(), alg.begin(),
-                           [](unsigned char c) { return static_cast<char>(std::tolower(c)); });
-
-            if (alg == "ffbidx")
-                indexing_algorithm = IndexingAlgorithmEnum::FFBIDX;
-            else if (alg == "fft")
+            case 'z':
+                ref_mtz = optarg;
+                break;
+            case 'F':
                 indexing_algorithm = IndexingAlgorithmEnum::FFT;
-            else if (alg == "fftw")
-                indexing_algorithm = IndexingAlgorithmEnum::FFTW;
-            else if (alg == "auto")
-                indexing_algorithm = IndexingAlgorithmEnum::Auto;
-            else if (alg == "none")
-                indexing_algorithm = IndexingAlgorithmEnum::None;
-            else {
-                logger.Error("Invalid indexing algorithm: {}", alg);
-                print_usage();
-                exit(EXIT_FAILURE);
-            }
-            break;
-        }
-        case 'r': {
-            std::string alg = optarg ? optarg : "";
-            std::transform(alg.begin(), alg.end(), alg.begin(),
-                           [](unsigned char c) { return static_cast<char>(std::tolower(c)); });
-            if (alg == "none")
-                refinement_algorithm = GeomRefinementAlgorithmEnum::None;
-            else if (alg == "beam_and_lattice")
-                refinement_algorithm = GeomRefinementAlgorithmEnum::BeamCenter;
-            else if (alg == "orientation")
-                refinement_algorithm = GeomRefinementAlgorithmEnum::OrientationOnly;
-            else if (alg == "pixelrefine")
-                refinement_algorithm = GeomRefinementAlgorithmEnum::PixelRefine;
-            else {
-                logger.Error("Invalid geom refinement algorithm: {}", alg);
-                print_usage();
-                exit(EXIT_FAILURE);
-            }
-            break;
-        }
-        case 'C': {
-            auto uc = parse_unit_cell_arg(optarg);
-            if (!uc.has_value()) {
-                logger.Error(
-                    "Invalid unit cell. Expected: \"a,b,c,alpha,beta,gamma\" (6 floats, comma-separated, no spaces). Got: {}",
-                    optarg ? optarg : "<null>");
-                print_usage();
-                exit(EXIT_FAILURE);
-            }
-            fixed_reference_unit_cell = uc;
-            logger.Info(
-                "Fixed reference unit cell set: a={:.3f} b={:.3f} c={:.3f} alpha={:.3f} beta={:.3f} gamma={:.3f}",
-                uc->a, uc->b, uc->c, uc->alpha, uc->beta, uc->gamma);
-            break;
-        }
-        case 'z':
-            ref_mtz = optarg;
-            break;
-        case 'F':
-            indexing_algorithm = IndexingAlgorithmEnum::FFT;
-            break;
-        case 'A':
-            anomalous_mode = true;
-            break;
-        case 'B':
-            refine_bfactor = true;
-            break;
-        case 'w':
-            refine_wedge = true;
-            if (optarg)
-                wedge_for_scaling = std::stod(optarg);
-            break;
-        case 'S':
-            space_group_number = atoi(optarg);
-            break;
-        case 'P':
-            if (strcmp(optarg, "unity") == 0)
-                partiality_model = PartialityModel::Unity;
-            else if (strcmp(optarg, "fixed") == 0)
-                partiality_model = PartialityModel::Fixed;
-            else if (strcmp(optarg, "rot") == 0)
-                partiality_model = PartialityModel::Rotation;
-            else {
-                logger.Error("Invalid partiality mode: {}", optarg);
-                print_usage();
-                exit(EXIT_FAILURE);
-            }
-            break;
-        case OPT_SPOT_SIGMA:
-            sigma_spot_finding = atof(optarg);
-            logger.Info("Noise threshold level for spot finding set to {:.2f} sigma", sigma_spot_finding);
-            break;
-        case OPT_SPOT_THRESHOLD:
-            photon_count_threshold_spot_finding = atoi(optarg);
-            logger.Info("Photon-count threshold level for spot finding set to {:d}",
-                        photon_count_threshold_spot_finding);
-            break;
-        case OPT_SPOT_RESOLUTION:
-            d_min_spot_finding = atof(optarg);
-            logger.Info("High resolution limit for spot finding set to {:.2f} A", d_min_spot_finding);
-            break;
-        case OPT_MAX_SPOTS:
-            max_spot_count_override = atoll(optarg);
-            logger.Info("Max spot count overridden to {}", max_spot_count_override.value());
-            break;
-        case 'M':
-            run_scaling = true;
-            break;
-        case OPT_MIN_PARTIALITY:
-            min_partiality = std::stod(optarg);
-            break;
-        case OPT_MIN_IMAGE_CC:
-            min_image_cc = std::stod(optarg);
-            break;
-        case OPT_SCALING_HIGH_RESOLUTION:
-            d_min_scale_merge = atof(optarg);
-            break;
-        case OPT_SCALING_OUTPUT:
-            if (strcmp(optarg, "mtz") == 0) {
-                intensity_format = IntensityFormat::MTZ;
-            } else if (strcmp(optarg, "cif") == 0) {
-                intensity_format = IntensityFormat::mmCIF;
-            } else if (strcmp(optarg, "txt") == 0) {
-                intensity_format = IntensityFormat::Text;
-            } else {
-                logger.Error("Invalid intensity format: {}", optarg);
-                exit(EXIT_FAILURE);
-            }
-            break;
-        case OPT_SCALING_ITERATIONS:
-            scaling_iter = atoi(optarg);
-            if (scaling_iter <= 0) {
-                logger.Error("Invalid scaling iteration count: {}", scaling_iter);
-                exit(EXIT_FAILURE);
-            }
-            break;
-        case OPT_BANDWIDTH:
-            bandwidth_fwhm = atof(optarg);
-            if (!(bandwidth_fwhm.value() >= 0.0f)) {
-                logger.Error("Invalid bandwidth: {}", optarg);
-                exit(EXIT_FAILURE);
-            }
-            break;
+                break;
+            case 'A':
+                anomalous_mode = true;
+                break;
+            case 'B':
+                refine_bfactor = true;
+                break;
+            case 'w':
+                refine_wedge = true;
+                if (optarg)
+                    wedge_for_scaling = std::stod(optarg);
+                break;
+            case 'S':
+                space_group_number = atoi(optarg);
+                break;
+            case 'P':
+                if (strcmp(optarg, "unity") == 0)
+                    partiality_model = PartialityModel::Unity;
+                else if (strcmp(optarg, "fixed") == 0)
+                    partiality_model = PartialityModel::Fixed;
+                else if (strcmp(optarg, "rot") == 0)
+                    partiality_model = PartialityModel::Rotation;
+                else {
+                    logger.Error("Invalid partiality mode: {}", optarg);
+                    print_usage();
+                    exit(EXIT_FAILURE);
+                }
+                break;
+            case OPT_SPOT_SIGMA:
+                sigma_spot_finding = atof(optarg);
+                logger.Info("Noise threshold level for spot finding set to {:.2f} sigma", sigma_spot_finding);
+                break;
+            case OPT_SPOT_THRESHOLD:
+                photon_count_threshold_spot_finding = atoi(optarg);
+                logger.Info("Photon-count threshold level for spot finding set to {:d}",
+                            photon_count_threshold_spot_finding);
+                break;
+            case OPT_SPOT_RESOLUTION:
+                d_min_spot_finding = atof(optarg);
+                logger.Info("High resolution limit for spot finding set to {:.2f} A", d_min_spot_finding);
+                break;
+            case OPT_MAX_SPOTS:
+                max_spot_count_override = atoll(optarg);
+                logger.Info("Max spot count overridden to {}", max_spot_count_override.value());
+                break;
+            case 'M':
+                run_scaling = true;
+                break;
+            case OPT_MIN_PARTIALITY:
+                min_partiality = std::stod(optarg);
+                break;
+            case OPT_MIN_IMAGE_CC:
+                min_image_cc = std::stod(optarg);
+                break;
+            case OPT_SCALING_HIGH_RESOLUTION:
+                d_min_scale_merge = atof(optarg);
+                break;
+            case OPT_SCALING_OUTPUT:
+                if (strcmp(optarg, "mtz") == 0) {
+                    intensity_format = IntensityFormat::MTZ;
+                } else if (strcmp(optarg, "cif") == 0) {
+                    intensity_format = IntensityFormat::mmCIF;
+                } else if (strcmp(optarg, "txt") == 0) {
+                    intensity_format = IntensityFormat::Text;
+                } else {
+                    logger.Error("Invalid intensity format: {}", optarg);
+                    exit(EXIT_FAILURE);
+                }
+                break;
+            case OPT_SCALING_ITERATIONS:
+                scaling_iter = atoi(optarg);
+                if (scaling_iter <= 0) {
+                    logger.Error("Invalid scaling iteration count: {}", scaling_iter);
+                    exit(EXIT_FAILURE);
+                }
+                break;
+            case OPT_BANDWIDTH:
+                bandwidth_fwhm = atof(optarg);
+                if (!(bandwidth_fwhm.value() >= 0.0f)) {
+                    logger.Error("Invalid bandwidth: {}", optarg);
+                    exit(EXIT_FAILURE);
+                }
+                break;
 
-        default:
-            print_usage();
-            exit(EXIT_FAILURE);
+            default:
+                print_usage();
+                exit(EXIT_FAILURE);
+        }
     }
-}
 
     if (optind != argc - 1) {
         logger.Error("Input file not specified");
@@ -603,7 +661,7 @@ int main(int argc, char **argv) {
 
     if (images_to_process <= 0) {
         logger.Warning("No images to process (Start: {}, End: {} Stride: {}, Total: {})", start_image, end_image,
-            image_stride, total_images_in_file);
+                       image_stride, total_images_in_file);
         return 0;
     }
 
@@ -668,6 +726,17 @@ int main(int argc, char **argv) {
 
     experiment.ImportScalingSettings(scaling_settings);
 
+    BraggIntegrationSettings bragg_integration_settings;
+    bragg_integration_settings.UseAzimProfile(use_azim_for_integration);
+    experiment.ImportBraggIntegrationSettings(bragg_integration_settings);
+
+    AzimuthalIntegrationSettings azimuthal_integration_settings;
+    if (azim_q_spacing)
+        azimuthal_integration_settings.QSpacing_recipA(azim_q_spacing.value());
+    azimuthal_integration_settings.QRange_recipA(azim_q_min.value_or(azimuthal_integration_settings.GetLowQ_recipA()),
+                                                 azim_q_max.value_or(azimuthal_integration_settings.GetHighQ_recipA()));
+    experiment.ImportAzimuthalIntegrationSettings(azimuthal_integration_settings);
+
     SpotFindingSettings spot_settings;
     spot_settings.enable = true;
     spot_settings.indexing = true;
@@ -680,9 +749,6 @@ int main(int argc, char **argv) {
     // Initialize Analysis Components
     PixelMask pixel_mask = dataset->pixel_mask;
 
-    // If dataset has a mask you wish to use, you might need to load it into pixel_mask here
-    // e.g. pixel_mask.LoadUserMask(dataset->pixel_mask, ...);
-
     AzimuthalIntegrationMapping mapping(experiment, pixel_mask);
     IndexerThreadPool indexer_pool(experiment.GetIndexingSettings());
 
@@ -768,7 +834,7 @@ int main(int argc, char **argv) {
 
                     analysis.Analyze(msg, profile, first_pass_spot_settings);
                 }
-               indexer.AddImageToRotationIndexer(msg);
+                indexer.AddImageToRotationIndexer(msg);
             } catch (const std::exception &e) {
                 logger.Warning("Failed to add image {} to first-pass rotation indexing: {}", image_idx, e.what());
             }
@@ -814,7 +880,8 @@ int main(int argc, char **argv) {
             msg.image = img->image;
             msg.number = image_ordinal;
             msg.original_number = image_idx;
-            if (dataset->efficiency.size() > image_idx) msg.image_collection_efficiency = dataset->efficiency[image_idx];
+            if (dataset->efficiency.size() > image_idx)
+                msg.image_collection_efficiency = dataset->efficiency[image_idx];
 
             total_uncompressed_bytes += msg.image.GetUncompressedSize();
 
@@ -921,8 +988,8 @@ int main(int argc, char **argv) {
     if (consensus_cell) {
         logger.Info("Consensus unit cell found in {:.2f} ms", consensus_duration * 1e3);
         logger.Info("UC:         a={:.2f} b={:.2f} c={:.2f} alpha={:.2f} beta={:.2f} gamma={:.2f}",
-                        consensus_cell->a, consensus_cell->b, consensus_cell->c,
-                        consensus_cell->alpha, consensus_cell->beta, consensus_cell->gamma);
+                    consensus_cell->a, consensus_cell->b, consensus_cell->c,
+                    consensus_cell->alpha, consensus_cell->beta, consensus_cell->gamma);
     } else
         logger.Info("Consensus unit cell not found - calculation tool {:.2f} ms", consensus_duration * 1e3);
     end_msg.unit_cell = consensus_cell;
@@ -930,9 +997,8 @@ int main(int argc, char **argv) {
     if (end_msg.indexing_rate.has_value()
         && end_msg.indexing_rate > 0
         && (run_scaling || !reference_data.empty())) {
-
         const bool pixel_refine_path =
-            (refinement_algorithm == GeomRefinementAlgorithmEnum::PixelRefine);
+                (refinement_algorithm == GeomRefinementAlgorithmEnum::PixelRefine);
 
         // Scaling is only the classical, no-reference ScaleOnTheFly post-pass.
         // - With a reference (classical path): per-image live scaling already ran.
@@ -973,11 +1039,12 @@ int main(int argc, char **argv) {
         MergeOnTheFly merge_engine(experiment);
         if (consensus_cell.has_value())
             merge_engine.ReferenceCell(*consensus_cell);
-        for (auto &i : indexer.GetIntegrationOutcome())
+        for (auto &i: indexer.GetIntegrationOutcome())
             merge_engine.AddImage(i);
 
         auto merged_reflections = merge_engine.ExportReflections();
-        auto merged_statistics = merge_engine.MergeStats(merged_reflections, indexer.GetIntegrationOutcome(), reference_data);
+        auto merged_statistics = merge_engine.MergeStats(merged_reflections, indexer.GetIntegrationOutcome(),
+                                                         reference_data);
 
         auto merge_end = std::chrono::steady_clock::now();
         double merge_time = std::chrono::duration<double>(merge_end - merge_start).count();
@@ -1031,17 +1098,17 @@ int main(int argc, char **argv) {
 
     auto image_mean_time = plots.GetMeanProcessingTime();
     std::cout << fmt::format(
-        "Per-image time: (mean; milliseconds): decompress {:.2f} preprocess {:.2f} azint {:.2f} spot finding {:.2f} indexing {:.2f} refinement {:.2f} indexing analysis {:.2f} prediction {:.2f} integration {:.2f} scaling {:.2f} total {:.2f}",
-        image_mean_time.compression * 1e3,
-        image_mean_time.preprocessing * 1e3,
-        image_mean_time.azint * 1e3,
-        image_mean_time.spot_finding * 1e3,
-        image_mean_time.indexing * 1e3,
-        image_mean_time.refinement * 1e3,
-        image_mean_time.indexing_analysis * 1e3,
-        image_mean_time.bragg_prediction * 1e3,
-        image_mean_time.integration * 1e3,
-        image_mean_time.image_scale * 1e3,
-        image_mean_time.processing * 1e3)
-    << std::endl;
+                "Per-image time: (mean; milliseconds): decompress {:.2f} preprocess {:.2f} azint {:.2f} spot finding {:.2f} indexing {:.2f} refinement {:.2f} indexing analysis {:.2f} prediction {:.2f} integration {:.2f} scaling {:.2f} total {:.2f}",
+                image_mean_time.compression * 1e3,
+                image_mean_time.preprocessing * 1e3,
+                image_mean_time.azint * 1e3,
+                image_mean_time.spot_finding * 1e3,
+                image_mean_time.indexing * 1e3,
+                image_mean_time.refinement * 1e3,
+                image_mean_time.indexing_analysis * 1e3,
+                image_mean_time.bragg_prediction * 1e3,
+                image_mean_time.integration * 1e3,
+                image_mean_time.image_scale * 1e3,
+                image_mean_time.processing * 1e3)
+            << std::endl;
 }