Mercurial > public > atmospheric_lidar_ghk / file comparison

comparison: lidar_correction_ghk.py

lidar_correction_ghk.py

changeset 23
ef8a64173c96
parent 21
857c95060313
child 25
4c66b9ca23be
equal deleted inserted replaced
22:74ea89fcbeb2 23:ef8a64173c96
1 # -*- coding: utf-8 -*- 1 # -*- coding: utf-8 -*-
2 """ 2 """
3 Copyright 2016 Volker Freudenthaler 3 Copyright 2016, 2017 Volker Freudenthaler
4 4
5 Licensed under the EUPL, Version 1.1 only (the "Licence"). 5 Licensed under the EUPL, Version 1.1 only (the "Licence").
6 6
7 You may not use this work except in compliance with the Licence. 7 You may not use this work except in compliance with the Licence.
8 A copy of the licence is distributed with the code. Alternatively, you may obtain 8 A copy of the licence is distributed with the code. Alternatively, you may obtain
115 # Do you want to calculate the errors? If not, just the GHK-parameters are determined. 115 # Do you want to calculate the errors? If not, just the GHK-parameters are determined.
116 Error_Calc = True 116 Error_Calc = True
117 LID = "internal" 117 LID = "internal"
118 EID = "internal" 118 EID = "internal"
119 # --- IL Laser IL and +-Uncertainty 119 # --- IL Laser IL and +-Uncertainty
120 bL = 1. # degree of linear polarization; default 1 120 DOLP, dDOLP, nDOLP = 0.995, 0.005, 1 # degree of linear polarization; default 1
121 RotL, dRotL, nRotL = 0.0, 0.0, 1 # alpha; rotation of laser polarization in degrees; default 0 121 RotL, dRotL, nRotL = 0.0, 0.0, 1 # alpha; rotation of laser polarization in degrees; default 0
122 # --- ME Emitter and +-Uncertainty 122 # --- ME Emitter and +-Uncertainty
123 DiE, dDiE, nDiE = 0., 0.00, 1 # Diattenuation 123 DiE, dDiE, nDiE = 0., 0.00, 1 # Diattenuation
124 TiE = 1. # Unpolarized transmittance 124 TiE = 1. # Unpolarized transmittance
125 RetE, dRetE, nRetE = 0., 180.0, 0 # Retardance in degrees 125 RetE, dRetE, nRetE = 0., 180.0, 0 # Retardance in degrees
202 # end of initial definition of variables 202 # end of initial definition of variables
203 # ******************************************************************************************************************************* 203 # *******************************************************************************************************************************
204 204
205 # --- Read actual lidar system parameters from optic_input.py (should be in sub-directory 'system_settings') 205 # --- Read actual lidar system parameters from optic_input.py (should be in sub-directory 'system_settings')
206 InputFile = 'optic_input_example_lidar.py' 206 InputFile = 'optic_input_example_lidar.py'
207 InputFile = 'optic_input_ver8c_POLIS_532_Sep2016.py'
208 InputFile = 'optic_input_Bertha_Saltrace_532.py'
207 209
208 ''' 210 '''
209 print("From ", dname) 211 print("From ", dname)
210 print("Running ", fname) 212 print("Running ", fname)
211 print("Reading input file ", InputFile, " for") 213 print("Reading input file ", InputFile, " for")
224 # Y = -1 226 # Y = -1
225 # LocC = 4 #location of calibrator: 1 = behind laser; 2 = behind emitter; 3 = before receiver; 4 = before PBS 227 # LocC = 4 #location of calibrator: 1 = behind laser; 2 = behind emitter; 3 = before receiver; 4 = before PBS
226 ##TypeC = 6 Don't change the TypeC here 228 ##TypeC = 6 Don't change the TypeC here
227 # RotationErrorEpsilonForNormalMeasurements = True 229 # RotationErrorEpsilonForNormalMeasurements = True
228 # LDRCal = 0.25 230 # LDRCal = 0.25
229 # bL = 0.8
230 ## --- Errors 231 ## --- Errors
232 DOLP0, dDOLP, nDOLP = DOLP, dDOLP, nDOLP
231 RotL0, dRotL, nRotL = RotL, dRotL, nRotL 233 RotL0, dRotL, nRotL = RotL, dRotL, nRotL
232 234
233 DiE0, dDiE, nDiE = DiE, dDiE, nDiE 235 DiE0, dDiE, nDiE = DiE, dDiE, nDiE
234 RetE0, dRetE, nRetE = RetE, dRetE, nRetE 236 RetE0, dRetE, nRetE = RetE, dRetE, nRetE
235 RotE0, dRotE, nRotE = RotE, dRotE, nRotE 237 RotE0, dRotE, nRotE = RotE, dRotE, nRotE
272 DiR = (RP - RS) / (RP + RS) 274 DiR = (RP - RS) / (RP + RS)
273 ZiR = (1. - DiR ** 2) ** 0.5 275 ZiR = (1. - DiR ** 2) ** 0.5
274 276
275 277
276 # --- this subroutine is for the calculation with certain fixed parameters ----------------------------------------------------- 278 # --- this subroutine is for the calculation with certain fixed parameters -----------------------------------------------------
277 def Calc(RotL, RotE, RetE, DiE, RotO, RetO, DiO, RotC, RetC, DiC, TP, TS, RP, RS, ERaT, RotaT, RetT, ERaR, RotaR, RetR, 279 def Calc(DOLP, RotL, RotE, RetE, DiE, RotO, RetO, DiO, RotC, RetC, DiC, TP, TS, RP, RS, ERaT, RotaT, RetT, ERaR, RotaR, RetR,
278 LDRCal): 280 LDRCal):
279 # ---- Do the calculations of bra-ket vectors 281 # ---- Do the calculations of bra-ket vectors
280 h = -1. if TypeC == 2 else 1 282 h = -1. if TypeC == 2 else 1
281 # from input file: assumed LDRCal for calibration measurements 283 # from input file: assumed LDRCal for calibration measurements
282 aCal = (1. - LDRCal) / (1 + LDRCal) 284 aCal = (1. - LDRCal) / (1 + LDRCal)
294 S2ab = np.sin(np.deg2rad(2 * RotL - 2 * RotE)) 296 S2ab = np.sin(np.deg2rad(2 * RotL - 2 * RotE))
295 C2ab = np.cos(np.deg2rad(2 * RotL - 2 * RotE)) 297 C2ab = np.cos(np.deg2rad(2 * RotL - 2 * RotE))
296 S2g = np.sin(np.deg2rad(2 * RotO)) 298 S2g = np.sin(np.deg2rad(2 * RotO))
297 C2g = np.cos(np.deg2rad(2 * RotO)) 299 C2g = np.cos(np.deg2rad(2 * RotO))
298 300
299 # Laser with Degree of linear polarization DOLP = bL 301 # Laser with Degree of linear polarization DOLP
300 IinL = 1. 302 IinL = 1.
301 QinL = bL 303 QinL = DOLP
302 UinL = 0. 304 UinL = 0.
303 VinL = (1. - bL ** 2) ** 0.5 305 VinL = (1. - DOLP ** 2) ** 0.5
304 306
305 # Stokes Input Vector rotation Eq. E.4 307 # Stokes Input Vector rotation Eq. E.4
306 A = C2a * QinL - S2a * UinL 308 A = C2a * QinL - S2a * UinL
307 B = S2a * QinL + C2a * UinL 309 B = S2a * QinL + C2a * UinL
308 # Stokes Input Vector rotation Eq. E.9 310 # Stokes Input Vector rotation Eq. E.9
378 S2aT = np.sin(np.deg2rad(h * 2 * RotaT)) 380 S2aT = np.sin(np.deg2rad(h * 2 * RotaT))
379 C2aT = np.cos(np.deg2rad(2 * RotaT)) 381 C2aT = np.cos(np.deg2rad(2 * RotaT))
380 S2aR = np.sin(np.deg2rad(h * 2 * RotaR)) 382 S2aR = np.sin(np.deg2rad(h * 2 * RotaR))
381 C2aR = np.cos(np.deg2rad(2 * RotaR)) 383 C2aR = np.cos(np.deg2rad(2 * RotaR))
382 384
383 # Aanalyzer As before the PBS Eq. D.5 385 # Analyzer As before the PBS Eq. D.5
384 ATP1 = (1 + C2aT * DaT * DiT) 386 ATP1 = (1 + C2aT * DaT * DiT)
385 ATP2 = Y * (DiT + C2aT * DaT) 387 ATP2 = Y * (DiT + C2aT * DaT)
386 ATP3 = Y * S2aT * DaT * ZiT * CosT 388 ATP3 = Y * S2aT * DaT * ZiT * CosT
387 ATP4 = S2aT * DaT * ZiT * SinT 389 ATP4 = S2aT * DaT * ZiT * SinT
388 ATP = np.array([ATP1, ATP2, ATP3, ATP4]) 390 ATP = np.array([ATP1, ATP2, ATP3, ATP4])
451 # parameters for calibration with aCal 453 # parameters for calibration with aCal
452 AT = ATP1 * IinP + h * ATP4 * VinP 454 AT = ATP1 * IinP + h * ATP4 * VinP
453 BT = ATP3e * QinP - h * ATP2e * UinP 455 BT = ATP3e * QinP - h * ATP2e * UinP
454 AR = ARP1 * IinP + h * ARP4 * VinP 456 AR = ARP1 * IinP + h * ARP4 * VinP
455 BR = ARP3e * QinP - h * ARP2e * UinP 457 BR = ARP3e * QinP - h * ARP2e * UinP
456 # Correction paremeters for normal measurements; they are independent of LDR 458 # Correction parameters for normal measurements; they are independent of LDR
457 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 459 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
458 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 460 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
459 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 461 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
460 GT = np.dot(ATP, IS1) 462 GT = np.dot(ATP, IS1)
461 GR = np.dot(ARP, IS1) 463 GR = np.dot(ARP, IS1)
472 # parameters for calibration with aCal 474 # parameters for calibration with aCal
473 AT = ATP1 * IinP + ATP3e * UinPe + ZiC * CosC * (ATP2e * QinPe + ATP4 * VinP) 475 AT = ATP1 * IinP + ATP3e * UinPe + ZiC * CosC * (ATP2e * QinPe + ATP4 * VinP)
474 BT = DiC * (ATP1 * UinPe + ATP3e * IinP) - ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe) 476 BT = DiC * (ATP1 * UinPe + ATP3e * IinP) - ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe)
475 AR = ARP1 * IinP + ARP3e * UinPe + ZiC * CosC * (ARP2e * QinPe + ARP4 * VinP) 477 AR = ARP1 * IinP + ARP3e * UinPe + ZiC * CosC * (ARP2e * QinPe + ARP4 * VinP)
476 BR = DiC * (ARP1 * UinPe + ARP3e * IinP) - ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe) 478 BR = DiC * (ARP1 * UinPe + ARP3e * IinP) - ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe)
477 # Correction paremeters for normal measurements; they are independent of LDR 479 # Correction parameters for normal measurements; they are independent of LDR
478 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 480 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
479 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 481 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
480 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 482 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
481 GT = np.dot(ATP, IS1) 483 GT = np.dot(ATP, IS1)
482 GR = np.dot(ARP, IS1) 484 GR = np.dot(ARP, IS1)
499 ATP2e * UinPe + ATP3e * QinPe) - sqr05 * ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe) 501 ATP2e * UinPe + ATP3e * QinPe) - sqr05 * ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe)
500 AR = ARP1 * IinP + sqr05 * DiC * (ARP1 * QinPe + ARP2e * IinP) + (1 - 0.5 * WiC) * ( 502 AR = ARP1 * IinP + sqr05 * DiC * (ARP1 * QinPe + ARP2e * IinP) + (1 - 0.5 * WiC) * (
501 ARP2e * QinPe + ARP3e * UinPe) + ZiC * (sqr05 * SinC * (ARP3e * VinP - ARP4 * UinPe) + ARP4 * CosC * VinP) 503 ARP2e * QinPe + ARP3e * UinPe) + ZiC * (sqr05 * SinC * (ARP3e * VinP - ARP4 * UinPe) + ARP4 * CosC * VinP)
502 BR = sqr05 * DiC * (ARP1 * UinPe + ARP3e * IinP) + 0.5 * WiC * ( 504 BR = sqr05 * DiC * (ARP1 * UinPe + ARP3e * IinP) + 0.5 * WiC * (
503 ARP2e * UinPe + ARP3e * QinPe) - sqr05 * ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe) 505 ARP2e * UinPe + ARP3e * QinPe) - sqr05 * ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe)
504 # Correction paremeters for normal measurements; they are independent of LDR 506 # Correction parameters for normal measurements; they are independent of LDR
505 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 507 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
506 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 508 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
507 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 509 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
508 GT = np.dot(ATP, IS1) 510 GT = np.dot(ATP, IS1)
509 GR = np.dot(ARP, IS1) 511 GR = np.dot(ARP, IS1)
587 # parameters for calibration with aCal 589 # parameters for calibration with aCal
588 AT = ATF1 * IinF + ATF4 * h * VinF 590 AT = ATF1 * IinF + ATF4 * h * VinF
589 BT = ATF3e * QinF - ATF2e * h * UinF 591 BT = ATF3e * QinF - ATF2e * h * UinF
590 AR = ARF1 * IinF + ARF4 * h * VinF 592 AR = ARF1 * IinF + ARF4 * h * VinF
591 BR = ARF3e * QinF - ARF2e * h * UinF 593 BR = ARF3e * QinF - ARF2e * h * UinF
592 # Correction paremeters for normal measurements; they are independent of LDR 594 # Correction parameters for normal measurements; they are independent of LDR
593 if (not RotationErrorEpsilonForNormalMeasurements): 595 if (not RotationErrorEpsilonForNormalMeasurements):
594 GT = ATF1 * IinE + ATF4 * VinE 596 GT = ATF1 * IinE + ATF4 * VinE
595 GR = ARF1 * IinE + ARF4 * VinE 597 GR = ARF1 * IinE + ARF4 * VinE
596 HT = ATF2 * QinE - ATF3 * UinE - ATF4 * 2 * VinE 598 HT = ATF2 * QinE - ATF3 * UinE - ATF4 * 2 * VinE
597 HR = ARF2 * QinE - ARF3 * UinE - ARF4 * 2 * VinE 599 HR = ARF2 * QinE - ARF3 * UinE - ARF4 * 2 * VinE
607 AT = np.dot(ATFe, IF1e) 609 AT = np.dot(ATFe, IF1e)
608 AR = np.dot(ARFe, IF1e) 610 AR = np.dot(ARFe, IF1e)
609 BT = np.dot(ATFe, IF2e) 611 BT = np.dot(ATFe, IF2e)
610 BR = np.dot(ARFe, IF2e) 612 BR = np.dot(ARFe, IF2e)
611 613
612 # Correction paremeters for normal measurements; they are independent of LDR --- the same as for TypeC = 6 614 # Correction parameters for normal measurements; they are independent of LDR --- the same as for TypeC = 6
613 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 615 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
614 IS1 = np.array([IinE, 0, 0, VinE]) 616 IS1 = np.array([IinE, 0, 0, VinE])
615 IS2 = np.array([0, QinE, -UinE, -2 * VinE]) 617 IS2 = np.array([0, QinE, -UinE, -2 * VinE])
616 GT = np.dot(ATF, IS1) 618 GT = np.dot(ATF, IS1)
617 GR = np.dot(ARF, IS1) 619 GR = np.dot(ARF, IS1)
637 AT = np.dot(ATFe, IF1e) 639 AT = np.dot(ATFe, IF1e)
638 AR = np.dot(ARFe, IF1e) 640 AR = np.dot(ARFe, IF1e)
639 BT = np.dot(ATFe, IF2e) 641 BT = np.dot(ATFe, IF2e)
640 BR = np.dot(ARFe, IF2e) 642 BR = np.dot(ARFe, IF2e)
641 643
642 # Correction paremeters for normal measurements; they are independent of LDR 644 # Correction parameters for normal measurements; they are independent of LDR
643 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 645 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
644 # IS1 = np.array([IinE,0,0,VinE]) 646 # IS1 = np.array([IinE,0,0,VinE])
645 # IS2 = np.array([0,QinE,-UinE,-2*VinE]) 647 # IS2 = np.array([0,QinE,-UinE,-2*VinE])
646 IS1 = np.array([IinFo, 0, 0, VinFo]) 648 IS1 = np.array([IinFo, 0, 0, VinFo])
647 IS2 = np.array([0, QinFa, UinFa, VinFa]) 649 IS2 = np.array([0, QinFa, UinFa, VinFa])
722 AT = ATE1o * IinE + (ATE4o + aCal * ATE4a) * h * VinE 724 AT = ATE1o * IinE + (ATE4o + aCal * ATE4a) * h * VinE
723 BT = aCal * (ATE3ae * QinEe - ATE2ae * h * UinEe) 725 BT = aCal * (ATE3ae * QinEe - ATE2ae * h * UinEe)
724 AR = ARE1o * IinE + (ARE4o + aCal * ARE4a) * h * VinE 726 AR = ARE1o * IinE + (ARE4o + aCal * ARE4a) * h * VinE
725 BR = aCal * (ARE3ae * QinEe - ARE2ae * h * UinEe) 727 BR = aCal * (ARE3ae * QinEe - ARE2ae * h * UinEe)
726 728
727 # Correction paremeters for normal measurements; they are independent of LDR 729 # Correction parameters for normal measurements; they are independent of LDR
728 if (not RotationErrorEpsilonForNormalMeasurements): 730 if (not RotationErrorEpsilonForNormalMeasurements):
729 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F) 731 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F)
730 GT = ATE1o * IinE + ATE4o * h * VinE 732 GT = ATE1o * IinE + ATE4o * h * VinE
731 GR = ARE1o * IinE + ARE4o * h * VinE 733 GR = ARE1o * IinE + ARE4o * h * VinE
732 HT = ATE2a * QinE + ATE3a * h * UinEe + ATE4a * h * VinE 734 HT = ATE2a * QinE + ATE3a * h * UinEe + ATE4a * h * VinE
742 AT = ATE1 * IinE + ZiC * CosC * (ATE2e * QinEe + ATE4 * VinE) + ATE3e * UinEe 744 AT = ATE1 * IinE + ZiC * CosC * (ATE2e * QinEe + ATE4 * VinE) + ATE3e * UinEe
743 BT = DiC * (ATE1 * UinEe + ATE3e * IinE) + ZiC * SinC * (ATE4 * QinEe - ATE2e * VinE) 745 BT = DiC * (ATE1 * UinEe + ATE3e * IinE) + ZiC * SinC * (ATE4 * QinEe - ATE2e * VinE)
744 AR = ARE1 * IinE + ZiC * CosC * (ARE2e * QinEe + ARE4 * VinE) + ARE3e * UinEe 746 AR = ARE1 * IinE + ZiC * CosC * (ARE2e * QinEe + ARE4 * VinE) + ARE3e * UinEe
745 BR = DiC * (ARE1 * UinEe + ARE3e * IinE) + ZiC * SinC * (ARE4 * QinEe - ARE2e * VinE) 747 BR = DiC * (ARE1 * UinEe + ARE3e * IinE) + ZiC * SinC * (ARE4 * QinEe - ARE2e * VinE)
746 748
747 # Correction paremeters for normal measurements; they are independent of LDR 749 # Correction parameters for normal measurements; they are independent of LDR
748 if (not RotationErrorEpsilonForNormalMeasurements): 750 if (not RotationErrorEpsilonForNormalMeasurements):
749 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F) 751 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F)
750 GT = ATE1o * IinE + ATE4o * VinE 752 GT = ATE1o * IinE + ATE4o * VinE
751 GR = ARE1o * IinE + ARE4o * VinE 753 GR = ARE1o * IinE + ARE4o * VinE
752 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE 754 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE
773 AT = np.dot(ATEe, IE1e) 775 AT = np.dot(ATEe, IE1e)
774 AR = np.dot(AREe, IE1e) 776 AR = np.dot(AREe, IE1e)
775 BT = np.dot(ATEe, IE2e) 777 BT = np.dot(ATEe, IE2e)
776 BR = np.dot(AREe, IE2e) 778 BR = np.dot(AREe, IE2e)
777 779
778 # Correction paremeters for normal measurements; they are independent of LDR 780 # Correction parameters for normal measurements; they are independent of LDR
779 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 781 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
780 GT = ATE1o * IinE + ATE4o * VinE 782 GT = ATE1o * IinE + ATE4o * VinE
781 GR = ARE1o * IinE + ARE4o * VinE 783 GR = ARE1o * IinE + ARE4o * VinE
782 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE 784 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE
783 HR = ARE2a * QinE + ARE3a * UinE + ARE4a * VinE 785 HR = ARE2a * QinE + ARE3a * UinE + ARE4a * VinE
846 848
847 849
848 # ******************************************************************************************************************************* 850 # *******************************************************************************************************************************
849 851
850 # --- CALC truth with assumed true parameters from the input file 852 # --- CALC truth with assumed true parameters from the input file
851 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(RotL0, RotE0, RetE0, DiE0, RotO0, 853 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(DOLP0, RotL0, RotE0, RetE0, DiE0, RotO0,
852 RetO0, DiO0, RotC0, RetC0, DiC0, 854 RetO0, DiO0, RotC0, RetC0, DiC0,
853 TP0, TS0, RP0, RS0, ERaT0, 855 TP0, TS0, RP0, RS0, ERaT0,
854 RotaT0, RetT0, ERaR0, RotaR0, 856 RotaT0, RetT0, ERaR0, RotaR0,
855 RetR0, LDRCal0) 857 RetR0, LDRCal0)
856 858
861 print("Running ", fname) 863 print("Running ", fname)
862 print("Reading input file ", InputFile) # , " for Lidar system :", EID, ", ", LID) 864 print("Reading input file ", InputFile) # , " for Lidar system :", EID, ", ", LID)
863 print("for Lidar system: ", EID, ", ", LID) 865 print("for Lidar system: ", EID, ", ", LID)
864 # --- Print iput information********************************* 866 # --- Print iput information*********************************
865 print(" --- Input parameters: value ±error / ±steps ----------------------") 867 print(" --- Input parameters: value ±error / ±steps ----------------------")
866 print("{0:8} {1:8} {2:8.5f}; {3:8} {4:7.4f}±{5:7.4f}/{6:2d}".format("Laser: ", "DOLP = ", bL, 868 print("{0:5}{1:5} {2:6.4f}±{3:7.4f}/{4:2d}; {5:8} {6:8.4f}±{7:7.4f}/{8:2d}".format("Laser: ", "DOLP =", DOLP0, dDOLP, nDOLP,
867 " rotation alpha = ", RotL0, dRotL, 869 " Rotation alpha = ", RotL0, dRotL,
868 nRotL)) 870 nRotL))
869 print(" Diatt., Tunpol, Retard., Rotation (deg)") 871 print(" Diatt., Tunpol, Retard., Rotation (deg)")
870 print("{0:12} {1:7.4f}±{2:7.4f}/{8:2d}, {3:7.4f}, {4:3.0f}±{5:3.0f}/{9:2d}, {6:7.4f}±{7:7.4f}/{10:2d}".format( 872 print("{0:12} {1:7.4f}±{2:7.4f}/{8:2d}, {3:7.4f}, {4:3.0f}±{5:3.0f}/{9:2d}, {6:7.4f}±{7:7.4f}/{10:2d}".format(
871 "Emitter ", DiE0, dDiE, TiE, RetE0, dRetE, RotE0, dRotE, nDiE, nRetE, nRotE)) 873 "Emitter ", DiE0, dDiE, TiE, RetE0, dRetE, RotE0, dRotE, nDiE, nRetE, nRotE))
872 print("{0:12} {1:7.4f}±{2:7.4f}/{8:2d}, {3:7.4f}, {4:3.0f}±{5:3.0f}/{9:2d}, {6:7.4f}±{7:7.4f}/{10:2d}".format( 874 print("{0:12} {1:7.4f}±{2:7.4f}/{8:2d}, {3:7.4f}, {4:3.0f}±{5:3.0f}/{9:2d}, {6:7.4f}±{7:7.4f}/{10:2d}".format(
914 LDRsimxList = np.zeros(3) 916 LDRsimxList = np.zeros(3)
915 LDRCalList = 0.004, 0.2, 0.45 917 LDRCalList = 0.004, 0.2, 0.45
916 # The loop over LDRCalList is ony for checking whether and how much the LDR depends on the LDRCal during calibration and whether the corrections work. 918 # The loop over LDRCalList is ony for checking whether and how much the LDR depends on the LDRCal during calibration and whether the corrections work.
917 # Still with assumed true parameters in input file 919 # Still with assumed true parameters in input file
918 for i, LDRCal in enumerate(LDRCalList): 920 for i, LDRCal in enumerate(LDRCalList):
919 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(RotL0, RotE0, RetE0, 921 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(DOLP0, RotL0, RotE0, RetE0,
920 DiE0, RotO0, RetO0, 922 DiE0, RotO0, RetO0,
921 DiO0, RotC0, RetC0, 923 DiO0, RotC0, RetC0,
922 DiC0, TP0, TS0, RP0, 924 DiC0, TP0, TS0, RP0,
923 RS0, ERaT0, RotaT0, 925 RS0, ERaT0, RotaT0,
924 RetT0, ERaR0, RotaR0, 926 RetT0, ERaR0, RotaR0,
942 # The loop over LDRtrueList is ony for checking how much the uncorrected LDRsimx deviates from LDRtrue ... and whether the corrections work. 944 # The loop over LDRtrueList is ony for checking how much the uncorrected LDRsimx deviates from LDRtrue ... and whether the corrections work.
943 # LDRsimx = LDRsim = Ir / It or 1/LDRsim 945 # LDRsimx = LDRsim = Ir / It or 1/LDRsim
944 # Still with assumed true parameters in input file 946 # Still with assumed true parameters in input file
945 for i, LDRtrue in enumerate(LDRrange): 947 for i, LDRtrue in enumerate(LDRrange):
946 #for LDRtrue in LDRrange: 948 #for LDRtrue in LDRrange:
947 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(RotL0, RotE0, RetE0, 949 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(DOLP0, RotL0, RotE0, RetE0,
948 DiE0, RotO0, RetO0, 950 DiE0, RotO0, RetO0,
949 DiO0, RotC0, RetC0, 951 DiO0, RotC0, RetC0,
950 DiC0, TP0, TS0, RP0, 952 DiC0, TP0, TS0, RP0,
951 RS0, ERaT0, RotaT0, 953 RS0, ERaT0, RotaT0,
952 RetT0, ERaR0, RotaR0, 954 RetT0, ERaR0, RotaR0,
972 f.close 974 f.close
973 sys.stdout = old_target 975 sys.stdout = old_target
974 ''' 976 '''
975 if (Error_Calc): 977 if (Error_Calc):
976 # --- CALC again assumed truth with LDRCal0 and with assumed true parameters in input file to reset all 0-values 978 # --- CALC again assumed truth with LDRCal0 and with assumed true parameters in input file to reset all 0-values
977 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(RotL0, RotE0, RetE0, DiE0, 979 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim0 = Calc(DOLP0, RotL0, RotE0, RetE0, DiE0,
978 RotO0, RetO0, DiO0, RotC0, 980 RotO0, RetO0, DiO0, RotC0,
979 RetC0, DiC0, TP0, TS0, RP0, 981 RetC0, DiC0, TP0, TS0, RP0,
980 RS0, ERaT0, RotaT0, RetT0, 982 RS0, ERaT0, RotaT0, RetT0,
981 ERaR0, RotaR0, RetR0, 983 ERaR0, RotaR0, RetR0,
982 LDRCal0) 984 LDRCal0)
983 # --- Start Errors calculation with variable parameters ------------------------------------------------------------------ 985 # --- Start Errors calculation with variable parameters ------------------------------------------------------------------
984 986
985 iN = -1 987 iN = -1
986 N = ((nRotL * 2 + 1) * 988 N = ((nDOLP * 2 + 1) * (nRotL * 2 + 1) *
987 (nRotE * 2 + 1) * (nRetE * 2 + 1) * (nDiE * 2 + 1) * 989 (nRotE * 2 + 1) * (nRetE * 2 + 1) * (nDiE * 2 + 1) *
988 (nRotO * 2 + 1) * (nRetO * 2 + 1) * (nDiO * 2 + 1) * 990 (nRotO * 2 + 1) * (nRetO * 2 + 1) * (nDiO * 2 + 1) *
989 (nRotC * 2 + 1) * (nRetC * 2 + 1) * (nDiC * 2 + 1) * 991 (nRotC * 2 + 1) * (nRetC * 2 + 1) * (nDiC * 2 + 1) *
990 (nTP * 2 + 1) * (nTS * 2 + 1) * (nRP * 2 + 1) * (nRS * 2 + 1) * (nERaT * 2 + 1) * (nERaR * 2 + 1) * 992 (nTP * 2 + 1) * (nTS * 2 + 1) * (nRP * 2 + 1) * (nRS * 2 + 1) * (nERaT * 2 + 1) * (nERaR * 2 + 1) *
991 (nRotaT * 2 + 1) * (nRotaR * 2 + 1) * (nRetT * 2 + 1) * (nRetR * 2 + 1) * (nLDRCal * 2 + 1)) 993 (nRotaT * 2 + 1) * (nRotaR * 2 + 1) * (nRetT * 2 + 1) * (nRetR * 2 + 1) * (nLDRCal * 2 + 1))
1010 #LDRrange = 0.004, 0.02, 0.1, 0.3, 0.45 1012 #LDRrange = 0.004, 0.02, 0.1, 0.3, 0.45
1011 # aLDRsimx = np.zeros(N) 1013 # aLDRsimx = np.zeros(N)
1012 # aLDRsimx2 = np.zeros(N) 1014 # aLDRsimx2 = np.zeros(N)
1013 # aLDRcorr = np.zeros(N) 1015 # aLDRcorr = np.zeros(N)
1014 # aLDRcorr2 = np.zeros(N) 1016 # aLDRcorr2 = np.zeros(N)
1017 aDOLP = np.zeros(N)
1015 aERaT = np.zeros(N) 1018 aERaT = np.zeros(N)
1016 aERaR = np.zeros(N) 1019 aERaR = np.zeros(N)
1017 aRotaT = np.zeros(N) 1020 aRotaT = np.zeros(N)
1018 aRotaR = np.zeros(N) 1021 aRotaR = np.zeros(N)
1019 aRetT = np.zeros(N) 1022 aRetT = np.zeros(N)
1054 for iLDRCal in range(-nLDRCal, nLDRCal + 1): 1057 for iLDRCal in range(-nLDRCal, nLDRCal + 1):
1055 # from input file: assumed LDRCal for calibration measurements 1058 # from input file: assumed LDRCal for calibration measurements
1056 LDRCal = LDRCal0 1059 LDRCal = LDRCal0
1057 if nLDRCal > 0: LDRCal = LDRCal0 + iLDRCal * dLDRCal / nLDRCal 1060 if nLDRCal > 0: LDRCal = LDRCal0 + iLDRCal * dLDRCal / nLDRCal
1058 1061
1059 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim = Calc(RotL0, RotE0, RetE0, 1062 GT0, HT0, GR0, HR0, K0, Eta0, LDRsimx, LDRCorr, DTa0, DRa0, TTa0, TRa0, F11sim = Calc(DOLP0, RotL0, RotE0, RetE0,
1060 DiE0, RotO0, RetO0, DiO0, 1063 DiE0, RotO0, RetO0, DiO0,
1061 RotC0, RetC0, DiC0, TP0, 1064 RotC0, RetC0, DiC0, TP0,
1062 TS0, RP0, RS0, ERaT0, 1065 TS0, RP0, RS0, ERaT0,
1063 RotaT0, RetT0, ERaR0, 1066 RotaT0, RetT0, ERaR0,
1064 RotaR0, RetR0, LDRCal) 1067 RotaR0, RetR0, LDRCal)
1065 aCal = (1. - LDRCal) / (1 + LDRCal) 1068 aCal = (1. - LDRCal) / (1 + LDRCal)
1066 for iRotL, iRotE, iRetE, iDiE \ 1069 for iDOLP, iRotL, iRotE, iRetE, iDiE \
1067 in [(iRotL, iRotE, iRetE, iDiE) 1070 in [(iDOLP, iRotL, iRotE, iRetE, iDiE)
1071 for iDOLP in range(-nDOLP, nDOLP + 1)
1068 for iRotL in range(-nRotL, nRotL + 1) 1072 for iRotL in range(-nRotL, nRotL + 1)
1069 for iRotE in range(-nRotE, nRotE + 1) 1073 for iRotE in range(-nRotE, nRotE + 1)
1070 for iRetE in range(-nRetE, nRetE + 1) 1074 for iRetE in range(-nRetE, nRetE + 1)
1071 for iDiE in range(-nDiE, nDiE + 1)]: 1075 for iDiE in range(-nDiE, nDiE + 1)]:
1072 1076
1077 if nDOLP > 0: DOLP = DOLP0 + iDOLP * dDOLP / nDOLP
1073 if nRotL > 0: RotL = RotL0 + iRotL * dRotL / nRotL 1078 if nRotL > 0: RotL = RotL0 + iRotL * dRotL / nRotL
1074 if nRotE > 0: RotE = RotE0 + iRotE * dRotE / nRotE 1079 if nRotE > 0: RotE = RotE0 + iRotE * dRotE / nRotE
1075 if nRetE > 0: RetE = RetE0 + iRetE * dRetE / nRetE 1080 if nRetE > 0: RetE = RetE0 + iRetE * dRetE / nRetE
1076 if nDiE > 0: DiE = DiE0 + iDiE * dDiE / nDiE 1081 if nDiE > 0: DiE = DiE0 + iDiE * dDiE / nDiE
1077 1082
1082 S2b = np.sin(2 * np.deg2rad(RotE)) 1087 S2b = np.sin(2 * np.deg2rad(RotE))
1083 C2b = np.cos(2 * np.deg2rad(RotE)) 1088 C2b = np.cos(2 * np.deg2rad(RotE))
1084 S2ab = np.sin(np.deg2rad(2 * RotL - 2 * RotE)) 1089 S2ab = np.sin(np.deg2rad(2 * RotL - 2 * RotE))
1085 C2ab = np.cos(np.deg2rad(2 * RotL - 2 * RotE)) 1090 C2ab = np.cos(np.deg2rad(2 * RotL - 2 * RotE))
1086 1091
1087 # Laser with Degree of linear polarization DOLP = bL 1092 # Laser with Degree of linear polarization DOLP
1088 IinL = 1. 1093 IinL = 1.
1089 QinL = bL 1094 QinL = DOLP
1090 UinL = 0. 1095 UinL = 0.
1091 VinL = (1. - bL ** 2) ** 0.5 1096 VinL = (1. - DOLP ** 2) ** 0.5
1092 1097
1093 # Stokes Input Vector rotation Eq. E.4 1098 # Stokes Input Vector rotation Eq. E.4
1094 A = C2a * QinL - S2a * UinL 1099 A = C2a * QinL - S2a * UinL
1095 B = S2a * QinL + C2a * UinL 1100 B = S2a * QinL + C2a * UinL
1096 # Stokes Input Vector rotation Eq. E.9 1101 # Stokes Input Vector rotation Eq. E.9
1274 # parameters for calibration with aCal 1279 # parameters for calibration with aCal
1275 AT = ATP1 * IinP + h * ATP4 * VinP 1280 AT = ATP1 * IinP + h * ATP4 * VinP
1276 BT = ATP3e * QinP - h * ATP2e * UinP 1281 BT = ATP3e * QinP - h * ATP2e * UinP
1277 AR = ARP1 * IinP + h * ARP4 * VinP 1282 AR = ARP1 * IinP + h * ARP4 * VinP
1278 BR = ARP3e * QinP - h * ARP2e * UinP 1283 BR = ARP3e * QinP - h * ARP2e * UinP
1279 # Correction paremeters for normal measurements; they are independent of LDR 1284 # Correction parameters for normal measurements; they are independent of LDR
1280 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1285 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1281 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 1286 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
1282 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 1287 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
1283 GT = np.dot(ATP, IS1) 1288 GT = np.dot(ATP, IS1)
1284 GR = np.dot(ARP, IS1) 1289 GR = np.dot(ARP, IS1)
1295 # parameters for calibration with aCal 1300 # parameters for calibration with aCal
1296 AT = ATP1 * IinP + ATP3e * UinPe + ZiC * CosC * (ATP2e * QinPe + ATP4 * VinP) 1301 AT = ATP1 * IinP + ATP3e * UinPe + ZiC * CosC * (ATP2e * QinPe + ATP4 * VinP)
1297 BT = DiC * (ATP1 * UinPe + ATP3e * IinP) - ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe) 1302 BT = DiC * (ATP1 * UinPe + ATP3e * IinP) - ZiC * SinC * (ATP2e * VinP - ATP4 * QinPe)
1298 AR = ARP1 * IinP + ARP3e * UinPe + ZiC * CosC * (ARP2e * QinPe + ARP4 * VinP) 1303 AR = ARP1 * IinP + ARP3e * UinPe + ZiC * CosC * (ARP2e * QinPe + ARP4 * VinP)
1299 BR = DiC * (ARP1 * UinPe + ARP3e * IinP) - ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe) 1304 BR = DiC * (ARP1 * UinPe + ARP3e * IinP) - ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe)
1300 # Correction paremeters for normal measurements; they are independent of LDR 1305 # Correction parameters for normal measurements; they are independent of LDR
1301 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1306 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1302 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 1307 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
1303 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 1308 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
1304 GT = np.dot(ATP, IS1) 1309 GT = np.dot(ATP, IS1)
1305 GR = np.dot(ARP, IS1) 1310 GR = np.dot(ARP, IS1)
1326 AR = ARP1 * IinP + sqr05 * DiC * (ARP1 * QinPe + ARP2e * IinP) + (1 - 0.5 * WiC) * ( 1331 AR = ARP1 * IinP + sqr05 * DiC * (ARP1 * QinPe + ARP2e * IinP) + (1 - 0.5 * WiC) * (
1327 ARP2e * QinPe + ARP3e * UinPe) + ZiC * ( 1332 ARP2e * QinPe + ARP3e * UinPe) + ZiC * (
1328 sqr05 * SinC * (ARP3e * VinP - ARP4 * UinPe) + ARP4 * CosC * VinP) 1333 sqr05 * SinC * (ARP3e * VinP - ARP4 * UinPe) + ARP4 * CosC * VinP)
1329 BR = sqr05 * DiC * (ARP1 * UinPe + ARP3e * IinP) + 0.5 * WiC * ( 1334 BR = sqr05 * DiC * (ARP1 * UinPe + ARP3e * IinP) + 0.5 * WiC * (
1330 ARP2e * UinPe + ARP3e * QinPe) - sqr05 * ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe) 1335 ARP2e * UinPe + ARP3e * QinPe) - sqr05 * ZiC * SinC * (ARP2e * VinP - ARP4 * QinPe)
1331 # Correction paremeters for normal measurements; they are independent of LDR 1336 # Correction parameters for normal measurements; they are independent of LDR
1332 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1337 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1333 IS1 = np.array([IinPo, QinPo, UinPo, VinPo]) 1338 IS1 = np.array([IinPo, QinPo, UinPo, VinPo])
1334 IS2 = np.array([IinPa, QinPa, UinPa, VinPa]) 1339 IS2 = np.array([IinPa, QinPa, UinPa, VinPa])
1335 GT = np.dot(ATP, IS1) 1340 GT = np.dot(ATP, IS1)
1336 GR = np.dot(ARP, IS1) 1341 GR = np.dot(ARP, IS1)
1416 AT = ATF1 * IinF + ATF4 * h * VinF 1421 AT = ATF1 * IinF + ATF4 * h * VinF
1417 BT = ATF3e * QinF - ATF2e * h * UinF 1422 BT = ATF3e * QinF - ATF2e * h * UinF
1418 AR = ARF1 * IinF + ARF4 * h * VinF 1423 AR = ARF1 * IinF + ARF4 * h * VinF
1419 BR = ARF3e * QinF - ARF2e * h * UinF 1424 BR = ARF3e * QinF - ARF2e * h * UinF
1420 1425
1421 # Correction paremeters for normal measurements; they are independent of LDR 1426 # Correction parameters for normal measurements; they are independent of LDR
1422 if (not RotationErrorEpsilonForNormalMeasurements): 1427 if (not RotationErrorEpsilonForNormalMeasurements):
1423 GT = ATF1 * IinE + ATF4 * VinE 1428 GT = ATF1 * IinE + ATF4 * VinE
1424 GR = ARF1 * IinE + ARF4 * VinE 1429 GR = ARF1 * IinE + ARF4 * VinE
1425 HT = ATF2 * QinE - ATF3 * UinE - ATF4 * 2 * VinE 1430 HT = ATF2 * QinE - ATF3 * UinE - ATF4 * 2 * VinE
1426 HR = ARF2 * QinE - ARF3 * UinE - ARF4 * 2 * VinE 1431 HR = ARF2 * QinE - ARF3 * UinE - ARF4 * 2 * VinE
1438 AT = np.dot(ATFe, IF1e) 1443 AT = np.dot(ATFe, IF1e)
1439 AR = np.dot(ARFe, IF1e) 1444 AR = np.dot(ARFe, IF1e)
1440 BT = np.dot(ATFe, IF2e) 1445 BT = np.dot(ATFe, IF2e)
1441 BR = np.dot(ARFe, IF2e) 1446 BR = np.dot(ARFe, IF2e)
1442 1447
1443 # Correction paremeters for normal measurements; they are independent of LDR --- the same as for TypeC = 6 1448 # Correction parameters for normal measurements; they are independent of LDR --- the same as for TypeC = 6
1444 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1449 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1445 IS1 = np.array([IinE, 0, 0, VinE]) 1450 IS1 = np.array([IinE, 0, 0, VinE])
1446 IS2 = np.array([0, QinE, -UinE, -2 * VinE]) 1451 IS2 = np.array([0, QinE, -UinE, -2 * VinE])
1447 1452
1448 GT = np.dot(ATF, IS1) 1453 GT = np.dot(ATF, IS1)
1472 AT = np.dot(ATFe, IF1e) 1477 AT = np.dot(ATFe, IF1e)
1473 AR = np.dot(ARFe, IF1e) 1478 AR = np.dot(ARFe, IF1e)
1474 BT = np.dot(ATFe, IF2e) 1479 BT = np.dot(ATFe, IF2e)
1475 BR = np.dot(ARFe, IF2e) 1480 BR = np.dot(ARFe, IF2e)
1476 1481
1477 # Correction paremeters for normal measurements; they are independent of LDR 1482 # Correction parameters for normal measurements; they are independent of LDR
1478 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1483 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1479 # IS1 = np.array([IinE,0,0,VinE]) 1484 # IS1 = np.array([IinE,0,0,VinE])
1480 # IS2 = np.array([0,QinE,-UinE,-2*VinE]) 1485 # IS2 = np.array([0,QinE,-UinE,-2*VinE])
1481 IS1 = np.array([IinFo, 0, 0, VinFo]) 1486 IS1 = np.array([IinFo, 0, 0, VinFo])
1482 IS2 = np.array([0, QinFa, UinFa, VinFa]) 1487 IS2 = np.array([0, QinFa, UinFa, VinFa])
1560 AT = ATE1o * IinE + (ATE4o + aCal * ATE4a) * h * VinE 1565 AT = ATE1o * IinE + (ATE4o + aCal * ATE4a) * h * VinE
1561 BT = aCal * (ATE3ae * QinEe - ATE2ae * h * UinEe) 1566 BT = aCal * (ATE3ae * QinEe - ATE2ae * h * UinEe)
1562 AR = ARE1o * IinE + (ARE4o + aCal * ARE4a) * h * VinE 1567 AR = ARE1o * IinE + (ARE4o + aCal * ARE4a) * h * VinE
1563 BR = aCal * (ARE3ae * QinEe - ARE2ae * h * UinEe) 1568 BR = aCal * (ARE3ae * QinEe - ARE2ae * h * UinEe)
1564 1569
1565 # Correction paremeters for normal measurements; they are independent of LDR 1570 # Correction parameters for normal measurements; they are independent of LDR
1566 if (not RotationErrorEpsilonForNormalMeasurements): 1571 if (not RotationErrorEpsilonForNormalMeasurements):
1567 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F) 1572 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F)
1568 GT = ATE1o * IinE + ATE4o * h * VinE 1573 GT = ATE1o * IinE + ATE4o * h * VinE
1569 GR = ARE1o * IinE + ARE4o * h * VinE 1574 GR = ARE1o * IinE + ARE4o * h * VinE
1570 HT = ATE2a * QinE + ATE3a * h * UinEe + ATE4a * h * VinE 1575 HT = ATE2a * QinE + ATE3a * h * UinEe + ATE4a * h * VinE
1580 AT = ATE1 * IinE + ZiC * CosC * (ATE2e * QinEe + ATE4 * VinE) + ATE3e * UinEe 1585 AT = ATE1 * IinE + ZiC * CosC * (ATE2e * QinEe + ATE4 * VinE) + ATE3e * UinEe
1581 BT = DiC * (ATE1 * UinEe + ATE3e * IinE) + ZiC * SinC * (ATE4 * QinEe - ATE2e * VinE) 1586 BT = DiC * (ATE1 * UinEe + ATE3e * IinE) + ZiC * SinC * (ATE4 * QinEe - ATE2e * VinE)
1582 AR = ARE1 * IinE + ZiC * CosC * (ARE2e * QinEe + ARE4 * VinE) + ARE3e * UinEe 1587 AR = ARE1 * IinE + ZiC * CosC * (ARE2e * QinEe + ARE4 * VinE) + ARE3e * UinEe
1583 BR = DiC * (ARE1 * UinEe + ARE3e * IinE) + ZiC * SinC * (ARE4 * QinEe - ARE2e * VinE) 1588 BR = DiC * (ARE1 * UinEe + ARE3e * IinE) + ZiC * SinC * (ARE4 * QinEe - ARE2e * VinE)
1584 1589
1585 # Correction paremeters for normal measurements; they are independent of LDR 1590 # Correction parameters for normal measurements; they are independent of LDR
1586 if (not RotationErrorEpsilonForNormalMeasurements): 1591 if (not RotationErrorEpsilonForNormalMeasurements):
1587 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F) 1592 # Stokes Input Vector before receiver optics Eq. E.19 (after atmosphere F)
1588 GT = ATE1o * IinE + ATE4o * VinE 1593 GT = ATE1o * IinE + ATE4o * VinE
1589 GR = ARE1o * IinE + ARE4o * VinE 1594 GR = ARE1o * IinE + ARE4o * VinE
1590 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE 1595 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE
1611 AT = np.dot(ATEe, IE1e) 1616 AT = np.dot(ATEe, IE1e)
1612 AR = np.dot(AREe, IE1e) 1617 AR = np.dot(AREe, IE1e)
1613 BT = np.dot(ATEe, IE2e) 1618 BT = np.dot(ATEe, IE2e)
1614 BR = np.dot(AREe, IE2e) 1619 BR = np.dot(AREe, IE2e)
1615 1620
1616 # Correction paremeters for normal measurements; they are independent of LDR 1621 # Correction parameters for normal measurements; they are independent of LDR
1617 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out 1622 if (not RotationErrorEpsilonForNormalMeasurements): # calibrator taken out
1618 GT = ATE1o * IinE + ATE4o * VinE 1623 GT = ATE1o * IinE + ATE4o * VinE
1619 GR = ARE1o * IinE + ARE4o * VinE 1624 GR = ARE1o * IinE + ARE4o * VinE
1620 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE 1625 HT = ATE2a * QinE + ATE3a * UinE + ATE4a * VinE
1621 HR = ARE2a * QinE + ARE3a * UinE + ARE4a * VinE 1626 HR = ARE2a * QinE + ARE3a * UinE + ARE4a * VinE
1702 aX[2, iN] = HR 1707 aX[2, iN] = HR
1703 aX[3, iN] = HT 1708 aX[3, iN] = HT
1704 aX[4, iN] = K 1709 aX[4, iN] = K
1705 1710
1706 aLDRCal[iN] = iLDRCal 1711 aLDRCal[iN] = iLDRCal
1712 aDOLP[iN] = iDOLP
1707 aERaT[iN] = iERaT 1713 aERaT[iN] = iERaT
1708 aERaR[iN] = iERaR 1714 aERaR[iN] = iERaR
1709 aRotaT[iN] = iRotaT 1715 aRotaT[iN] = iRotaT
1710 aRotaR[iN] = iRotaR 1716 aRotaR[iN] = iRotaR
1711 aRetT[iN] = iRetT 1717 aRetT[iN] = iRetT
1784 # fig.savefig(LID + '_' + aVar + '.png', dpi=150, bbox_inches='tight', pad_inches=0) 1790 # fig.savefig(LID + '_' + aVar + '.png', dpi=150, bbox_inches='tight', pad_inches=0)
1785 # plt.close 1791 # plt.close
1786 return 1792 return
1787 1793
1788 1794
1795 if (nDOLP > 0): PlotSubHist("DOLP", aDOLP, DOLP0, dDOLP, iDOLP, nDOLP)
1789 if (nRotL > 0): PlotSubHist("RotL", aRotL, RotL0, dRotL, iRotL, nRotL) 1796 if (nRotL > 0): PlotSubHist("RotL", aRotL, RotL0, dRotL, iRotL, nRotL)
1790 if (nRetE > 0): PlotSubHist("RetE", aRetE, RetE0, dRetE, iRetE, nRetE) 1797 if (nRetE > 0): PlotSubHist("RetE", aRetE, RetE0, dRetE, iRetE, nRetE)
1791 if (nRotE > 0): PlotSubHist("RotE", aRotE, RotE0, dRotE, iRotE, nRotE) 1798 if (nRotE > 0): PlotSubHist("RotE", aRotE, RotE0, dRotE, iRotE, nRotE)
1792 if (nDiE > 0): PlotSubHist("DiE", aDiE, DiE0, dDiE, iDiE, nDiE) 1799 if (nDiE > 0): PlotSubHist("DiE", aDiE, DiE0, dDiE, iDiE, nDiE)
1793 if (nRetO > 0): PlotSubHist("RetO", aRetO, RetO0, dRetO, iRetO, nRetO) 1800 if (nRetO > 0): PlotSubHist("RetO", aRetO, RetO0, dRetO, iRetO, nRetO)
1861 #plt.show() 1868 #plt.show()
1862 #fig.savefig(LID + '_' + aVar + '.png', dpi=150, bbox_inches='tight', pad_inches=0) 1869 #fig.savefig(LID + '_' + aVar + '.png', dpi=150, bbox_inches='tight', pad_inches=0)
1863 #plt.close 1870 #plt.close
1864 return 1871 return
1865 1872
1873 if (nDOLP > 0): PlotSubHistF11("DOLP", aDOLP, DOLP0, dDOLP, iDOLP, nDOLP)
1866 if (nRotL > 0): PlotSubHistF11("RotL", aRotL, RotL0, dRotL, iRotL, nRotL) 1874 if (nRotL > 0): PlotSubHistF11("RotL", aRotL, RotL0, dRotL, iRotL, nRotL)
1867 if (nRetE > 0): PlotSubHistF11("RetE", aRetE, RetE0, dRetE, iRetE, nRetE) 1875 if (nRetE > 0): PlotSubHistF11("RetE", aRetE, RetE0, dRetE, iRetE, nRetE)
1868 if (nRotE > 0): PlotSubHistF11("RotE", aRotE, RotE0, dRotE, iRotE, nRotE) 1876 if (nRotE > 0): PlotSubHistF11("RotE", aRotE, RotE0, dRotE, iRotE, nRotE)
1869 if (nDiE > 0): PlotSubHistF11("DiE", aDiE, DiE0, dDiE, iDiE, nDiE) 1877 if (nDiE > 0): PlotSubHistF11("DiE", aDiE, DiE0, dDiE, iDiE, nDiE)
1870 if (nRetO > 0): PlotSubHistF11("RetO", aRetO, RetO0, dRetO, iRetO, nRetO) 1878 if (nRetO > 0): PlotSubHistF11("RetO", aRetO, RetO0, dRetO, iRetO, nRetO)

EARLINET Repository: Lidar polarisation GHK-parameter and error calculation

mercurial