Changeset 38136

Timestamp:

Apr 16, 2015, 4:19:35 PM (11 years ago)

Author:

eugene

Message:

add code to calculate Mxx,Mxy,Myy -> E1,E2 for each object; calculate mean R7 values

File:

• branches/eam_branches/ipp-20150405/Ohana/src/dvolens/src/update_objects_catalog.c (modified) (7 diffs)

branches/eam_branches/ipp-20150405/Ohana/src/dvolens/src/update_objects_catalog.c

@@ -11 +11 @@
   int Nsecfilt = GetPhotcodeNsecfilt ();
   
+  float *Mxx_obj = NULL;
+  float *Mxy_obj = NULL;
+  float *Myy_obj = NULL;
+
+  ALLOCATE (Mxx_obj, float, Nsecfilt);
+  ALLOCATE (Mxy_obj, float, Nsecfilt);
+  ALLOCATE (Myy_obj, float, Nsecfilt);
+
   // I think we have already allocated lensobj
   REALLOCATE (catalog->lensobj, Lensobj, catalog->Naverage * Nsecfilt);
@@ -27 +35 @@
     for (j = 0; j < Nsecfilt; j++) {
       dvo_lensobj_init (&lensobj[j], TRUE); // init accumulated values to 0
+      Mxx_obj[j] = 0.0;
+      Mxy_obj[j] = 0.0;
+      Myy_obj[j] = 0.0;
     }
 
     // reset the average values
-
     off_t Loff = average->lensingOffset;
     off_t Moff = average->measureOffset;
@@ -67 +77 @@
       if (measure->psfQFperf < 0.85) continue; 
       
-
+      Mxx_obj[Nsec] += measure->Mxx;
+      Mxy_obj[Nsec] += measure->Mxy;
+      Myy_obj[Nsec] += measure->Myy;
+      
       lensobj[Nsec].X11_sm_obj += lensing->X11_sm_obj;
       lensobj[Nsec].X12_sm_obj += lensing->X12_sm_obj;
@@ -101 +114 @@
       lensobj[Nsec].sF_ApR6 += SQ(lensing->sF_ApR6);
       lensobj[Nsec].fF_ApR6 += lensing->fF_ApR6;
+
+      lensobj[Nsec]. F_ApR7 += lensing-> F_ApR7;
+      lensobj[Nsec].dF_ApR7 += SQ(lensing->dF_ApR7);
+      lensobj[Nsec].sF_ApR7 += SQ(lensing->sF_ApR7);
+      lensobj[Nsec].fF_ApR7 += lensing->fF_ApR7;
 
       lensobj[Nsec].Nmeas   ++;
@@ -139 +157 @@
       lensobj[j].sF_ApR5 /= Nmeas;
       lensobj[j].fF_ApR5 /= Nmeas;
+
       lensobj[j]. F_ApR6 /= Nmeas;
       lensobj[j].dF_ApR6 /= Nmeas;
@@ -144 +163 @@
       lensobj[j].fF_ApR6 /= Nmeas;
 
+      lensobj[j]. F_ApR7 /= Nmeas;
+      lensobj[j].dF_ApR7 /= Nmeas;
+      lensobj[j].sF_ApR7 /= Nmeas;
+      lensobj[j].fF_ApR7 /= Nmeas;
+
       lensobj[j].dF_ApR5 = sqrt(lensobj[j].dF_ApR5);
       lensobj[j].sF_ApR5 = sqrt(lensobj[j].sF_ApR5);
       lensobj[j].dF_ApR6 = sqrt(lensobj[j].dF_ApR6);
       lensobj[j].sF_ApR6 = sqrt(lensobj[j].sF_ApR6);
+      lensobj[j].dF_ApR7 = sqrt(lensobj[j].dF_ApR7);
+      lensobj[j].sF_ApR7 = sqrt(lensobj[j].sF_ApR7);
 
       // somewhere in here I should calculate the values of E1, E2, and gamma
+
+      float e0 = Mxx_obj[j] + Myy_obj[j];
+      lensobj[j].E1 = (Mxx_obj[j] - Myy_obj[j]) / e0;
+      lensobj[j].E2 = 2.0 * Mxy_obj[j] / e0;
     }
 
@@ -159 +189 @@
   }
 
+  free (Mxx_obj);
+  free (Mxy_obj);
+  free (Myy_obj);
+
   catalog->Nlensobj = Nlensobj;
   catalog->Nlensobj_disk = 0;