dustem_check_data.pro 20.7 KB
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417
PRO dustem_check_data,data_sed,data_ext,sed,ext,polext,polsed,polfrac,psi_em,qsed,used,qext,uext,psi_ext,fpolext

;WRITE HELP SECTION . URGENT
;NOTA BENE: LargeP has not been debiasied yet.       

;+
; NAME:
;    dustem_set_data
; PURPOSE:
;    ## adapting the data structures to the format used by the fitting algorithm 
;    Initializes the dustem data structure (!dustem_data) used in the fitting
;    Initializes the !dustem_show structure which is either a copy or a more detailed version of !dustem_data 
; CATEGORY:
;    Dustem
; CALLING SEQUENCE:
;    dustem_check_data,data_sed,data_ext,sed,ext,polext,polsed,polfrac,psi_em,qsed,used,qext,uext,psi_ext,fpolext,[/help]
; INPUTS:
;    data_sed : Emission dustemwrap structure
;    data_ext : Extinction dustemwrap structure
; OPTIONAL INPUT PARAMETERS:
;    None
; OUTPUTS:
;    ## adapting the data structures to the format used by the fitting algorithm 
;    sed : emission formatted structure 
;    ext : extinction formatted structure
;    polext : polarized extinction formatted structure
;    polsed : polarized emission formatted structure 
;    polfrac : fractional polarized emission formatted structure
;    psi_em  : emission polarization angle formatted structure 
;    qsed : stokesQ emission formatted structure
;    used : stokesU emission formatted structure
;    qext : stokesQ extinction formatted structure
;    uext : stokesU extinction formatted structure
;    psi_ext : extinction polarization angle formatted structure
;    fpolext : fractional extinction polarization formatted structure
; OPTIONAL OUTPUT PARAMETERS:
;    None
; ACCEPTED KEY-WORDS:
;    help      = If set print this help
; COMMON BLOCKS:
;    None
; SIDE EFFECTS:
;    Structures 
; RESTRICTIONS:
;    The dustem idl wrapper must be installed
; PROCEDURE:
;    None
; EXAMPLES
;    dustem_check_data,data_sed,data_ext,sed,ext,polext,polsed,polfrac,psi_em,qsed,used,qext,uext,psi_ext,fpolext,[/help]
;-


IF keyword_set(help) THEN BEGIN
  doc_library,'dustem_check_data'
  goto,the_end
ENDIF

