1 Modulation Format Independent Joint Polarization and Phase Tracking for Coherent Receivers Cristian B. Czegledi, Erik Agrell, Senior Member, IEEE, Magnus Karlsson, Senior Member, IEEE, Fellow, OSA, and Pontus Johannisson, Member, OSA Abstract—The state of polarization and the carrier phase noise and drift of the state of polarization (SOP), which have drift dynamically during transmission in a random fashion in to be dynamically tracked in the receiver [1], [2]. 6 coherent optical fibercommunications. Thetypical digital signal 1 The phase and SOP tracking are important DSP blocks at processingsolutiontomitigatetheseimpairmentsconsistsoftwo 0 the receiver and are different from the chromatic dispersion separateblocksthattrackeachphenomenonindependently.Such 2 algorithms have been developed without taking into account compensation, which can be set once and then forgotten due n mathematical models describing the impairments. We study a to its static behavior. The SOP drift has its origin in the a blind, model-based tracking algorithm to compensate for these imperfectionsofthemanufacturingprocessofthefibercables, J impairments. The algorithm dynamically recovers the carrier mechanical/thermalstress on the deployed fibers, splices, etc. 2 phase and state of polarization jointly for an arbitrary modu- Duetotheserandomvariations,theSOPchangesdynamically 2 lation format. Simulation results show the effectiveness of the proposed algorithm, having a fast convergence rate and an in time and along the fiber, which makes it difficult to fully ] excellent tolerance to phase noise and dynamic drift of the compensate for. The phase noise originates from the finite T polarization. The computational complexity of the algorithm is coherence length of the transmitter and receiver lasers and I lowercomparedtostate-of-the-artalgorithmsatsimilarorbetter it drifts in time as a Wiener random walk. Despite the fact s. performance,whichmakesitastrongcandidateforfutureoptical that the SOP drift and the phase noise arise due to different c systems. [ hardware imperfections, they can be modeled jointly as dy- Index Terms—Coherent optical fiber communication, model- namicrotationsoftheopticalfield[3].Adeterministicorstatic 1 based, phase noise, phase recovery, polarization demultiplexing, behavior of these phenomena would be straightforward to v polarization drift, polarization recovery. 7 resolve,butwhentheimpairmentsdriftrandomly,thereceiver 2 must adjust dynamically to track the phenomena. 9 I. INTRODUCTION The commonDSP solution for SOP trackingis done in the 5 Jonesspaceusingtheconstantmodulusalgorithm(CMA)[2], 0 DIGITAL signal processing (DSP) enables spectrally effi- initially developed for two-dimensional modulation formats 1. cient communications based on coherent transmission. [4], or modified versions of it to accommodate for various 0 Contrary to traditional optical transmission links that are modulation formats, such as the multiple modulus algorithm 6 based on intensity-modulation and direct-detection, coherent (MMA) [5], [6] or the polarization-switched (PS-)CMA [7]. 1 transmissions carry the information in both the intensity and : Alternatively, the polarization demultiplexing can be done in v phase of the optical field, in both polarizations, and benefit the Stokes space [8], [9], which in addition also aligns the i from improvedsensitivities, higher-ordermodulation formats, X phase of the two polarizations, thus improving the phase and digital impairment mitigation. Polarization-multiplexed r tracking by enabling joint phase estimation over the two a quadraturephase-shiftkeying(PM-QPSK) introducedfor 100 polarizations. In general, the phase tracking is performed Gb/s transmission has been widely deployedand reachedma- independently of the SOP tracking, using algorithms such turity. Recently, 200 Gb/s transceivers have been made com- as the Viterbi–Viterbi algorithm [10] or the blind phase- mercially available based on 16-ary polarization-multiplexed search algorithm (BPS) [11], which treat each polarization quadrature amplitude modulation (PM-16-QAM) and it is independently. expected that in the near future, higher-order PM-M-QAM Recently, Louchet et al. proposed the Kabsch algorithm modulation formats will become a necessity for higher data [12], which addresses the two impairments jointly in both rates. However, the improved spectral efficiency comes at the polarizations in the real four-dimensional (4D) space and costofareducedtolerancetoimpairmentssuchaslaserphase accommodates any modulation format. In general, joint es- timation leads to better performance, and it is expected that C. B. Czegledi and E. Agrell are with the Department of Signals and future transceivers will benefit from improved performance Systems,ChalmersUniversityofTechnology,SE-41296Gothenburg,Sweden (e-mail: [email protected]). from such integrations of different DSP blocks [13]. M. Karlsson and P. Johannisson are with the Department of Microtech- However, very few algorithms present in the literature take nology and Nanoscience, Chalmers University of Technology, SE-412 96 into consideration analytical models describing the behavior Gothenburg, Sweden. This work was supported by the Swedish research council (VR) and of the impairments. Model-based algorithms have a restricted performedwithintheFiberOpticCommunicationsResearchCenter(FORCE) flexibility and therefore fewer degrees of freedom (DOFs) to atChalmers.ThesimulationswereperformedonresourcesatChalmersCentre adjust. The DOFs of model-based algorithms are restricted forComputationalScienceandEngineering(C3SE)providedbytheSwedish National Infrastructure forComputing (SNIC). to only the ones covered by the impairment to compensate 2 for,thusresultingin a moreefficientimpairmentcancellation, The received discrete symbols, at distance L, are obtained rather than scanning over a larger domain in order to find the from the receivedelectric field E(L,t) after matched filtering optimal setup. and sampling Inthispaper,weproposeamodel-basedalgorithmtojointly ∞ recover the carrier phase and SOP for arbitrary modulation rk = E(L,t)p∗(t kT)dt. (3) − formats. The design of the algorithm is based on a channel Z−∞ model (described in Section II) that can emulate temporal The discrete transmitted symbols u are drawn indepen- k stochastic polarization and phase drifts, and has been suc- dently from a finite constellation C = c ,c ,...,c with 1 2 M { } cessfully validated with data measured on installed fibers equal probability. The average energy per symbol is the [3]. The algorithm (described in Section III) uses a non- average of u 2 and in this case equals k k k dataaideddecision-directedarchitecture,hencezerooverhead, M and operates jointly on both polarizations. The performance E = 1 c 2. (4) s k of the algorithm is investigated in Section IV by comparing M k k k=1 it with state-of-the-art algorithms for different modulation X Assuming that the chromatic dispersion has been success- formats, whereas the complexity is evaluated in Section V. fully compensated for and polarization-dependent losses and Theproposedalgorithmperformssimilarlyorbetterthanstate- polarizationmodedispersionarenegligible,thepropagationof of-the-art algorithms and provides a good trade-off between the opticalfield can be described by a unitary 2 2 complex- complexity and performance regardless of the modulation × valued Jones matrix J [14, p. 18]. The received symbol format. High performance and fast convergence rate at low k r C2, in the presenceofopticalamplifiernoise, SOPdrift, complexity, for any modulation format, make the algorithm k ∈ and laser phase noise, can be related to the input u as a strong candidate for future elastic optical systems, where k the modulation format can be changed dynamically during rk =e−iφkJkuk+nk, (5) transmissiontoaccommodateforvariouschannelandnetwork conditions. where i = √ 1, φk models the carrier phase noise, and n C2 deno−tes the additive noise, which is represented Thefollowingnotationconventionsareusedthroughoutthe k ∈ paper: column vectors are denoted by bold lower case (e.g., by two independent complex circular zero-mean Gaussian u) and matrices by bold upper case (e.g., U), except a few randomvariableswith varianceN0/2 per realdimension,i.e., specific cases, for literature consistency reasons, denoted by E[nknHk]=N0I2 [15]. smallGreekletterssuchasthePaulimatricesσσσ ,the4Dbasis The phase noise is modeled as a Wiener process [11], [16] i matrices ρρρ , λλλ , and the electric field Jones vector E. Trans- • i i φ =φ +φ , (6) position is written as uT, conjugation as u , and conjugate k+1 k k ∗ transposeasuH.Then nidentitymatrixiswrittenasI and whereφ• is theinnovationofthe phasenoise.Theinnovation n k theexpectationoperator×asE[].Thedotoperationααα σ~σσ should φ• is a random variable drawn independently at each time k · · be interpreted as a linear combination of the three matrices instance k from a zero-mean Gaussian distribution forming the tensor σ~σσ = (σσσ ,σσσ ,σσσ ). Multiplication of a 1 2 3 φ• (0,σ2), (7) matrix with the tensor σ~σσ result in a tensor with element-wise k ∼N ν multiplications, e.g., Uσ~σσ = (Uσσσ1,Uσσσ2,Uσσσ3). The absolute where the variance σ2 = 2π∆νT, ∆ν is the sum of the ν value is denoted by and the Euclidean norm by . linewidths of the transmitter and receiver lasers. The initial |·| k·k phase φ is modeled as a random variable uniformly dis- 0 tributed in the interval [0,2π). II. DISCRETE-TIMECHANNEL MODEL The time evolution of the SOP drift can be emulated by modeling J as a sequence of random Jones matrices [3] The coherent optical signal has two quadratures in two k polarizations and can be described by a Jones vector J =J(ααα• )J , (8) k+1 k k E(z,t)= Ex(z,t) , (1) whereJ(ααα•k)istheSOPinnovationmatrix(cf.φ•k in(6)).The Ey(z,t) matrix function J(ααα) is defined using the matrix exponential (cid:18) (cid:19) [17, p. 165] parameterized by three DOFs ααα [18] as at propagationdistance z andtime t, where E (E ) is the x- x y polarized(y-polarized)electricfield,representedasacomplex J(ααα)=exp( iααα σ~σσ), (9) − · baseband signal. The linearly modulated transmitted electric where ααα=(α ,α ,α ) is a three-componentreal vector and field into the transmission medium is 1 2 3 σ~σσ = (σσσ ,σσσ ,σσσ ) is a tensor of the Pauli spin matrices [14, 1 2 3 eq. (2.5.19)] E(0,t)= u p(t kT), (2) k − k 1 0 0 1 0 i X σσσ1 = 0 1 , σσσ2 = 1 0 , σσσ3 = i −0 . (10) where uk C2 are the information symbols for k Z, T (cid:18) − (cid:19) (cid:18) (cid:19) (cid:18) (cid:19) ∈ ∈ is the symbol (baud) interval, and p(t) is a real-valued pulse The vector ααα can be expressed as a product ααα=θa of its shape. length θ = ααα and the unit vector a = (a ,a ,a ), which 1 2 3 k k 3 joint equalization CD and/or Frequency Symbol Front-end Carrier phase nonlinearity offset SOPrecovery estimation correction recovery compensation estimation and decoding Figure1. Receiver blockdiagram withelementary DSPmodules. represents its direction on the unit sphere. Based on this on the Jones formalism rely on complex two-dimensional decompositionofααα, (9)can be rewritteninto anexplicitform vectors and matrices that have four DOFs. This description issufficientforwavepropagationsinceitcancoveranylinear J(ααα)=I cosθ ia σ~σσsinθ. (11) 2 − · phenomenon that can arise during photon propagation. The Since the transformation J(ααα) is unitary, the inverse can Stokes description is preferred in some situations since the be foundby the conjugatetranspose operationor negatingthe Stokesvectorsare observablequantitiesandcan bevisualized argument, J(ααα) 1 =J(ααα)H =J( ααα). as points on a three-dimensional sphere, called the Poincare´ − − TherandomnatureoftheSOPdriftisemulatedbydrawing sphere. In this case, the channel matrix J is replaced by a k the three innovation parameters ααα• of the innovation J(ααα• ) 3 3 Mueller matrix M with three DOFs that models only k k × k independently from a zero-mean real Gaussian distribution at the changesof the SOP. The Stokes description cannotmodel each time instance k absolute phase shifts, therefore it is immune to phase noise. The real 4D formalism models the channel behavior using ααα• (0,σ2I ), (12) k ∼N p 3 4 4 real rotation matrices that have six DOFs, which span × where σ2 = 2π∆pT and ∆p is referred as the polarization over a richer space than the Jones (four DOFs) or Mueller p linewidth, which quantifies the speed of the SOP drift, analo- (three DOFs) matrices can. However, only four of them are gously to the linewidth describing the phase noise, cf. (7). physically realizable for propagating photons and the other The initial state of the channel J = J(ααα ) is formed two can be synthesized using DSP [20]. 0 0 from the vector ααα = θa, which is identified from the In the remainder of this section, we will provide the 0 unitvector(cosθ,a sinθ,a sinθ,a sinθ)T =g/ g where derivation of the proposed tracking algorithm using the Jones 1 2 3 k k g (0,I ). This ensures that J u is uniformly distributed description, after which the equivalent algorithms using the 4 0 ∼N over all possible SOP for a fixed u [3], [19]. Stokes and the real 4D descriptions are given. The details of The phase noise and the SOP drift can be combined into a the channel model descriptions using the Stokes and real 4D single operation Hk =e−iφkJk and (5) can be rewritten as descriptions are omitted and can be found in [3]. r =H u +n . (13) k k k k A. Jones Description The update of H can be expressed analogously to (8) as k The considered DSP setup is shown in Fig. 1, where H =H(φ• ,ααα• )H , (14) we combine the SOP drift and carrier phase tracking, i.e., k+1 k k k estimation of H , into a single block after the frequency where the phase innovation φ• and the random vector ααα• k k k offset compensation, which can be done in this case using are defined as (7) and (12), respectively. The matrix function spectrum-basedmethods[21],[22].ConsideringthatH does H(φ,ααα) can be expressed as k not change significantly over a symbol duration, we estimate H(φ,ααα)=e−iφJ(ααα)=exp( i(ααα σ~σσ+φI2)) (15) the transmitted symbol from Hˆ−k1rk based on a previously − · calculated estimate of the channel matrix Hˆ using the mini- =(cosφ isinφ)(I cosθ ia σ~σσsinθ), (16) k 2 − − · mum Euclidean distance criterion whichcombinestheeffectsofbothphasenoiseandSOPdrift. 2 uˆk =arcgmCin Hˆ−k1rk−c . (17) III. POLARIZATION AND PHASETRACKING ALGORITHM ∈ (cid:13) (cid:13) In order to successfully decode the data at the receiver, Thereafter Hˆk is updated as (cid:13)(cid:13) (cid:13)(cid:13) the channel matrix Hk (or, equivalently, φk and Jk) needs Hˆ =H(φˆ ,αααˆ )Hˆ , (18) to be estimated and tracked during transmission such that k+1 k k k it is possible to accurately estimate uk from the received whereφˆk andαααˆk areestimatesofφ•k andααα•k. Theseestimates sampler .Thissectionprovidesadescriptionoftheproposed are calculated such that the error function k algorithm to estimate H , first using the Jones formalism, k 1 2 thereafteralternativesusingtheStokesandreal4Dformalisms e = H(φ,ααα)Hˆ − r uˆ , (19) k k k− k are given. (cid:13) (cid:13) (cid:13)(cid:16) (cid:17) (cid:13) TheJonesdescriptioncanbe replacedby theStokesor real is minimized with(cid:13)respect to φ and ααα, i.e., (cid:13) (cid:13) (cid:13) 4D descriptions, which can provide benefits in some situa- [φˆ ,αααˆ ]=argmine . (20) tions [3]. The analytics describing wave propagation based k k k φ,ααα 4 This can be achieved by computing φˆ and αααˆ using the k k 2π gradient descent method [23, p. 466] 3 φˆk = µφRe ∂ek φ=0 (21) π2 − ∂φ(cid:12)ααα=[0,0,0]T! (cid:12) =−2µφRe i(Hˆ(cid:12)(cid:12)−k1rk−uˆk)HHˆ−k1rk , (22) 2π 7 (cid:16) (cid:17) 1.6 1.7 1.8 1.9 2 2.1 ×104 αααˆ = µ Re e (23) k − ααα ∇ααα k(cid:12)αααφ==0[0,0,0]T! π Σφ•k, Σααα•k - true (cid:12) 2 k k =−2µαααRe i(Hˆ−k(cid:12)(cid:12)1rk−uˆk)HHˆ−k1σ~σσrk , (24) π4 Σkφˆk, Σkαααˆk - estimate (cid:16) (cid:17) α1 where µ and µ are positive tracking step sizes of the phase 0 φ ααα and of the SOP parameters, respectively,which determine the −π4 α2 speed of the algorithm’s convergence, the tracking accuracy, −π 2 φ and the rate at which changes in the channel can be tracked. −3π Thederivationsof(22)and(24)canbefoundintheAppendix. 4 α3 Both innovation parameters φ• and ααα• have zero mean by −π k k (7) and (12); therefore the partial derivatives in (21) and (23) 0 1 2 3 4 5 6 7 8 9 10 are evaluated at φ = 0,ααα = [0,0,0]T, which results in no Symbolindex k ×104 preferred direction of φˆ and αααˆ . Evaluating the gradient at k k Figure 2. Tracked channel parameters using the proposed algorithm with non-zero values could compensate for constant offsets; e.g., ∆ν=1MHzand∆p=1kHzat28GbaudPM-16-QAMareshown.Ascan using φ=0 could compensate for frequency offsets. be seen, the algorithm has excellent tracking capabilities without exhibiting 6 It is important to note that H(φ,ααα) is a many-to-one cycleslips. function. Therefore φˆ and αααˆ are not necessarily equal to k k φ• ,ααα• and can have different values, but resulting in same k k B. Stokes Description matrix H(φ,ααα). Since the typical drift time of the SOP is slower ( 1 ms In this case, the propagation of the electric field can be ∼ or larger) [3], [24], [25] than the drift of the phase noise ( modeled as [3] ∼ 1 µs) [11], we chose the tracking steps µ ,µ of the two phenomenadifferently.Forthesamereasons,φtheαααupdateofthe srk =Mksuk +snk, (25) SOPestimateαααˆ canbedonelessoftenthantheupdateofthe absolute phase kestimate φˆ , which will result in a decreased where suk =uHkσ~σσuk and srk =rHkσ~σσrk are the corresponding k Stokes vectors of u and r [14, eq. (2.5.26)]. The noise DSP complexity.In this case, φˆ can be calculated using (22) k k at every time instance k and αααˆk using (24) only at every P term is snk = (Hkuk)Hσ~σσnk + nHkσ~σσHkuk + nHkσ~σσnk and k can be identified by equating terms after substituting (13) in symbols, otherwise should be set to [0,0,0]T. Fig. 2 shows an example of the algorithm’s tracking capa- srk = rHkσ~σσrk. As can be noted, snk is signal dependent and (31) is not an additive noise model, opposed to (13). The bility, where we compare the sum of the innovations φ• ,ααα• k k channel matrix M modeling the evolution of the SOP can with their estimates φˆ ,αααˆ obtained using the proposed al- k k k be expressed using a 3 3 Mueller matrix defined as [18] gorithm. Even though, as mentioned above, φˆ ,αααˆ do not × k k have to follow φ•k,ααα•k to obtain a good estimate of the M(ααα)=exp(2 (ααα)), (26) channel matrix H , the algorithm manages to obtain similar K k parameters efficiently without exhibiting cycle slips over the where (ααα) denotes the cross-productoperator [20, eq. (11)] 105simulatedsymbols.Perhapsdifferentvaluesofφˆ ,αααˆ may K k k be obtained if the initial values φˆ ,αααˆ are not set to be the 0 ααα ααα 0 0 − 3 2 same as φ0,ααα0. Note that the plotted parameters are just for K(ααα)= ααα3 0 −ααα1. (27) demonstration purposes and do not reflect the behavior of ααα ααα 0 2 1 H or Hˆ , since H = H(Σφ• ,Σααα• ) because in general − k k k 6 k k k k The inverse can be obtained as M(ααα) 1 = M(ααα)T = H(φ +φ ,ααα +ααα )=H(φ ,ααα )H(φ ,ααα ). − 1 2 1 2 6 1 1 2 2 M( ααα). The polarization transformations introduced by M Inaprestudyforthiswork[26],weinvestigatedasimilaral- − k canbeseenasrotationsofthePoincare´ spherearoundtheunit gorithmthattracksboththecarrierphaseandtheSOPjointly. vector a by an angle 2θ. The main difference between the two algorithms consists in the updating rule of Hˆ , which in (18) is done as a matrix AnalogouslytotheJonesdescription,thealgorithmdecides k first which was the transmitted Stokes vector based on the multiplication. However, the method in [26] was based on minimum Euclidean distance criterion a different system model that did not reflect random SOP drift accurately; therefore it is not suited for installed fiber 2 transmissions. ˆsuk =arcgmCin Mˆ −k1srk −sc , (28) ∈ (cid:13) (cid:13) (cid:13) (cid:13) (cid:13) (cid:13) 5 using the inverted estimate of Mk, where1 sc = cHσ~σσc. where, analogously with (21)–(24) and (30), Thereafter, Mˆ is updated analogously to (18) as k Mˆ =M(αααˆ )Mˆ . (29) φˆk =2µφ(Rˆ−k1rk−uˆk)TRˆ−k1λλλ1rk, (37) k+1 k k αααˆk =−2µααα(Rˆ−k1rk−uˆk)TRˆ−k1~ρρρrk. (38) Analogously with (20)–(24), it can be shown that the optimal In the above description, only four DOFs, i.e., the scalars αααˆ =[αˆ ,αˆ ,αˆ ]T are computed as k k,1 k,2 k,3 φ, α , α , and α corresponding to λλλ , ρρρ , ρρρ , and ρρρ , 1 2 3 1 1 2 3 αˆk,i =4µααα(Mˆ −k1srk −ˆsuk)TMˆ −k1K(ei)srk (30) rsepsopnedcttiovetlhye, ocfartrhieermphatarsiexaRn(dφt,hαααe)SwOePredurisfet.d,Thwehiocthhecrotrwreo- for i=1,2,3, using the gradient descent algorithm such that DOFs, i.e., the scalars corresponding to λλλ ,λλλ in (32), can 2 3 the Euclidean distance in the Stokes space is minimized. The beusedtocorrectcertaintransmitterand/orreceiverhardware vectors e form the standard basis in R3. i imperfections, which cannot be done using Jones or Stokes Note that the error function that minimizes the Euclidean formalisms, such as 90 I/Q error or the time skew between I ◦ distance is not optimum in this case since sn includes not and Q [29]. k onlynoise,butsignal–noiseinteraction.Furthermore,thenoise Thealgorithmpresentedinthissectionisfullyequivalentto termnHkσ~σσnk isnon-Gaussian.Bettermetrics[8]thattakeinto the one in Section III-A and will have the same performance, accountthenon-Gaussiandistributionofthenoisecanbeused, but not to the one in Section III-B. The latter may have a but it is outside the scope of this work. differentperformanceduetothesuboptimalerrorfunctionand a separate solution to mitigate the phase noise is required in C. 4D Real Description the latter case since it will not be covered by Mˆ . k Inthe4Dformalism,thephaseandSOPdriftsarecombined and modeled using a 4 4 real orthogonal matrix R [20], IV. RESULTS k × [27], [28] as We evaluated the achievable performance of the proposed vrk =Rkvuk +vnk, (31) recovery algorithm numerically. The details of the simulation setuparedescribedinSectionIV-A,whereasinSectionsIV-B where vz for any z = [z1,z2]T C2 is defined as to IV-E various performance metrics of the algorithm are ∈ [Re(z ),Im(z ),Re(z ),Im(z )]T.ThechannelmatrixR can 1 1 2 2 k evaluated. be expressed using the matrix function [20] R(φ,ααα)=exp((φ,0,0) λ~λλ ααα ~ρρρ), (32) A. Simulation Setup · − · where~ρρρ=(ρρρ ,ρρρ ,ρρρ ) and ~λλλ=(λλλ ,λλλ ,λλλ ) are six2 constant The considered modulation formats are PS-QPSK, PM- 1 2 3 1 2 3 QPSK, PM-16-QAM, PM-64-QAM, and PM-256-QAM at a basis matrices [20, eqs. (20)–(25)] symbol rate of 28 Gbaud. The performance is quantified by 0 1 0 0 0 0 0 1 countingthe numberof erroneous(4D) symbolsto obtain the − − 1 0 0 0 0 0 1 0 symbol error rate (SER) at the receiver for various setups in ρρρ = , ρρρ = , (33) 1 0 0 0 1 2 0 1 0 0 thepresenceoflaserphasenoise,SOPdrift,andadditivewhite − 0 0 1 0 1 0 0 0 Gaussian noise (AWGN). The latter is quantified through the − 0 0 1 0 0 1 0 0 signal-to-noise (SNR) ratio defined as SNR = Es/N0. To 0 0 0 1 1 0 0 0 remove phase ambiguities differential coding [30, Sec. 2.6.1] ρρρ3 = 1 0 0 0, λλλ1 =−0 0 0 1. (34) was employed independently in each polarization. Note that −0 1 0 0 0 0 1 0 the presented results when only AWGN is considered (for − − comparison reasons) still imply differential coding. The pre- The inverse can be obtained as R(φ,ααα)−1 = R(φ,ααα)T = sented results are obtained using the Jones description of the R( φ, ααα). algorithm from Section III-A. − − Analogously to the Jones/Stokes description, the algorithm The proposed algorithm was implemented such that both decides first which was the transmitted symbol based on the (22) and (24) were calculated for every k, i.e., P = 1. minimum Euclidean distance criterion The tracking step size µ and µ were chosen for each set φ ααα 2 of system parameters according to the heuristically obtained vˆuk =arcgmCin Rˆ−k1vrk −vc , (35) relations ∈ (cid:13) (cid:13) √∆νTc usingtheinvertedestimateof(cid:13)(cid:13)Rk.Theestima(cid:13)(cid:13)teofthechannel µφ = E , (39) matrix Rˆ is updated analogously to (18) as s k √∆pTc Rˆk+1 =R(φˆk,αααˆk)Rˆk, (36) µααα = E , (40) s 1Note that for constellations with rotational symmetry, more than one where c is a constant given in Table I, which depends on the constellation point will correspond to the same Stokes vector; e.g., the PM- modulationformat.Theconstantcwasoptimizedtoensurethe QAMmodulation formathas afour-fold rotational symmetry, therefore four bestperformanceinthesteady-stateregimeatSER=10 3.In distinct constellation pointscwillcorrespond tothesameStokesvector sc. − someapplications,thelinewidthparametersoftheopticallink 2Thematricesλλλ2andλλλ3arenotshownsincetheydonotinfluenceR(φ,ααα) in(32). maybeunknownatthereceiver,andthereforeitisimpossible 6 TableI ALGORITHMPARAMETERS,ACHIEVABLEPERFORMANCE,ANDHARDWARECOMPLEXITY Algorithm Max. tol. Max.tol.∆p Max. tol. Max.tol.∆ν Operations Comparisons Memory parameters ∆p·T at28Gbaud ∆ν·T at28Gbaud units Kabsch NKab=31 0.34·10−4 0.95MHz 0.91·10−4 2.55MHz 132 7 0.5 NBPS=13 PS-QPSK PS-CMA+BPS PBPS=32 0.33·10−4 0.93MHz 6.67·10−4 18.68MHz 1320 267 856 µCMA=0.04/Es2 Proposedalg. c=27 3.20·10−4 8.96MHz 11.5·10−4 32.2MHz 346 7 8 Kabsch NKab=16 0.17·10−4 0.47MHz 0.79·10−4 2.23MHz 171 4 1 NBPS=19 PM-QPSK CMA+BPS PBPS=32 0.37·10−4 1.04MHz 6.98·10−4 19.54MHz 2936 304 1264 µCMA=0.16/Es2 Proposedalg. c=64 1.17·10−4 3.28MHz 9.06·10−4 25.37MHz 346 4 8 Kabsch NKab=16 0.44·10−5 122.1kHz 0.15·10−4 0.412MHz 171 12 1 NBPS=19 PM-16-QAM MMA+BPS PBPS=32 0.14·10−5 37.8kHz 1.48·10−4 4.14MHz 2940 840 1264 µMMA=0.04/Es2 Proposedalg. c=400 2.48·10−5 694.4kHz 1.35·10−4 3.78MHz 346 12 8 Kabsch NKab=16 0.66·10−6 18.4kHz 0.32·10−5 90.2kHz 171 28 1 NBPS=19 PM-64-QAM MMA+BPS PBPS=64 0.11·10−6 3.1kHz 3.06·10−5 856.8kHz 5768 3806 2480 µMMA=0.035/Es2 Proposedalg. c=2352 4.59·10−6 128.5kHz 2.54·10−5 711.2kHz 346 28 8 Kabsch NKab=16 0.19·10−6 5.3kHz 0.74·10−6 20.6kHz 171 60 1 NBPS=19 PM-256-QAM MMA+BPS PBPS=64 8.42·10−9 0.2kHz – – 5814 8058 2480 µMMA =0.017/Es2 Proposedalg. c=6084 1.22·10−6 34.1kHz 6.30·10−6 176.4kHz 346 60 8 to accurately compute µ and µ . To overcomethis problem, H isunknowntothereceiver.Theachievedresultsareshown φ ααα 0 the ∆p and ∆ν parameters should be overestimated to be on inFig.3,whereeachrowcorrespondstooneofSectionsIV-B– the safe side. Typically, overestimating these parameters does IV-Eandeachcolumncorrespondstoamodulationformat.We not lead to considerable degradation in the performance. omitted results obtained for PS-QPSK and PM-QPSK in the figureforspace reasons,and these aresummarizedin Table I. For comparison, results obtained by the Kabsch algo- rithm [12] and combinations of the (PS-)CMA/MMA [2], [6], [7] and BPS [11] algorithms are presented. The Kabsch B. Polarization Noise Tolerance algorithm operates simultaneously on both polarizations in a In this section, the ability to track time-varying SOP drift decision-directedblock-wisefashion,usingarectangularfilter is evaluated. To measure the polarization sensitivity, the po- of size equalto N =16 [12],exceptfor PS-QPSK, where Kab larization linewidth ∆p is varied while keeping the other it was set to NKab = 31 since we observed better results parameters fixed. The SNR is set such that SER = 10 3 − with a longer block size. The BPS algorithm uses a sliding- is achieved in an AWGN scenario and the accumulated laser window technique in each polarization independently since a linewidth ∆ν is chosen such that ∆ν T = 3.6 10 5, 3.6 − relative phase offset between the two polarization may occur. 10 5, 0.36 10 5, 0.18 10 5, 0.04 ·10 5{ (co·rresponding· − − − − The length of the window was set to N =19 [11], except · · · } BPS to ∆ν = 1000, 1000, 100, 50, 10 kHz at 28 Gbaud) for forPS-QPSK, whereforthe sameconsiderationsasabovethe { } PS-QPSK, PM-QPSK, PM-16-QAM, PM-64-QAM, and PM- filterlengthwassettoN =13.Thenumberoftestphases BPS 256-QAM, respectively. of the BPS algorithm was set to P = 32 for PS-QPSK, BPS The results of the simulation are shown in the top row PM-QPSK, PM-16-QAM and P = 64 for PM-64-QAM BPS of Fig. 3, where the SER is plotted versus the polarization and PM-256-QAM [11]. The convergence parameter of the linewidthtimesthesymboltime.Ascanbeseen,theproposed (PS-)CMA/MMA algorithm µ was optimized to tolerate CMA algorithm offers the best tolerance to SOP drift at any ∆p T themostpolarizationnoise andthe chosenvaluesare listed in · forallmodulationformats.HightolerancetoSOPdriftenables Table I. theuseofhigh-ordermodulationformatsevenonrapidlyvary- In Sections IV-B to IV-D, the SOP and phase tracking ingchannels,suchasaerialfibers.Bothcompetitoralgorithms capabilities of the algorithms are evaluated and therefore the showerrorfloorshigherthantheproposedalgorithm,possibly initial channel matrix H is considered to be known to the due to the block-wise operation of the Kabsch algorithm and 0 receiver, i.e., Hˆ = H . Section IV-E presents results on the finitenumberoftestphasesoftheBPSalgorithm,respectively. 0 0 convergencerate of the algorithms, where the channel matrix Perhaps the performanceof both competitor algorithmscould 7 MMA+BPS * Kabsch ut Proposedalgorithm AWGN 100 PM-16-QAM * *ut *ut 100 PM-64-QAM * *ut *ut *ut 100 PMbb-256-QAM * *ut *ut *ut *ut ∆ν·T=0.36·10−5 ∆ν·T=0.18·10−5 ∆ν·T=0.04·10−5 10−1 Es/N0=21.28dB 10−1 Es/N0=27.48dB 10−1 Es/N0=33.55dB R E ut ut S * * * ut 10−2 10−2 10−2 * * * 10−3*ut *ut bbbbbb *ut *ut *ut ut ut 10−3*ut *ut bbbbbb *ut *ut ut ut 10−3*ut *ut bbbb *ut ut ut 10−8 10−7 10−6 10−5 10−4 10−3 10−8 10−7 10−6 10−5 10−4 10−3 10−8 10−7 10−6 10−5 10−4 10−3 ∆p·T ∆p·T ∆p·T 100 PM-16-QAM * *ut *ut 100 PM-64-QAM * *ut *ut *ut 100 PM-256-Qbb AM * *ut *ut *ut *ut ∆p·T=3.57·10−8 * ∆p·T=3.57·10−8 ∆p·T=3.57·10−8 10−1 Es/N0=21.28dB 10−1 Es/N0=27.48dB 10−1 Es/N0=33.55dB R ut * ut SE10−2 10−2 10−2 * ut * * * 10−3*ut *ut *ut *ut *utbbbbbb ut ut 10−3*ut *ut *ut *ut bbbbbb ut ut 10−3*ut *ut *ut bbbb ut ut 10−8 10−7 10−6 10−5 10−4 10−3 10−2 10−8 10−7 10−6 10−5 10−4 10−3 10−2 10−8 10−7 10−6 10−5 10−4 10−3 10−2 ∆ν·T ∆ν·T ∆ν·T 101−010***ututut *ut *ut ut* ut* 101−010***ututut *ut *ut ut* ut* 101−010***ututut *ut *ut ut* ut* bb R10−2 ut * ut * 10−2 ut * ut * 10−2 ut * ut * SE10−3 bbbbbbut * 10−3 utbbbbbb * 10−3 utbbbb * 10−4 P∆pM·-T16=-Q3.A57M·10−8 ut * 10−4 P∆pM·-T64=-Q3.A57M·10−8 ut * 10−4 P∆pM·-T25=63-Q.57A·M10−8 ut * 10−5 ∆ν·T=0.36·10−5 10−5 ∆ν·T=0.18·10−5 utut ** 10−5 ∆ν·T=0.04·10−5 utut ** utut 10−612 14 16 18 20 22 24ut *** 26 10−618 20 22 24 26 28 30 ut *32 10−624 26 28 30 32 34 36 ut *38 Es/N0 Es/N0 Es/N0 100 PM-16-QAM 100 PM-64-QAM 100 PM-256-QAM ∆p·T=3.57·10−8 ∆p·T=3.57·10−8 ∆p·T=3.57·10−8 ∆ν·T=0.36·10−5 ∆ν·T=0.18·10−5 ∆ν·T=0.04·10−5 10−1 Es/N0=21.28dB 10−1 Es/N0=27.48dB 10−1 Es/N0=33.55dB QPSK16-QAM64-QAM R E S 10−2 10−2 10−2 QPSK 16-QAM µφ=(27), µααα=(28) µφ=(27), µααα=(28) µφ=(27), µααα=(28) 10−3 µφ=µααα=0.1/Es 10−3 µφ=µααα=0.1/Es 10−3 µφ=µααα=0.1/Es 0 1000 2000 3000 4000 5000 0 1000 2000 3000 4000 5000 0 1000 2000 3000 4000 5000 Symbolindexk Symbolindexk Symbolindexk Figure3. Theachievable performanceofthethreetrackingschemesforPM-16-QAM,PM-64-QAM,andPM-256-QAMisshown.Eachcolumncorresponds toamodulation format, whereas the rows present different performance metrics. Thepolarization-noise tolerance is showninthe firstrow byplotting SER versus ∆p·T. Thepenalty compared to the AWGNcurve atlow ∆p·T is due tothe applied phase noise. Thetolerance to phase noise is plotted onthe secondrow,where∆ν·T isvaried, whereasthenoisesensitivity isshowninthethirdrowbyvaryingtheSNR.Theroundmarkersshowninthefirstthree rows correspond to thesame channel conditions, i.e., the same∆p·T, ∆ν·T,and SNR.Theconvergence rate is compared onthe fourth row, where the SERisplottedversusthesymbolindexk. be improved if the block size/tracking step/number of test seen in Table I, the proposed algorithm performs the best in phaseswouldbeoptimizedforeachset ofsystem parameters. all scenarios, being able to tolerate up to nine times more However, this is outside the scope of the paper. The perfor- polarization noise. mance of the MMA-based algorithms degrades significantly as the size of the constellation increases, providing the worst C. Phase-Noise Tolerance performance for PM-256-QAM. Another important feature of the considered tracking algo- Table I summarizes the maximum-tolerable polarization rithms is the tolerance to phase noise. The phase-noise sensi- linewidthtimesthesymboltimesuchthatthereceiverrequires tivity is evaluated by varying ∆ν and keeping the rest of the an extra 1 dB SNR to achieve SER = 10 3 compared parameters fixed. The SNR is set such that a SER=10 3 is − − to the case without laser and polarization noise. As can be achievedinan AWGN scenarioandthe polarizationlinewidth 8 issetsuchthat∆p T =3.57 10 8(correspondingto∆p=1 stepsizesgivefasterconvergencebutlessaccuratetrackingin − · · kHz at 28 Gbaud). steady state. Fig. 3, second row, showsresults obtainedby the evaluated Fig. 3, bottom row, shows that the proposed algorithm algorithms, where the SER is plotted versus the accumulated provides the fastest convergence compared to the other two laserlinewidthtimesthesymboltime.Theproposedalgorithm algorithms and it approachesthe AWGN scenario the closest. has similar performance compared to BPS, which is known In the consideredscenarios, the convergencetakesup to 2500 to have one of the best phase-noise sensitivities for QAM symbols, depending on the constellation, corresponding to modulation formats [31]. 100 ps. ∼ Table I summarizes the maximum tolerable laser linewidth The convergence rate of the three algorithms can be im- times the symbol time for an extra 1 dB SNR to achieve proved by parallelization [33], where the computations are SER=10 3 comparedtothecasewithoutlaserandpolariza- performedforanumberofdifferentinitialmatricesinparallel − tion noise. As can be seen in Table I, the proposed algorithm and the best candidate is chosen at a later stage. performs the best in all scenarios, except for PM-16-QAM and PM-64-QAM where it has a slightly worse phase-noise V. HARDWARE COMPLEXITY sensitivity, being able to tolerate laser linewidths up to 32 The assessment of the hardware complexity of the three MHz. considered algorithms is done by comparing the number of required real3 operations (additions and multiplications), comparisons, and memory units. Although this approach of D. Additive-Noise Sensitivity quantifyingcomplexitydoes not accurately measure the algo- Fig. 3, third row, shows the SER as a function of SNR for rithm’sefficiencyandisaroughapproximation,itisastarting fixed ∆p and ∆ν shown in the figure. The Kabsch algorithm point. has the biggest penalty compared to the AWGN scenario, up Forpedagogicalreasons,in(18),(29),and(36),theestimate 1.2 dB for PM-256-QAM, while the proposed algorithm has of the channel matrix is computed (instead of the estimate the smallest penalty for all three modulation formats with of the inverse), which is then inverted and used throughout a maximum of 0.4 dB for PM-256-QAM. The performance the rest of the algorithm. The inversion step is unnecessary of the MMA+BPS combination is in between the two other and can be skipped by estimating the inverse of the channel algorithms, but is the worst at low SNR. matrixalreadyin(18),(29),and(36)bynegatingφˆ ,αααˆ ,thus k k reducing the required number of operations. In addition, cal- culationsthatoccursmultipletimesinthedecodingprocedure, E. Convergence Rate e.g.,Hˆ 1r in (17), (22), (24), are computedonly once.This An important quality of non-data aided algorithms is the −k k holds for all three algorithms. convergence rate based on blind data. To compare the con- In Table I, the obtained results corresponding to one pro- vergence rates of the three algorithms, additional numerical cessed 4D symbol are listed. In this evaluation, the required simulationshave been performed,in which the initial channel number of operations was minimized by using memory units matrix H was not known by the algorithms and generated 0 instead. This choice highly affects the BPS algorithm, where by the method described in Section II. Note that the channel instead of recomputingpreviouslydone calculationsfor every matrixwasdynamicandupdatedaccordingto(14).Inorderto processed symbol, the values are stored in the memory. compare the algorithms in a fair manner, the initial estimates Therefore, the results are different from the ones presented ofthechannelmatrixweresettobetheidentity,i.e.,Hˆ =I 0 2 in [34], where no memory was used. or I4 for the Kabsch algorithm. As can be seen in Table I, the Kabsch algorithm4 has The convergencerates for the differentmodulation formats the lowest complexity since it operates in a block-wise are shown in Fig. 3, bottom row, where the SER evolution is fashion processing the entire block at a time. This reduces plotted versus the symbol index k. The results were obtained the complexity significantly, which depends on the block byaveragingover4 105realizations,whereineachrealization length, i.e., a longer block will require less computing power · theinitialchannelmatrixwasrandomlygenerated.Toimprove per symbol, but also degrades the performance since it as- the convergence rate of the MMA algorithm, the training sumes that the channel does not change during the block. regime has been split into several stages with increasing The (PS-)CMA/MMA+BPS combination requires the highest numbers of target radii from QPSK to the final constellation computational effort and it increases with the constellation [32],asisshowninthefigures.TheperformanceoftheKabsch size and the size of the sliding-window over which the BPS algorithm could be improved by splitting the training regime algorithmoperates.Thecomplexityoftheproposedalgorithm into different stages with different window lengths NKab, but was evaluated using the Jones description, which is the same it is outside the scope of this work. The convergence rate of with respect to different modulation formats and it is in the proposed algorithm is improved by dividing the training between the other two algorithms. The overall complexity regimeintotwostages,eachwithdifferenttrackingparameters µ ,µ . The first stage, k =1,...,2000, uses the same value 3Thecomplexmultiplications havebeenconverted torealoperations such φ ααα for both step sizes µ = µ = 0.1/E , which is bigger than that acomplex multiplication requires fourreal multiplications and tworeal φ ααα s summations,thussixoperations intotal. theonesusedinthesecondstage,k >2000,wherethevalues 4WeconsideredtheGolub–Reinsch method[35,p.493]forsingularvalue were obtained using (39) and (40). It was observed that high decomposition required bythealgorithm. 9 Algorithm 1: Proposed algorithm The description shown in Algorithm 1 uses the Jones Input: r , G , k formalism. This can be interchanged with the Stokes or real k k Output: uˆk, Gk+1 4D formalisms by replacing the Jones matrix Gk with Mˆ −k1 1 G0 =I2 //initialize the channel matrix (Section III-B) or Rˆ−k1 (Section III-C), and φˆk, αααˆk should 2 uˆ =argmin G r c 2 //decide the be calculated using (30) or (37)–(38), respectively. In this k c∈C k k k− k symbol case, the Jones vectors rk, uˆk become Stokes vectors or 3 if k 2000 then //set the tracking steps real 4D vectors. The Jones description of the algorithm is ≤ 4 µ =µ =0.1/E //convergence stage fully equivalentwith the real 4D descriptionand they achieve φ ααα s the same performance. However, the latter involves more 5 else computationssincemultiplyingtwo4 4realmatricesrequires 6 µ = √∆νTc //tracking stage × φ Es moreoperationsthanmultiplyingtwo2 2complexmatrices. 7 µααα = √∆EpsTc Nevertheless, the real 4D description×can be modified to account for hardware imperfections. On the other hand, the 8 φˆ = 2µ Re i(G r uˆ )HG r k − φ k k− k k k Stokesdescriptionisimmunetoabsolutephaseshiftsandonly 9 αααˆ = 2µ Re(cid:16)i(G r uˆ )HG σ~σσr(cid:17) trackstheSOP,requiringaseparatesolutionforabsolutephase k − ααα k k− k k k 10 G =G H(cid:16)( φˆ , αααˆ ) (cid:17)//update the tracking. This can be beneficial in some situations since fast k+1 k k k − − oscillations of the phase noise will not affect the algorithm, channel matrix andthereforethe SOP tracking.Performance-wise,the Stokes description will behave differently from the other two as it relies on a suboptimal error function (28). of the algorithm can be reduced if the SOP tracking is not updated at every symbol in (18), as discussed in Section VII. DISCUSSION AND CONCLUSIONS III-A. This is a reasonable approach since the SOP does not drift as fast as the phase noise and it can be considered We have proposeda model-basedalgorithmto jointly track constant over a limited observation time in some scenarios. random polarization and phase drifts. The algorithm uses In this case, if αααˆ is computed every P symbols and φˆ the gradient descent optimization algorithm in a decision- k k (which compensates for the phase noise) is updated at every directed architecture processing one symbol at a time. The symbol, the number of required operations is reduced from achievable performance of the algorithm was evaluated by 346 (Table I) to 203 + 143/P. Of course in this scenario meansof numericalsimulationsandcomparedto state-of-the- the performance will degrade, thus resulting in less tolerable artalgorithms.Results show the effectivenessof the proposed polarization noise, but still sufficient for most installed fiber algorithm, having a fast convergence rate and an excellent links. tolerancetophasenoiseanddynamicdriftsofthepolarization, in particular when using high-order modulation formats. At similarorbetterperformance,thecomputationalcomplexityof VI. ALGORITHMIC SUMMARY ouralgorithmislowercomparedtostate-of-the-artalgorithms. The increased phase and polarization noise tolerance for This section provides an easily implementable form of the high-order modulation formats at low complexity makes the proposed algorithm without requiring knowledge about the algorithm a strong candidate for future optical systems. details of the derivations. The proposed algorithm is tested in a somewhat idealized The algorithm is summarized in Algorithm 1, where for scenario assuming that chromatic dispersion, polarization- ease of notation we denote Gk = Hˆ−k1 to be the inverted dependentlosses,andpolarizationmodedispersionarenonex- estimate of the channel matrix at time k. The algorithm istent.Ingeneral,thepolarizationdemultiplexingisperformed receivesasinputsthereceivedsymbolr ,thepreviousinverse k using finite impulse response filters in conjunction with the estimate of the channel matrix G , and the symbol index k, k CMA or MMA algorithms, which, in addition to SOP re- and outputs the decided symbol uˆ and the updated matrix k covery,performchannelequalization.This equalization is not G . The step sizes, which should be precomputed,are higher k covered by the proposed algorithm, but it can be overcame during convergence (k 2000) than tracking (k > 2000). ≤ by coupling the algorithm with a separate equalization stage In the tracking stage, they are computed based on the laser as in [36] to compensate for the slowly-varying transfer linewidth∆ν andthepolarizationlinewidth∆p,whichshould characteristics of the channel. beoverestimatediftheyareunknowntothereceiver.ThePauli matrices σ~σσ are given in (10), C is the set of constellation points, E is the average symbol energy (4), T is the symbol APPENDIX s time, and H() is computed according to (15). Inthisappendix,thegradientoftheerrorfunction(19)with · To decrease the computational effort, for slowly-varying respect to φ,ααα is derived as channels,theupdateofGkcanbedonelessfrequentlythanfor 1 2 every received symbol. Moreover, since the SOP drift varies e = H(φ,ααα)Hˆ − r uˆ (41) ∇ααα k ∇ααα k k− k slower than the phase noise, the SOP update of G can be (cid:13) (cid:13) k (cid:13)(cid:16) (cid:17) (cid:13) H donelessoftenbynotcalculatingαααˆk inline9ateachiteration =Re (cid:13)(cid:13)2 H(φ,ααα)Hˆ −1r uˆ(cid:13)(cid:13) but setting it to αααˆk =[0,0,0]T instead. (cid:18)(cid:16) k(cid:17) k− k(cid:19) 10 H(φ,ααα)Hˆ −1r uˆ (42) [10] A. Viterbi and A. Viterbi, “Nonlinear estimation of PSK-modulated ·∇ααα(cid:18)(cid:16) k(cid:17) k− k(cid:19)! cTarrarniesractpihoansseonwiItnhfoarmppaltiicoantioTnhetoory,buvroslt. 2d9ig,itnaol.t4r,anpspm.i5ss4i3o–n5,”51I,EJEuEl. H 1983. =2Re (cid:18)Hˆ−k1H(−φ,−ααα)rk−uˆk(cid:19) Hˆ−k1 [11] Tre.cPeifvaeur,cSo.nHceopfftmwainthn,feaenddfoRr.wNarode´,ca“rHriaerrdrweacroev-eerfyficfioerntMc-oQhAerMentcodnigstietal-l lations,”JournalofLightwaveTechnology, vol.27,no.8,pp.989–999, Apr.2009. αααH( φ, ααα)rk , (43) [12] H. Louchet, K. Kuzmin, and A. Richter, “Joint carrier-phase and ·∇ − − ! polarization rotation recovery forarbitrary signal constellations,” IEEE Photonics Technology Letters, vol.26,no.9,pp.922–924, May2014. where (42) follows because, for any y Cn and x Rm, x y 2 = x(yHy)=2Re(yH xy). T∈hepartialderi∈vatives [13] DO.. EA..AMgaozrzeir,o“,DMes.igAn.trCadaes-toriflflso´na,ndA.chAalgleunigrrees, iMn.prRac.ticHaulecdoah,eraenndt o∇f Hk(kφ, ∇ααα)=eiφJ( ααα) are,∇from (11) and (15), opticaltransceiver implementations,” JournalofLightwaveTechnology, − − − toappear, 2016. ∂H( φ, ααα) [14] J.N.Damask,PolarizationOpticsinTelecommunications. NewYork, − − =ieiφJ( ααα), (44) NY:Springer, 2005. ∂φ − [15] E.AgrellandM.Karlsson,“Power-efficient modulation formats inco- ∂H( φ, ααα) ∂J( ααα) herenttransmissionsystems,”JournalofLightwaveTechnology,vol.27, − − =eiφ − (45) no.22,pp.5115–5126, Nov.2009. ∂α ∂α i i [16] M. Tur, B.Moslehi, and J.Goodman, “Theory oflaser phase noise in ∂ I cosθ+ia σ~σσsinθ recirculating fiber-optic delay lines,” JournalofLightwave Technology, 2 =eiφ · (46) vol.3,no.1,pp.20–31,Feb.1985. (cid:16) ∂α (cid:17) [17] R.Bellman,IntroductiontoMatrixAnalysis. NewYork,NY:McGraw- i σσσ a Hill,1960. =eiφ I a sinθ+i i ia σ~σσ sinθ [18] J.P.GordonandH.Kogelnik, “PMDfundamentals: polarization mode 2 i − θ − θ · dispersion in optical fibers,” Proceedings of the National Academy of (cid:18) (cid:16) (cid:17) SciencesoftheUnitedStatesofAmerica,vol.97,no.9,pp.4541–4550, +ia σ~σσa cosθ , (47) Apr.2000. i · [19] M. Karlsson, C. B. Czegledi, and E. Agrell, “Coherent transmission (cid:19) channels as 4d rotations,” in Proc. Signal Processing in Photonic where (47) follows because ∂θ/∂α =α / ααα =a . i i i Communication (SPPcom),Boston,MA,Jul.2015,p.SpM3E.2. k k Evaluating(44)and(47)atφ=0,ααα=[0,0,0]T resultsin [20] M. Karlsson, “Four-dimensional rotations in coherent optical commu- nications,” Journal ofLightwave Technology, vol.32, no.6,pp.1246– ∂H( φ, ααα) ∂H( φ, ααα) 1257,Mar.2014. − − =iI , − − =iσσσ , (48) ∂φ 2 ∂α i [21] H. Sun and K.-T. Wu, “A novel dispersion and PMD tolerant clock i phasedetectorforcoherenttransmissionsystems,”inProc.OpticalFiber asθ 0foranydirectiona.Thegradientoftheerrorfunction Communication Conference (OFC), Los Angeles, CA, Mar. 2011, p. canb→eobtainedbysubstituting(48)andH(0,0)=I in(43), OMJ4. 2 [22] T.Nakagawa,K.Ishihara,T.Kobayashi,R.Kudo,M.Matsui,Y.Taka- which is then substituted in (21) and (23) to obtain (22) and tori,andM.Mizoguchi,“Wide-rangeandfast-tracking frequencyoffset (24), respectively. estimator for optical coherent receivers,” in Proc. European Confer- ence on Optical Communication (ECOC), Torino, Italy, Sept. 2010, p. We.7.A.2. REFERENCES [23] S.BoydandL.Vandenberghe, Convex Optimization. NewYork,NY: CambridgeUniversity Press,2004. [1] R. Kim, M. O’Sullivan, Kuang-Tsan Wu, Han Sun, A. Awadalla, D. Krause, and C. Laperle, “Performance of dual-polarization QPSK [24] K.Ogaki,M.Nakada,Y.Nagao,andK.Nishijima,“Fluctuation differ- foropticaltransportsystems,”JournalofLightwaveTechnology,vol.27, encesintheprincipalstatesofpolarization inaerialandburiedcables,” no.16,pp.3546–3559, Aug.2009. inProc.OpticalFiberCommunicationConference(OFC),Atlanta,GA, [2] S.J.Savory,“Digitalcoherentopticalreceivers:algorithmsandsubsys- Mar.2003,p.MF13. tems,”IEEEJournalofSelectedTopicsinQuantumElectronics,vol.16, [25] P. Krummrich, E.-D. Schmidt, W. Weiershausen, and A. Mattheus, no.5,pp.1164–1179, Sept.–Oct.2010. “Fieldtrialresultsonstatisticsoffastpolarization changesinlonghaul [3] C.B.Czegledi,M.Karlsson,E.Agrell,andP.Johannisson,“Polarization WDM transmission systems,” in Proc. Optical Fiber Communication drift channel modelforcoherent fibre-optic systems,”NatureScientific Conference (OFC),Anaheim,CA,Mar.2005,p.OThT6. Reports,toappear, 2016. [26] C. B. Czegledi, E. Agrell, and M. Karlsson, “Symbol-by-symbol joint [4] D. N. Godard, “Self-recovering equalization and carrier tracking in polarization andphasetrackingincoherentreceivers,” inProc.Optical two-dimensional data communication systems,” IEEE Transactions on FiberCommunicationConference(OFC),LosAngeles,CA,Mar.2015, Communications, vol.28,no.11,pp.1867–1875,Nov.1980. p.W1E.3. [5] J. Yang, J.-J. Werner, and G. A. Dumont, “The multimodulus blind [27] S.Betti, F.Curti, G.DeMarchis, andE.Iannone, “A novelmultilevel equalization anditsgeneralized algorithms,” IEEEJournal onSelected coherent optical system:4-quadrature signaling,” Journal ofLightwave AreasinCommunications, vol.20,no.5,pp.997–1015, Jun.2002. Technology, vol.9,no.4,pp.514–523,Apr.1991. [6] H. Louchet, K. Kuzmin, and A. Richter, “Improved DSP algorithms [28] R. Cusani, E. Iannone, A. M. Salonico, and M. Todaro, “An efficient for coherent 16-QAM transmission,” in Proc. European Conference multilevelcoherentopticalsystem:M-4Q-QAM,”JournalofLightwave onOptical Communication (ECOC),Brussels,Belgium, Sept. 2008,p. Technology, vol.10,no.6,pp.777–786, Jun.1992. Tu.1.E.6. [29] D.E.Crivelli, M.R.Hueda, H.S.Carrer, M.delBarco, R.R.Lo´pez, [7] P.Johannisson,M.Sjo¨din,M.Karlsson,H.Wymeersch,E.Agrell,and P. Gianni, J. Finochietto, N. Swenson, P. Voois, and O. E. Agazzi, P.A.Andrekson,“Modifiedconstantmodulusalgorithmforpolarization- “Architectureofasingle-chip50Gb/sDP-QPSK/BPSKtransceiverwith switched QPSK,” Optics Express, vol. 19, no. 8, pp. 7734–7741, Apr. electronicdispersioncompensationforcoherentopticalchannels,”IEEE 2011. TransactionsonCircuitsandSystemsI:RegularPapers,vol.61,no.4, [8] M.Visintin,G.Bosco,P.Poggiolini,andF.Forghieri,“Adaptivedigital pp.1012–1025, Apr.2014. equalization in optical coherent receivers with Stokes-space update [30] M. Seimetz, High-Order Modulation for Optical Fiber Transmission. algorithm,”JournalofLightwaveTechnology,vol.32,no.24,pp.4759– Heidelberg, Germany:Springer, 2009. 4767,Dec.2014. [31] X. Zhou, “Efficient clock and carrier recovery algorithms for single- [9] N. J. Muga and A. N. Pinto, “Adaptive 3-D Stokes space-based po- carriercoherentopticalsystems:asystematicreviewonchallenges and larization demultiplexing algorithm,”JournalofLightwaveTechnology, recentprogress,”IEEESignalProcessingMagazine, vol.31,no.2,pp. vol.32,no.19,pp.3290–3298,Oct.2014. 35–45,Mar.2014.