Roll back IQU=0 for pol > 1, test multiple smoothings with nan handling

src/FOC_reduction.py

@@ -102,11 +102,11 @@ def main():
     smoothing_FWHM = 1.00           #If None, no smoothing is done
     smoothing_scale = 'arcsec'       #pixel or arcsec
     # Rotation
-    rotate_stokes = False           #rotation to North convention can give erroneous results
+    rotate_stokes = True           #rotation to North convention can give erroneous results
     rotate_data = False              #rotation to North convention can give erroneous results
     # Polarization map output
     figname = '3C405_FOC'         #target/intrument name
-    figtype = '_combine_FWHM100'    #additionnal informations
+    figtype = '_combine_FWHM100_rot'    #additionnal informations
     SNRp_cut = 3    #P measurments with SNR>3
     SNRi_cut = 30   #I measurments with SNR>30, which implies an uncertainty in P of 4.7%.
     step_vec = 1    #plot all vectors in the array. if step_vec = 2, then every other vector will be plotted

src/lib/reduction.py

@@ -1018,8 +1018,8 @@ def compute_Stokes(data_array, error_array, headers, FWHM=None,
             plt.show()
 
             #I_stokes[mask]=0.
-            Q_stokes[mask]=0.
+            #Q_stokes[mask]=0.
-            U_stokes[mask]=0.
+            #U_stokes[mask]=0.
 
         #Stokes covariance matrix
         Stokes_cov = np.zeros((3,3,I_stokes.shape[0],I_stokes.shape[1]))
@@ -1427,4 +1427,20 @@ def rotate2_Stokes(I_stokes, Q_stokes, U_stokes, Stokes_cov, headers, ang):
 
         new_headers.append(new_header)
 
+    # Nan handling :
+    fmax = np.finfo(np.float64).max
+
+    new_I_stokes[np.isnan(new_I_stokes)] = 0.
+    new_Q_stokes[new_I_stokes == 0.] = 0.
+    new_U_stokes[new_I_stokes == 0.] = 0.
+    new_Q_stokes[np.isnan(new_Q_stokes)] = 0.
+    new_U_stokes[np.isnan(new_U_stokes)] = 0.
+    new_Stokes_cov[np.isnan(new_Stokes_cov)] = fmax
+    P[np.isnan(P)] = 0.
+    s_P[np.isnan(s_P)] = fmax
+    s_PA[np.isnan(s_PA)] = fmax
+    debiased_P[np.isnan(debiased_P)] = 0.
+    s_P_P[np.isnan(s_P_P)] = fmax
+    s_PA_P[np.isnan(s_PA_P)] = fmax
+
     return new_I_stokes, new_Q_stokes, new_U_stokes, new_Stokes_cov, P, debiased_P, s_P, s_P_P, PA, s_PA, s_PA_P, new_headers