;CHECKING THE FORMAT OF THE SED STRUCTURE
IF KEYWORD_SET(data_sed) THEN BEGIN
            
    ans=''
    wavs=data_sed.wave
    ;=== Impose central wavelengths for photometric channels
    ind=where(data_sed.filter NE 'SPECTRUM',count)
    IF count NE 0 THEN BEGIN
        wavs(ind)=dustem_filter2wav(data_sed(ind).filter)
    ENDIF
    ;=== Initializing observation structure
    obs_EM={instru_names:data_sed.instru,filt_names:data_sed.filter,wav:wavs, values:data_sed.StokesI,sigma:sqrt(abs(data_sed.sigmaII))} 
    ;INITIALIZING OUTPUT STRUCTURES
    sed = obs_EM
    qsed = sed
    used = sed
    
    ;Test if there is StokesI data. If not move to polsed testing.
    ind = where(data_sed.StokesI NE la_undef() and finite(data_sed.StokesI) EQ 1 and data_sed.sigmaII NE la_undef() and finite(data_sed.sigmaII) EQ 1,N_sed) ;Testing if the data and errors are set
    IF N_sed EQ 0 THEN BEGIN
        sed = ptr_new()
        goto, qzone
    ENDIF ELSE BEGIN 
    ind = where(finite(data_sed.StokesI) EQ 0 or finite(data_sed.sigmaII) EQ 0, countnoki) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnoki NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in StokesI data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        read, ans, prompt='Would you like automatic format filtering? (Y/N):'
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF 
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerri=where((data_sed.sigmaII EQ la_undef() and data_sed.StokesI NE la_undef()) or (data_sed.sigmaII NE la_undef() and data_sed.StokesI EQ la_undef()),ctdaerri)
    IF ctdaerri NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testi=where(data_sed.sigmaII NE 0. and finite(data_sed.sigmaII) EQ 1 and data_sed.sigmaII NE la_undef() and data_sed.StokesI NE la_undef() and finite(data_sed.StokesI) EQ 1,ctesti)                
    IF ctesti NE 0 THEN BEGIN 
        instru_names = sed.instru_names(testi)
        filt_names = sed.filt_names(testi)
        wav = sed.wav(testi)
        values = sed.values(testi)
        sigma = sed.sigma(testi)
        sed={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF    
    ENDELSE
    qzone:
    ans=''
    qsed.values=data_sed.StokesQ
    qsed.sigma=sqrt(abs(data_sed.sigmaQQ))
    ;Testig if there is StokesQ data. If not move to used data
    ind = where(data_sed.StokesQ NE la_undef() and finite(data_sed.StokesQ) EQ 1 and data_sed.sigmaQQ NE la_undef() and finite(data_sed.sigmaQQ) EQ 1,NQ)
    IF NQ EQ 0 THEN BEGIN
        qsed = ptr_new()
        goto, uzone
    ENDIF ELSE BEGIN
    ind = where(finite(data_sed.StokesQ) EQ 0 or finite(data_sed.sigmaQQ) EQ 0, countnokq) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnokq NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in StokesQ data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        read, ans, prompt='Would you like automatic format filtering? (Y/N):' 
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerrq=where((data_sed.sigmaQQ EQ la_undef() and data_sed.StokesQ NE la_undef()) or (data_sed.sigmaQQ NE la_undef() and data_sed.StokesQ EQ la_undef()),ctdaerrq)
    IF ctdaerrq NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testq=where(data_sed.sigmaQQ NE 0. and finite(data_sed.sigmaQQ) EQ 1 and data_sed.sigmaQQ NE la_undef() and data_sed.StokesQ NE la_undef() and finite(data_sed.StokesQ) EQ 1,ctestq)                
    IF ctestq NE 0 THEN BEGIN 
        instru_names = qsed.instru_names(testq)
        filt_names = qsed.filt_names(testq)
        wav = qsed.wav(testq)
        values = qsed.values(testq)
        sigma = qsed.sigma(testq)
        qsed={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF
    ENDELSE
    uzone:
    ans=''
    used.values=data_sed.StokesU
    used.sigma=sqrt(abs(data_sed.sigmaUU))
    ;Testig if there is StokesU data. If not return empty value and end loop.
    ind = where(data_sed.StokesU NE la_undef() and finite(data_sed.StokesU) EQ 1 and data_sed.sigmaUU NE la_undef() and finite(data_sed.sigmaUU) EQ 1,NU)
    IF NU EQ 0 THEN BEGIN
        used = ptr_new()
        goto, psiemzone
    ENDIF ELSE BEGIN
    ind = where(finite(data_sed.StokesU) EQ 0 or finite(data_sed.sigmaUU) EQ 0, countnoku) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnoku NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in StokesU data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        read, ans, prompt='Would you like automatic format filtering? (Y/N):' 
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF  
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerru=where((data_sed.sigmaUU EQ la_undef() and data_sed.StokesU NE la_undef()) or (data_sed.sigmaUU NE la_undef() and data_sed.StokesU EQ la_undef()),ctdaerru)
    IF ctdaerru NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testu=where(data_sed.sigmaUU NE 0. and finite(data_sed.sigmaUU) EQ 1 and data_sed.sigmaUU NE la_undef() and data_sed.StokesU NE la_undef() and finite(data_sed.StokesU) EQ 1,ctestu)                
    IF ctestu NE 0 THEN BEGIN 
        instru_names = used.instru_names(testu)
        filt_names = used.filt_names(testu)
        wav = used.wav(testu)
        values = used.values(testu)
        sigma = used.sigma(testu)
        used={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF
    ENDELSE
    psiemzone:
    
    if nu eq 0 and nq eq 0 then goto, end_sed
    degtorad = !pi/180
    
    ;Generation of POLFRAC,POLSED and PSI data by default   
    ;Testing on the presence of covariances otherwise an error is returned
    ;The testing happens here because covariances are needed for the error estimation of P,p and PSI.   
     
    tstcovu = where((data_sed(testu).sigmaIU) EQ la_undef() or (data_sed(testu).sigmaQU) EQ la_undef(), ctcovu)
    if ctcovu NE 0  then message, 'Covariance data is missing in StokesU data. Data, errors and covariances do not match. Aborting...'
    tstcovq = where((data_sed(testq).sigmaIQ) EQ la_undef() or (data_sed(testq).sigmaQU) EQ la_undef(), ctcovq)
    if ctcovq NE 0  then message, 'Covariance data is missing in StokesQ data. Data, errors and covariances do not match. Aborting...'
    
    ;Matching both filter and spectrum points 
    
    match2,used.wav,qsed.wav,qpts,upts
    
    ;removing -1 values from the indices (qpts & upts)
    indnotq = where(qpts eq -1, cnotq)
    if cnotq ne 0 then remove, indnotq, qpts
    indnotu = where(upts eq -1, cnotu)
    if cnotu ne 0 then remove, indnotu, upts
    
    match2, sed.wav, qsed.wav(qpts), qtmp, sedpts
    indnoti = where(sedpts eq -1, cnoti)
    if cnoti ne 0 then remove, indnoti, sedpts
    
    instru_names = used.instru_names(upts)
    filt_names = used.filt_names(upts)
    wav = used.wav(upts)
    values = 0.5*atan(used.values(upts),qsed.values(qpts))/degtorad
    polar_variance_iqu2ippsi,sed.values(sedpts),qsed.values(qpts),used.values(upts),(data_sed(testi).sigmaII)(sedpts),(data_sed(testq).sigmaQQ)(qpts),(data_sed(testu).sigmaUU)(upts),(data_sed(testi).sigmaIQ)(sedpts),(data_sed(testi).sigmaIU)(sedpts),(data_sed(testi).sigmaQU)(sedpts),variance_smallp,variance_psi,variance_largep
    sigma = sqrt(variance_psi)
    psi_em={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    
    goto, pzone ; the variances for polfrac and polsed were already computed. 
        
    pzone:
    
    polsed = psi_em
    polsed.values = sqrt(used.values(upts)^2+qsed.values(qpts)^2) 
    polsed.sigma = sqrt(variance_largep)
    goto, fpzone
    
    fpzone:
    
    polfrac = psi_em
    polfrac.values = polsed.values/sed.values(sedpts) 
    polfrac.sigma = sqrt(variance_smallp)
    
    goto, end_sed
    end_sed:
    
ENDIF


;EXTINCTION HASN'T BEEN TREATED YET 


;SMALLL ISSUE .. POLEXTFRAC is not present. Let's create it for completness

IF KEYWORD_SET(data_ext) THEN begin
    
    ans=''
    wavs=data_ext.wave
    ;=== Impose central wavelengths for photometric channels
    ind=where(data_ext.filter NE 'SPECTRUM',count)
    IF count NE 0 THEN BEGIN
        wavs(ind)=dustem_filter2wav(data_ext(ind).filter)
    ENDIF
    ;=== define observations structure 
    obs_EXT={instru_names:data_ext.instru,filt_names:data_ext.filter,wav:wavs, values:data_ext.EXT_I,sigma:sqrt(abs(data_ext.sigextII))} 
    ;INITIALIZING OUTPUT STRUCTURES
    ext = obs_EXT
    qext = ext
    uext = ext
    
    ind = where(data_ext.EXT_I NE la_undef() and finite(data_ext.EXT_I) EQ 1 and data_ext.sigextII NE la_undef() and finite(data_ext.sigextII) EQ 1,N_ext) ;Testing if the data and errors are set
    IF N_ext EQ 0 THEN BEGIN
        ext = ptr_new()
        goto, qexzone
    ENDIF ELSE BEGIN
    ind = where(finite(data_ext.EXT_I) EQ 0 or finite(data_ext.sigextII) EQ 0, countnokex) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnokex NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in EXT_I data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        ans=''
        read, ans, prompt='Would you like automatic format filtering? (Y/N):' 
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF      
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $ 
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerrex=where((data_ext.sigextII EQ la_undef() and data_ext.EXT_I NE la_undef()) or (data_ext.sigextII NE la_undef() and data_ext.EXT_I EQ la_undef()),ctdaerrex)
    IF ctdaerrex NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testex=where(data_ext.sigextII NE 0. and finite(data_ext.sigextII) EQ 1 and data_ext.sigextII NE la_undef() and data_ext.EXT_I NE la_undef() and finite(data_ext.EXT_I) EQ 1,ctestex)                
    IF ctestex NE 0 THEN BEGIN 
        instru_names = ext.instru_names(testex)
        filt_names = ext.filt_names(testex)
        wav = ext.wav(testex)
        values = ext.values(testex)
        sigma = ext.sigma(testex)
        ext={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF 
    ENDELSE
    
    qexzone:
    ans=''
    qext.values=data_ext.EXT_Q
    qext.sigma=sqrt(abs(data_ext.sigextQQ))
    ;Testig if there is EXT_P data. If not move to EXT_U testing
    ind = where(data_ext.EXT_Q NE la_undef() and finite(data_ext.EXT_Q) EQ 1 and data_ext.sigextQQ NE la_undef() and finite(data_ext.sigextQQ) EQ 1,NQ_ext) ;Testing if the data and errors are set
    ;dustem_set_data.prostop   
    IF NQ_ext EQ 0 THEN BEGIN
        qext = ptr_new()
        goto, uexzone
    ENDIF ELSE BEGIN
    ind = where(finite(data_ext.EXT_Q) EQ 0 or finite(data_ext.sigextQQ) EQ 0, countnokqex) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnokqex NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in EXT_Q data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        ans=''
        read, ans, prompt='Would you like automatic format filtering? (Y/N):' 
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerrqex=where((data_ext.sigextQQ EQ la_undef() and data_ext.EXT_Q NE la_undef()) or (data_ext.sigextQQ NE la_undef() and data_ext.EXT_Q EQ la_undef()),ctdaerrqex)
    IF ctdaerrqex NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testqex=where(data_ext.sigextQQ NE 0. and finite(data_ext.sigextQQ) EQ 1 and data_ext.sigextQQ NE la_undef() and data_ext.EXT_Q NE la_undef() and finite(data_ext.EXT_Q) EQ 1,ctestqex)                
    IF ctestqex NE 0 THEN BEGIN 
        instru_names = qext.instru_names(testqex)
        filt_names = qext.filt_names(testqex)
        wav = qext.wav(testqex)
        values = qext.values(testqex)
        sigma = qext.sigma(testqex)
        qext={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF 
    ENDELSE
    
    uexzone:
    ans=''
    uext.values=data_ext.EXT_U
    uext.sigma=sqrt(abs(data_ext.sigextUU))
    ;Testig if there is EXT_P data. If not return empty value and end loop.
    ind = where(data_ext.EXT_U NE la_undef() and finite(data_ext.EXT_U) EQ 1 and data_ext.sigextUU NE la_undef() and finite(data_ext.sigextUU) EQ 1,NU_ext) ;Testing if the data and errors are set
    IF NU_ext EQ 0 THEN BEGIN
        uext = ptr_new()
        goto, psiexzone
    ENDIF ELSE BEGIN
    ind = where(finite(data_ext.EXT_U) EQ 0 or finite(data_ext.sigextUU) EQ 0, countnokuex) ;Format is not okay if there are nan values but the wrapper recognizes the la_undef() value 
    IF countnokuex NE 0 THEN BEGIN 
        message, "SED format isn't met. NaN value(s) detected in EXT_U data.",/continue
        message, 'Without automatic format fitering you cannot proceed.',/continue
        ans=''
        read, ans, prompt='Would you like automatic format filtering? (Y/N):' 
        IF strupcase(ans) EQ 'N' THEN BEGIN
            message,'You have chosen not to use automatic filtering. Please adjust the SED format.',/continue
            message, 'Aborting ...'
        ENDIF
    ENDIF
    IF strupcase(ans) EQ 'Y' THEN message, 'You have chosen automatic format filtering. Removal of undefined data points...',/continue $
    ELSE message, 'Data format is met. Removal of undefined data points...',/continue
    tdaerruex=where((data_ext.sigextUU EQ la_undef() and data_ext.EXT_U NE la_undef()) or (data_ext.sigextUU NE la_undef() and data_ext.EXT_U EQ la_undef()),ctdaerruex)
    IF ctdaerruex NE 0 THEN message, 'Data and errors do not match. Undefinded data detected. Aborting...'
    testuex=where(data_ext.sigextUU NE 0. and finite(data_ext.sigextUU) EQ 1 and data_ext.sigextUU NE la_undef() and data_ext.EXT_U NE la_undef() and finite(data_ext.EXT_U) EQ 1,ctestuex)                
    IF ctestuex NE 0 THEN BEGIN 
        instru_names = uext.instru_names(testuex)
        filt_names = uext.filt_names(testuex)
        wav = uext.wav(testuex)
        values = uext.values(testuex)
        sigma = uext.sigma(testuex)
        uext={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    ENDIF 
    ENDELSE
    
    psiexzone: ;Starting the analogy here. 
    
    if NU_EXT eq 0 and NQ_EXT eq 0 then goto, end_ext
    degtorad = !pi/180
    
    ;Generation of POLEXT and PSI_EXT data by default   
    ;Testing on the presence of covariances otherwise an error is returned
    ;The testing happens here because covariances are needed for the error estimation of EXT_P,PSI_EXT.   
    
    tstcovuex = where((data_ext(testuex).sigextIU) EQ la_undef() or (data_ext(testuex).sigextQU) EQ la_undef(), ctcovuex)
    if ctcovuex NE 0  then message, 'Covariance data is missing in the extinction data pertaining to the StokesU parameter. Data, errors and covariances do not match. Aborting...'
    tstcovqex = where((data_ext(testqex).sigextIQ) EQ la_undef() or (data_ext(testqex).sigextQU) EQ la_undef(), ctcovqex)
    if ctcovqex NE 0  then message, 'Covariance data is missing in the extinction pertaining to the StokesQ parameter. Data, errors and covariances do not match. Aborting...'
    
    ;Matching both filter and spectrum points 
    match2,uext.wav,qext.wav,qexpts,uexpts
    
    ;removing -1 values from the indices (qpts & upts)
    indnotqex = where(qexpts eq -1, cnotqex)
    if cnotqex ne 0 then remove, indnotqex, qexpts
    indnotuex = where(uexpts eq -1, cnotuex)
    if cnotuex ne 0 then remove, indnotuex, uexpts
    
    match2, ext.wav, qext.wav(qexpts), qextmp, extpts
    indnotx = where(extpts eq -1, cnotx)
    if cnotx ne 0 then remove, indnotx, extpts
    
    instru_names = uext.instru_names(uexpts)
    filt_names = uext.filt_names(uexpts)
    wav = uext.wav(uexpts)
    values = 0.5*atan(uext.values(uexpts),qext.values(qexpts))/degtorad
    polar_variance_iqu2ippsi,ext.values(extpts),qext.values(qexpts),uext.values(uexpts),(data_ext(testex).sigextII)(extpts),(data_ext(testqex).sigextQQ)(qexpts),(data_ext(testuex).sigextUU)(uexpts),(data_ext(testex).sigextIQ)(extpts),(data_ext(testex).sigextIU)(extpts),(data_ext(testex).sigextQU)(extpts),variance_smallp,variance_psi,variance_largep
    sigma = sqrt(variance_psi)
    psi_ext={instru_names:instru_names,filt_names:filt_names,wav:wav,values:values,sigma:sigma}
    
    goto, pexzone  ; the variances for polextPSI_EXT were already computed. 
        
    pexzone:
    
    polext = psi_ext
    polext.values = sqrt(uext.values(uexpts)^2+qext.values(qexpts)^2) 
    polext.sigma = sqrt(variance_largep)
    goto, fpexzone
    
    fpexzone:
    
    fpolext = psi_ext
    fpolext.values = polext.values/ext.values(extpts) 
    fpolext.sigma = sqrt(variance_smallp)
    
    goto, end_ext
    
    
    end_ext:
ENDIF

the_end:

;stop
; RETURN, 0.


END