Skip to content

Autocorrclean#1

Open
hmuellergoe wants to merge 19 commits intoARDG-NRAO:mainfrom
hmuellergoe:autocorrclean
Open

Autocorrclean#1
hmuellergoe wants to merge 19 commits intoARDG-NRAO:mainfrom
hmuellergoe:autocorrclean

Conversation

@hmuellergoe
Copy link

Major code in this request:
-> Implementation of autocorr-CLEAN

includes also a few, minor updates and some experimental code for other methods:

  • option to work with odd number of pixels in Asp-CLEAN
  • expose LBFGSB quantities (tolerances, exit criteria) to the top level
  • expose scales for Asp-CLEAN to the top level
  • add an option to never switch to Hogbom in Asp-CLEAN when fusedthreshold is set to a negative number
  • experimental code for approximating the psf as a sum of wavelets, to accelerate Asp-CLEAN LBFGS computation (that is like the Zhang version, but instead of a simple Gaussian approximating the main lobe, we use a few Gaussians approximating the main lobe and the first few sidelobes), not robustly tested yet
  • experimental alternative multifrequency Asp-CLEAN option: in which the Asp's are defined from a continuum image. Then the spectral imaging bypasses the iterative CLEAN part of Asp-CLEAN, simply reads in the continuum Asps, and for every spectral channel just solves the LBFGSB step (i.e. adapts the Asps locally to what is needed in a respective spectral channel), not robustly tested now
  • added options to save and load helper products

Hendrik Muller added 6 commits December 2, 2025 16:12
…-CLEAN (not robustly working yet), and solve an inconsistency in the normalization
…ents than for extended components, gains ~ 0.8 are fine for extended components, but could lead to divergence if the algorithms moves to point source model earlier than expected
@preshanth preshanth self-assigned this Feb 18, 2026
@preshanth
Copy link
Contributor

@hmuellergoe and @agawatw the building is failing on casacpp which means some of your changes have not propagated into this branch. Can you make sure the code compiles before we carry out a review please.

@preshanth
Copy link
Contributor

@hmuellergoe Once you fix the humbee tests we should be in a place to begin a proper review. @agawatw please detail any comments you have here.

@preshanth
Copy link
Contributor

@hmuellergoe the tests are still failing for both humbee and asp

Copy link
Contributor

@agawatw agawatw left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There should be better software architecture to implement the experimental features without touching the base AspClean class API for the additional parameters, like waveletXXX. The Asp code seems harder for maintenance now because of the new mfasp conditional loop.

Also, please address the following comments and the duplicated code for the fusedthreshold.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@hmuellergoe These added parameters seem fine but why do we need to expose lbfgsbepsf/lbfgsbepsx/lbfgsbepsg/lbfgsmaxit (LBFGSB optimization parameters of Asp-CLEAN)? Is it for the purpose of CG-Clean or autocorrelation, or just for better Asp-Clean performance if users know how to fine tune them?

itsUserLargestScale(largestScale),
itsMCsetup(true),
itsFusedThreshold(fusedThreshold),
itsScales(scales),
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why does "scales" need to be provided to Asp-Clean since the advantage of Asp-Clean is that it can automatically determined the optimal scale without users' input? Is it for computational efficiency? I see you added a new function of "loadInitScaleXfrs(itsScales);", and is the "scales" for that purpose?

If it is added for efficiency, have you measured the performance improvement to justify the added parameter?

// float initModelFlux = 0.0;

itsdimensionsareeven = (psfShape_p(0) == 2*(psfShape_p(0)/2));
if (itsFusedThreshold == -1.){
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Line 272 to 279 are duplicated functionality since this has already been implemented in AspClean. See the function switchedToHogbom(), where
if (itsFusedThreshold < 0) itsSwitchedToHogbom = false;
The tclean user guide has already mentioned that if fusedthrehold=-1, it never switches to hogbom.


if (itsFusedThreshold < 0)
{
os << LogIO::WARN << "Acceptable fusedthreshld values are >= 0. Changing fusedthreshold from " << itsFusedThreshold << " to -1." << LogIO::POST;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This warning message is wrong. See the other comment about the already existed functionality that fusedthreshold=-1 is valid which tells AspClean never switches to hogbom. In this if(itsFusedThreshold < 0) condition, when fusedthreshold=-1, it prints this warning message which is incorrect.

itsNumNoChange = 0;
}
if (!itsSwitchedToHogbom && itsNumNoChange >= 2)
if (!itsSwitchedToHogbom && itsNumNoChange >= 2 && itsstopMS)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

&& itsstopMS here is a duplicated addition as well since there is already !itsSwitchedToHogbom in this if condition.

os << "Failed to reach stopping threshold" << LogIO::POST;
}

