cabouman / mbircone

BSD 3-Clause "New" or "Revised" License
11 stars 9 forks source link

Restructure priorWeight_QGGMRF/priorWeight_proxMap flags. #63

Closed dyang37 closed 2 years ago

dyang37 commented 2 years ago

Currently in C code we have two flags to control whether recon() function is in qGGMRF mode or proximal map mode: In struct ReconParams:

float priorWeight_QGGMRF;       /* Prior mode: (-1: QGGMRF mode off; 1: QGGMRF mode on) */
float priorWeight_proxMap;      /* Prior mode: (-1: prox mode off; 1: prox mode on) */

There are three issues related to this:

  1. These flags need to be changed to int or bool instead of float.
  2. If it should be the case that only one of qGGMRF or prox modes may be turned on at the same time, then we may merge two flags into one bool flag.
  3. When calculating theta_1 and theta_2 values in ICD update, we have the update code as follows:
    float theta1, theta2;
    theta1 = icdInfo->theta1_f + reconParams->priorWeight_QGGMRF*icdInfo->theta1_p_QGGMRF + reconParams->priorWeight_proxMap*icdInfo->theta1_p_proxMap;
    theta2 = icdInfo->theta2_f + reconParams->priorWeight_QGGMRF*icdInfo->theta2_p_QGGMRF + reconParams->priorWeight_proxMap*icdInfo->theta2_p_proxMap;

    If we follow the definition that theta1_p_QGGMRF and priorWeight_proxMap is either -1 or 1, then this could be wrong since both qGGMRF and proxMap prior weights are added to theta values.

dyang37 commented 2 years ago

Regarding issue 3 the update rule of theta values mentioned above, I performed some investigation on the theta update rules.

When in qGGMRF mode, the variables icdInfo->theta1_p_proxMap and icdInfo->theta2_p_proxMap always have values of 0.

priorWeight_QGGMRF=1.0
theta_1_p_QGGMRF=1.4
theta_2_p_QGGMRF=1585.7
priorWeight_proxMap=-1.0
theta_1_p_proxMap=0.0
theta_2_p_proxMap=0.0

priorWeight_QGGMRF=1.0
theta_1_p_QGGMRF=102.3
theta_2_p_QGGMRF=1931.3
priorWeight_proxMap=-1.0
theta_1_p_proxMap=0.0
theta_2_p_proxMap=0.0

Therefore for update rules of theta values, the part reconParams->priorWeight_proxMap*icdInfo->theta1_p_proxMap would be 0 in qGGMRF mode even if priorWeight_proxMap is not 0.

Likewise, when in prox mode, icdInfo->theta1_p_qGGMRF and icdInfo->theta2_p_qGGMRF are all 0, and as a result, the part reconParams->priorWeight_QGGMRF*icdInfo->theta1_p_QGGMRF will be 0 (or approximately 0, given that both variables are float).

priorWeight_QGGMRF=-1.0
theta_1_p_QGGMRF=0.0
theta_2_p_QGGMRF=0.0
priorWeight_proxMap=1.0
theta_1_p_proxMap=2.9
theta_2_p_proxMap=242.0

priorWeight_QGGMRF=-1.0
theta_1_p_QGGMRF=0.0
theta_2_p_QGGMRF=0.0
priorWeight_proxMap=1.0
theta_1_p_proxMap=1.3
theta_2_p_proxMap=242.0

Therefore, our current code yields the correct results, even though the prior flags are poorly designed. Given that we are not planning to support the mode where both qGGMRF and proximal map priors are enabled at the same time, I propose the following changes:

  1. In python and C code, use a single bool variable prox_mode to replace the two float flags we have, where prox_mode=False corresponds to qGGMRF mode, and prox_mode=True corresponds to proximal map mode.
  2. In the code, whenever we need to check prior mode, replace if(reconParams->priorWeight_QGGMRF >= 0) with if(!prox_mode), and replace (reconParams->priorWeight_proxMap >= 0) with if(prox_mode).
  3. Update rule to theta values does not need to be changed, but if we are really uncomfortable with the idea of multiplying float with bool, we could modify the code to be an if-else statement.

(We don't have to change the code now, but at some point we need to redesign these flags)

cabouman commented 2 years ago

Let's discuss this. I think we need to remove the floating point flags and replace them with a single flag of prox_mode as you suggested.

dyang37 commented 2 years ago

Let's discuss this. I think we need to remove the floating point flags and replace them with a single flag of prox_mode as you suggested.

I'll go ahead and make the modification then, probably after adding comments to C code.

gbuzzard commented 2 years ago

I would prefer to change the update rule to be if-then to make the intention clearer.

dyang37 commented 2 years ago

I would prefer to change the update rule to be if-then to make the intention clearer.

Yeah that makes sense to me.

dyang37 commented 2 years ago

This is implemented and tested in branch c_comment. I tested with demo_mace3D_old.py. The recon results of both mace3D and qGGMRF are visually the same as before. We could probably PR the branch c_comment after resolving the previous PR regarding auto_sigma_prior_fix branch. One minor issue is that instead of bool I used int for the flag prox_mode, and the reason is that Cython interface does not recognize bool as a valid variable type. I saw that in svmbir we used char to define flags instead. I do not quite understand the logic behind that, but we could use char if that's the standard way of defining Boolean variables in Cython.

cabouman commented 2 years ago

Yes, I think it is better style to use char since that is just one byte. Yes, go ahead and finish this, and then do a PR.

dyang37 commented 2 years ago

We merged this into master last Friday.