Changeset 33963 for trunk/psModules/src/objects/pmSourcePhotometry.c
Legend:
- Unmodified
- Added
- Removed
-
trunk
-
trunk/psModules
- Property svn:mergeinfo changed
/branches/eam_branches/ipp-20111122/psModules merged: 33044,33085,33087,33096,33638 /branches/eam_branches/ipp-20120405/psModules (added) merged: 33947,33951
- Property svn:mergeinfo changed
-
trunk/psModules/src/objects/pmSourcePhotometry.c
r32998 r33963 899 899 } 900 900 901 double pmSourceModelWeight(const pmSource *Mi, int term, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)901 double pmSourceModelWeight(const pmSource *Mi, int term, const pmSourceFitVarMode fitVarMode, const float covarFactor, psImageMaskType maskVal) 902 902 { 903 903 PS_ASSERT_PTR_NON_NULL(Mi, NAN); 904 double flux = 0, wt = 0, factor = 0;904 double flux = 0, wt = 1.0, factor = 0; 905 905 906 906 const psImage *Pi = Mi->modelFlux; 907 907 assert (Pi != NULL); 908 const psImage *Wi = Mi->variance; 909 if (!unweighted_sum) { 910 assert (Wi != NULL); 911 } 908 909 const psImage *Wi = NULL; 910 switch (fitVarMode) { 911 case PM_SOURCE_PHOTFIT_CONST: 912 break; 913 case PM_SOURCE_PHOTFIT_IMAGE_VAR: 914 Wi = Mi->variance; 915 psAssert (Wi, "programming error"); 916 break; 917 case PM_SOURCE_PHOTFIT_MODEL_VAR: 918 Wi = Mi->modelVar; 919 psAssert (Wi, "programming error"); 920 break; 921 case PM_SOURCE_PHOTFIT_NONE: 922 psAbort("programming error"); 923 } 924 912 925 const psImage *Ti = Mi->maskObj; 913 926 assert (Ti != NULL); … … 915 928 for (int yi = 0; yi < Pi->numRows; yi++) { 916 929 for (int xi = 0; xi < Pi->numCols; xi++) { 917 if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi] & maskVal) 918 continue; 919 if (!unweighted_sum) { 930 if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi] & maskVal) continue; 931 if (fitVarMode != PM_SOURCE_PHOTFIT_CONST) { 920 932 wt = covarFactor * Wi->data.F32[yi][xi]; 921 if (wt == 0) 922 continue; 933 if (wt == 0) continue; 923 934 } 924 935 … … 937 948 } 938 949 939 if (unweighted_sum) { 940 flux += (factor * Pi->data.F32[yi][xi]); 941 } else { 942 flux += (factor * Pi->data.F32[yi][xi]) / wt; 943 } 950 // wt is 1.0 for CONST 951 flux += (factor * Pi->data.F32[yi][xi]) / wt; 944 952 } 945 953 } … … 947 955 } 948 956 949 double pmSourceModelDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)957 double pmSourceModelDotModel (const pmSource *Mi, const pmSource *Mj, const pmSourceFitVarMode fitVarMode, const float covarFactor, psImageMaskType maskVal) 950 958 { 951 959 PS_ASSERT_PTR_NON_NULL(Mi, NAN); … … 955 963 int xIs, xJs, yIs, yJs; 956 964 int xIe, yIe; 957 double flux, wt; 965 double flux; 966 double wt = 1.0; 958 967 959 968 const psImage *Pi = Mi->modelFlux; … … 962 971 assert (Pj != NULL); 963 972 964 const psImage *Wi = Mi->variance; 965 if (!unweighted_sum) { 966 assert (Wi != NULL); 967 } 973 const psImage *Wi = NULL; 974 switch (fitVarMode) { 975 case PM_SOURCE_PHOTFIT_CONST: 976 break; 977 case PM_SOURCE_PHOTFIT_IMAGE_VAR: 978 Wi = Mi->variance; 979 psAssert (Wi, "programming error"); 980 break; 981 case PM_SOURCE_PHOTFIT_MODEL_VAR: 982 Wi = Mi->modelVar; 983 psAssert (Wi, "programming error"); 984 break; 985 case PM_SOURCE_PHOTFIT_NONE: 986 psAbort("programming error"); 987 } 968 988 969 989 const psImage *Ti = Mi->maskObj; … … 995 1015 continue; 996 1016 997 // XXX skip the nonsense weight pixels? 998 if (unweighted_sum) { 999 flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]); 1000 } else { 1017 float value = (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]); 1018 switch (fitVarMode) { 1019 case PM_SOURCE_PHOTFIT_CONST: 1020 wt = 1.0; 1021 break; 1022 case PM_SOURCE_PHOTFIT_IMAGE_VAR: 1023 case PM_SOURCE_PHOTFIT_MODEL_VAR: 1001 1024 wt = covarFactor * Wi->data.F32[yi][xi]; 1002 if (wt > 0) { 1003 flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]) / wt; 1004 } 1005 } 1025 break; 1026 case PM_SOURCE_PHOTFIT_NONE: 1027 psAbort("programming error"); 1028 } 1029 // skip pixels with nonsense weight values 1030 if (wt <= 0) continue; 1031 1032 flux += value / wt; 1006 1033 } 1007 1034 } … … 1009 1036 } 1010 1037 1011 double pmSourceDataDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)1038 double pmSourceDataDotModel (const pmSource *Mi, const pmSource *Mj, const pmSourceFitVarMode fitVarMode, const float covarFactor, psImageMaskType maskVal) 1012 1039 { 1013 1040 PS_ASSERT_PTR_NON_NULL(Mi, NAN); … … 1017 1044 int xIs, xJs, yIs, yJs; 1018 1045 int xIe, yIe; 1019 double flux, wt; 1046 double flux; 1047 double wt = 1.0; 1020 1048 1021 1049 const psImage *Pi = Mi->pixels; … … 1024 1052 assert (Pj != NULL); 1025 1053 1026 const psImage *Wi = Mi->variance; 1027 if (!unweighted_sum) { 1028 assert (Wi != NULL); 1029 } 1054 const psImage *Wi = NULL; 1055 switch (fitVarMode) { 1056 case PM_SOURCE_PHOTFIT_CONST: 1057 break; 1058 case PM_SOURCE_PHOTFIT_IMAGE_VAR: 1059 Wi = Mi->variance; 1060 psAssert (Wi, "programming error"); 1061 break; 1062 case PM_SOURCE_PHOTFIT_MODEL_VAR: 1063 Wi = Mi->modelVar; 1064 psAssert (Wi, "programming error"); 1065 break; 1066 case PM_SOURCE_PHOTFIT_NONE: 1067 psAbort("programming error"); 1068 } 1030 1069 1031 1070 const psImage *Ti = Mi->maskObj; … … 1057 1096 continue; 1058 1097 1059 // XXX skip the nonsense weight pixels? 1060 if (unweighted_sum) { 1061 flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]); 1062 } else { 1098 float value = (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]); 1099 switch (fitVarMode) { 1100 case PM_SOURCE_PHOTFIT_CONST: 1101 wt = 1.0; 1102 break; 1103 case PM_SOURCE_PHOTFIT_IMAGE_VAR: 1104 case PM_SOURCE_PHOTFIT_MODEL_VAR: 1063 1105 wt = covarFactor * Wi->data.F32[yi][xi]; 1064 if (wt > 0) { 1065 flux += (Pi->data.F32[yi][xi] * Pj->data.F32[yj][xj]) / wt; 1066 } 1067 } 1106 break; 1107 case PM_SOURCE_PHOTFIT_NONE: 1108 psAbort("programming error"); 1109 } 1110 // skip pixels with nonsense weight values 1111 if (wt <= 0) continue; 1112 1113 flux += value / wt; 1114 1068 1115 } 1069 1116 }
Note:
See TracChangeset
for help on using the changeset viewer.