write_array(Optimums, std::string("./strengthoptimum"));
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If the added code in MatrixCleaner.cc is purely for debug information, we should not check it in. Or at least, not place them in MatrixCleaner.cc which is the base class of every deconvolver.

}

if(itsSwitchedToHogbom){
Optimums(ii) = 0.0;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think when it has been "switched to hogbom" (i.e. using only delta scale), the optimal strength should be set to 0 (i.e. Optimums(ii) = itsStrengthOptimum;.
Also, if the purpose of having these 3 additional variables
casacore::Vector<casacore::Float> Optimums(itsMaxNiter); casacore::Vector<casacore::Float> ScaleSizes(itsMaxNiter); casacore::Vector<casacore::IPosition> Positions(itsMaxNiter);
is purely for debug information, it is not necessary since there are already variables in AspClean holding this information already.

@sanbw
Copy link
Contributor

sanbw commented Mar 13, 2026 via email

@hmuellergoe
Copy link
Author

That would be probably helpful. However, I already refactored most of it today by implemeting the spectal and wavelet option simply as an inherited subclass of Asp-CLEAN, leaving the original Asp-CLEAN code untouched, looks most practical to me from a maintenance point of view.

Copy link
Contributor

@agawatw agawatw left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@hmuellergoe The recommended software architecture is as followed, which will eliminates all if(itsWaveletTrigger) checks and wavelet-specific parameters from the base implementation and facilitate future experimental work and better maintainability.

  1. Keep the base algorithm (AspClean) stable while improving it for extensibility.
  2. Move experimental behavior into a subclass (e.g. waveletAspClean)
  3. Use interchangeable strategies for scale models and normalization.

I'll work on item 1 and 3 above which will have usage examples of how to create a new scale model and normalization method in the subclass that inherits from AspClean.

For example, in terms of Item (3), there will be a new ScaleModel:
class ScaleModel {
public:
virtual ~ScaleModel() {}
virtual float value(
int i, int j,
int refi, int refj,
double scaleSize) const = 0;
};

In Asp, the gaussian scale model initialization will be changed to:
class GaussianScaleModel : public ScaleModel {
public:
double value(
int i, int j,
int refi, int refj,
float scaleSize) const override
{
return (1.0/(sqrt(2*M_PI)*scaleSize)) *
exp(-(pow(i-refi,2)+pow(j-refj,2)) *
0.5 / pow(scaleSize,2));
}
};

with a new function to pass this to the existing Asp function that sets up the scales:
void setScaleModel(std::unique_ptr model)
{
itsScaleModel = std::move(model);
}

Then, in the existing Asp's makeScaleImage function, the code will be changed from
if(itsWaveletTrigger == false){
iscale(i,j) = (1.0/(sqrt(2*M_PI)*scaleSize))exp(-(pow(i-refi,2) + pow(j-refj,2))0.5/pow(scaleSize,2)); // gaussian
}
else{
iscale(i,j) = (1.0/(2
M_PI
pow(scaleSize,2)))*exp(-(pow(i-refi,2) + pow(j-refj,2))*0.5/pow(scaleSize,2)); // wavelet
}

to

iscale(i,j) =
itsScaleModel->value(i,j,refi,refj,scaleSize);

This is just an example change. You'd need to implement your own
class WaveletScaleModel : public ScaleModel {}.

Also, the new WaveletAspClean would be something like:
class WaveletAspClean : public AspClean {
private:
double itsLbfgsEpsF;
double itsLbfgsEpsX;
double itsLbfgsEpsG;
int itsLbfgsMaxit;
casacore::Vector itsWaveletScales;
casacore::Vector itsWaveletAmps;
public:
WaveletAspClean()
{
setScaleModel(std::make_unique < WaveletScaleModel > ());
}
};

So there will be no wavelet parameters exist in AspClean.

Hendrik Muller added 2 commits March 15, 2026 19:28
…ation of ImageInterface to Asp-CLEAN such that building again
…te: I will need to bring back that code later to make the spectral option work
@preshanth
Copy link
Contributor

@hmuellergoe the code upstream has been updated a lot. It is time to rebase main onto your code and then come back and push here with --force. Then @agawatw and @preshanth will approve the PR and then I will grab this code and commit history and make sure it comes in from gitlab.

@agawatw
Copy link
Contributor

agawatw commented Mar 19, 2026

@hmuellergoe I wrote a README in synthesis/MeasurementEquations that describes how to derive from the AspClean base class with examples. Once you have updated your code accordingly, please check if you still get the same WaveletAsp and auto_CS results as before. We''ll then do the PR and I'll add unit tests for those in the ardg/main gitlab.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants