Public Member Functions | Private Member Functions | Private Attributes

ReduOptFiber Class Reference

The reduced optical fiber class. More...

#include <ocsReduOptFiber.hh>

List of all members.

Public Member Functions

 ReduOptFiber (string InFileName, ReduOptSignal *oReduOptSignalTmp, RanNumGen *oRanNumGenTmp)
 Class constructor.
 ReduOptFiber (const ReduOptFiber &oReduOptFiber)
 Copy constructor.
 ~ReduOptFiber ()
 Destructor.
void CopyReduOptFiber (const ReduOptFiber &oReduOptFiber)
 Called by the copy constructor.
void GetFiberSampleCompleteScatt (void)
 Gets a random realization of a fiber with PMD.
void SetScatteringMatrix (int IndexStep, double Theta, double Phi, double Psi)
 Sets the Euler angles of the fiber section indexed by IndexStep.
void SetLengthFiberAndQtStepsFiber (double Length, int numSteps)
void SetStepLengths (double *usrStepLength)
 Allows user to set the lengths of the fiber sections.
void SetAttenuation (double AttenuationFiber2)
double GetAttenuationFiber ()
double GetLengthFiber ()
void GetRandomStepLengths (double RelativeRangeLengths)
 Used to randomly assign the lengths of the fiber sections.
void SyncNumChannels (void)
void PropagateFiberCoupledNLS (double *PropagatedLength)
 Propagates the reduced dignal through the fiber.
void PropagateFiberCoupledNLS (void)
void WriteFileSteps (char *outFile)
 Writes the lengths of the fiber sections to a file.
void SetMeanDGD_SqrtLength (double MeanDGD_SqrtLength2)
void SetPMD_Fiber (double PMD_Fiber2)
 Sets the PMD strength as given by PMD_Fiber.
double RMS2meanDGD (void)
 Conversion factor from PMD_Fiber to MeanDGD_SqrtLength.
double GetDGD (double Freq)
double GetSndOrderDGD (double Freq)
 Returns the length of the second order polaroization dispersion vector.
void GetInPD_Vector (double Freq, double *InPD_Vector)
void GetOutPD_Vector (double Freq, double *OutPD_Vector)
void GetSndOrderOutPD_Vector (double Freq, double *SndOrderOutPD_Vector)
double GetExpectedDGD ()
double GetExpectedSndOrderDGD ()
double GetMeanDGD_SqrtLength (void)
 Returns the Mean DGD per square root length in ps/sqrt(km).

Private Member Functions

void AllocateMemory_SetLengthSteps (void)
void ReleaseMemory ()
 Releases allocated memory.

Private Attributes

ReduOptSignaloReduOptSignal
RanNumGenoRanNumGen
int NumChannels
 Number of channels of oReduOptSignal.
int GotFiberSample
 Flag that is set to 1 in call to GetFiberSampleCompleteScatt().
double LengthFiber
 The length of the fiber in meters.
int qtStepsFiber
 The number of sections in the fiber.
double PMD_Fiber
 The PMD of the fiber in s/sqrt(m).
double MeanDGD_SqrtLength
 The mean DGD per square root length of the fiber in s/sqrt(m).
double AttenuationFiber
 The fiber loss in dB/meter.
double NepperAttenuationFiber
 The fiber loss expressed on a log scale in base e = 2.718...
double CenterFreq
 The central frequency of the reduced optical signal.
double * psi
 Arrays for the Euler angles theta,phi,psi for each random rotation.
double * phi
double * theta
double * StepLength
 Array to store the lengths of the fiber sections.
double cvtheta2
double cvtheta
double DeltaZ
double dBeta_2dOmega
 The birefringence strength.
int TypeScattering
 The type of random rotation.

Detailed Description

The reduced optical fiber class.

This class is used to propagate the reduced optical signal through an optical fiber with PMD using the coarse step method applied to the Stokes parameters of the signal and noise in each channel. It also imparts a (lumped) fiber loss to the signal. The loss was included to make it simpler to compare the reduced model to a full time and frequency domain model. Note that chromatic dispersion and thr fiber's Kerr nonlinearity are not included in the reduced optical fiber.

Details of the algorithm can be found in Wang and Menyuk, JLT Vol 19, Num 4, April 2001 [WM 2001], page 488. The reader unfamiliar with the coarse step method is also encouraged to read Ref 8 of [WM 2001] and/or Curtis Menyuk's paper in the J. Engineering Math. Vol 36 pp 113-136, 1999. [JEM99]. The paper "Comparison of Polarization Mode Dispersion Emulators" by Ivan Lima et al. may also be helpful as this code is based on the full model code Ivan used to do the research in that paper.


Constructor & Destructor Documentation

ReduOptFiber::ReduOptFiber ( string  InFileName,
ReduOptSignal oReduOptSignalTmp,
RanNumGen oRanNumGenTmp 
)

Class constructor.

Parameters:
string InFileName: The name of the input parameter file.
ReduOptSignal *oReduOptSignalTmp: Pointer to the reduced optical signal object that we are going to propagate through the fiber.
RanNumGen *oRanNumGenTmp: Random number generator object used to obtain the random realizations of the fiber PMD

The important parameters read in from the input parameter file are

  1. LengthFiber: The length of the fiber in meters
  2. qtStepsFiber: The number of sections in the fiber.
  3. PMD_Fiber: The PMD of the fiber in ps/sqrt(km).
  4. MeanDGD_SqrtLength: The mean differential group delay per square root length in ps/sqrt(km).
  5. AttenuationFiber: The fiber loss in dB/meter.

You should set either PMD_Fiber or MeanDGD_SqrtLength. The code will use whichever of these is non-zero in the input parameter file.

The constructor also allocated memory for the Euler angles and computes derived parameters.

See also:
LengthFiber, qtStepsFiber, PMD_Fiber, MeanDGD_SqrtLength, oRanNumGen, oReduOptSignal.

References AttenuationFiber, CenterFreq, DeltaZ, ReduOptSignal::GetCenterFreq(), ReduOptSignal::GetNumChannels(), GotFiberSample, LengthFiber, log(), LogFile, LogFileSeparator(), LOWER_ONLY, MeanDGD_SqrtLength, NepperAttenuationFiber, NumChannels, oRanNumGen, oReduOptSignal, phi, PMD_Fiber, psi, qtStepsFiber, ReadDouble(), ReadInt(), SetMeanDGD_SqrtLength(), SetPMD_Fiber(), sqrt(), StepLength, theta, and TypeScattering.

ReduOptFiber::ReduOptFiber ( const ReduOptFiber oReduOptFiber  ) 
ReduOptFiber::~ReduOptFiber (  ) 

Destructor.

References ReleaseMemory().


Member Function Documentation

void ReduOptFiber::AllocateMemory_SetLengthSteps ( void   )  [private]

Allocates memory to store the Euler angles and the lengths of the fiber sections.

Also initializes dBeta_2dOmega.

See also:
void ReleaseMemory(),
ReduOptFiber(string ,ReduOptSignal*,RanNumGen *)

References dBeta_2dOmega, DeltaZ, GotFiberSample, LengthFiber, phi, PMD_Fiber, psi, qtStepsFiber, sqrt(), StepLength, and theta.

Referenced by SetLengthFiberAndQtStepsFiber().

void ReduOptFiber::CopyReduOptFiber ( const ReduOptFiber oReduOptFiber  ) 
double ReduOptFiber::GetAttenuationFiber (  )  [inline]

References AttenuationFiber.

double ReduOptFiber::GetDGD ( double  Freq  ) 

Returns the differential group delay (DGD) of the current random realization of the fiber at the given frequency.

The DGD is given in units of seconds.

Parameters:
double Freq: The frequency at which the DGD is to be computed. This frequency must be an absolute frequency, i.e. it is typically on the scale of THz. So for instance the DGD at the central frequency of the simulation window is given by GetDGD(oReduOptSignal->GetCenterFreq())

The method computes the polarization dispersion vector at the given frequency at the output of the fiber and returns its length, which is the DGD.

References GetOutPD_Vector(), sq(), and sqrt().

double ReduOptFiber::GetExpectedDGD ( void   )  [inline]

Returns the expected value of the DGD in seconds for a fiber of length LengthFiber and PMD = PMD_Fiber.

References LengthFiber, pi, PMD_Fiber, and sqrt().

Referenced by GetExpectedSndOrderDGD(), and GetSndOrderOutPD_Vector().

double ReduOptFiber::GetExpectedSndOrderDGD ( void   )  [inline]

Returns the expected value of the second DGD in seconds squared for a fiber of length LengthFiber and PMD = PMD_Fiber.

References GetExpectedDGD(), sq(), and sqrt().

void ReduOptFiber::GetFiberSampleCompleteScatt ( void   ) 

Gets a random realization of a fiber with PMD.

This method must be called by the user in an application program. It is called after the fiber constructor and before the propagation method. It should be called each time the user wants a new random realization of the fiber.

In the coarse-step method we divide the fiber into sections each of about 1 km in length. Before each section we randomly rotate the Stokes vectors on the Poincare sphere. This rotation is equivalent to a random rotation of the axes of birefringence of the fiber. The rotations are defined in terms of the three Euler angles theta,phi,psi. The purpose of this method is to randomly select the Euler angles for the rotations.

As explained in [WM 2001] page 488 the Euler angles are randomly chosen so that we get a random rotation chosen from a uniform probability distribution on the Poincare sphere.

The number of fiber sections is given by the class parameter qtStepsFiber. The Euler angles are stored in the arrays theta* phi* and psi*

See also:
qtStepsFiber, theta,phi,psi

References cvtheta, RanNumGen::GetRanNum(), GotFiberSample, oRanNumGen, phi, pi, PMD_Fiber, psi, qtStepsFiber, theta, and TypeScattering.

void ReduOptFiber::GetInPD_Vector ( double  Freq,
double *  InPD_Vector 
)

Computes the polarization dispersion vector at the given frequency at the fiber input

The polarization dispersion vector (PDV) at the input to this fiber is the Stokes vector that when transformed by the random realization of the fiber maps to the output PDV.

See also:
void GetOutPD_Vector((double Freq, double *OutPD_Vector)

References dBeta_2dOmega, GotFiberSample, phi, pi, psi, qtStepsFiber, RotatesAboutX(), RotatesAboutY(), RotatesAboutZ(), StepLength, theta, and TypeScattering.

double ReduOptFiber::GetLengthFiber (  )  [inline]

References LengthFiber.

double ReduOptFiber::GetMeanDGD_SqrtLength ( void   )  [inline]

Returns the Mean DGD per square root length in ps/sqrt(km).

References MeanDGD_SqrtLength, and sqrt().

void ReduOptFiber::GetOutPD_Vector ( double  Freq,
double *  OutPD_Vector 
)

Computes the polarization dispersion vector at the given frequency at the fiber output

Parameters:
Freq,: The absolute frequency at which we wish to compute the polarization dispersion vector (PDV).
OutPD_Vector,: A pointer to an array of length 3 which stores the PDV at output. We assume that memory for this array has been allocated by the user or by another method in this class.

To get the PDV vector for a system (or single fiber) you should set OutPD_Vector = {0,0,0} before calling this function for the first fiber in a system. Then you need to call the function once for each time the signal is passed through the fiber. The answer you get will of course depend on the random realization of the fiber that you choose.

This method is called by

See also:
double GetDGD(double Freq)

The method we use to compute the PDV is described in Section II of the paper "Comparison of Polarization Mode Dispersion Emulators" by Ivan Lima et al. which I think appeared in JLT in about 2001.

References dBeta_2dOmega, GotFiberSample, phi, pi, psi, qtStepsFiber, RotatesAboutX(), RotatesAboutY(), RotatesAboutZ(), StepLength, theta, and TypeScattering.

Referenced by GetDGD(), and GetSndOrderOutPD_Vector().

void ReduOptFiber::GetRandomStepLengths ( double  RelativeRangeLengths  ) 

Used to randomly assign the lengths of the fiber sections.

Currently not used by the reduced model

References DeltaZ, RanNumGen::GetRanNum(), LengthFiber, oRanNumGen, qtStepsFiber, and StepLength.

double ReduOptFiber::GetSndOrderDGD ( double  Freq  ) 

Returns the length of the second order polaroization dispersion vector.

See also:
double GetSndOrderOutPD_Vector(double Freq, double * SndOrderOutPD_Vector)

References GetSndOrderOutPD_Vector(), sq(), and sqrt().

void ReduOptFiber::GetSndOrderOutPD_Vector ( double  Freq,
double *  SndOrderOutPD_Vector 
)

Computes the second order polarization dispersion vector at the given frequency at the fiber output

Obtained by taking a numerical derivative of the output PDV with respect to frequency where the frequency step is given by the inverse of the expected DGD.

See also:
double GetExpectedDGD(void),
void GetOutPD_Vector(double Freq,double *OutPD_Vector)
Parameters:
Freq,: The absolute frequency at which we wish to compute the second order polarization dispersion vector (PDV2) at output.
OutPD_Vector,: A pointer to an array of length 3 which stores the PDV2 at output. We assume that memory for this array has been allocated by the user or by another method in this class.
You should set this vector to zero before obtaining the value for the first fiber

References GetExpectedDGD(), GetOutPD_Vector(), and pi.

Referenced by GetSndOrderDGD().

void ReduOptFiber::PropagateFiberCoupledNLS ( double *  PropagatedLengthFiber  ) 

Propagates the reduced dignal through the fiber.

Parameters:
double *PropagatedLengthFiber: Pointer to a double which keeps track of the length the signal has propagated.

You should call the constructor and void GetFiberSampleCompleteScatt(void) before clling this method.

Uses the coarse step method to model a fiber with PMD with random rotations given by the Euler angles.

Acts on the Stokes vectors of the signal and noise in each channel as in Equations (3,4a,4b) of [WM 2001] page 488.

Also includes fiber loss which is lumped at the end of the fiber. Depending on how the amplification is done you may wish to set the loss parameter to zero.

See also:
void GetFiberSampleCompleteScatt(void)

References dBeta_2dOmega, ReduOptSignal::GetFrequency(), GotFiberSample, LengthFiber, NepperAttenuationFiber, NumChannels, oReduOptSignal, phi, pi, psi, qtStepsFiber, RotatesAboutX(), RotatesAboutY(), RotatesAboutZ(), StepLength, ReduOptSignal::StokesNoise, ReduOptSignal::StokesSignal, theta, and TypeScattering.

void ReduOptFiber::PropagateFiberCoupledNLS ( void   )  [inline]
void ReduOptFiber::ReleaseMemory ( void   )  [private]

Releases allocated memory.

Releases memory allocated by void AllocateMemory_SetLengthSteps(void)
Called by the destructor ~ReduOptFiber()

References phi, psi, StepLength, and theta.

Referenced by SetLengthFiberAndQtStepsFiber(), and ~ReduOptFiber().

double ReduOptFiber::RMS2meanDGD ( void   ) 

Conversion factor from PMD_Fiber to MeanDGD_SqrtLength.

This conversion factor is different in short length (deterministic) regime and in the long length (stochastic) regime.

In the long length regime, we assume DGD is Maxwellian distributed and so the root mean square DGD is related to the mean DGD by a factor of RMS2meanDGD = sqrt(8./(3.*pi)). This conversion factor will be used is qtSteps > 1

Also, if the fiber consists only of one section, qtSteps=1, then we are in the short length regime and PMD is deterministic not stochastic and so it makes sense to define RMS2meanDGD = 1.

See also:
void SetPMD_Fiber(double PMD_Fiber2)
void SetMeanDGD_SqrtLength(double MeanDGD_SqrtLength2) OptFiber::rms2meanDGD(void)

References pi, qtStepsFiber, and sqrt().

Referenced by SetMeanDGD_SqrtLength(), and SetPMD_Fiber().

void ReduOptFiber::SetAttenuation ( double  AttenuationFiber2  ) 
void ReduOptFiber::SetLengthFiberAndQtStepsFiber ( double  Length,
int  numSteps 
)

Method enables user to change the length of fiber and number of sections in the fiber during the course of a simulation.

Sets

  1. LengthFiber = Length
  2. qtStepsFiber = numSteps
  3. Calls ReleaseMemory()
  4. Calls AllocateMemory_SetLengthSteps()
See also:
ReleaseMemory(), AllocateMemory_SetLengthSteps()

References AllocateMemory_SetLengthSteps(), LengthFiber, qtStepsFiber, and ReleaseMemory().

void ReduOptFiber::SetMeanDGD_SqrtLength ( double  MeanDGD_SqrtLength2  ) 

Sets the PMD strength as given by the Mean differential group delay per square root length.

Parameters:
MeanDGD_SqrtLength2,: The Mean differential group delay per square root length in units of ps/sqrt(km)

Initializes MeanDGD_SqrtLength,PMD_Fiber and dBeta_2dOmega. The parameter dBeta_2dOmega is used in void PropagateFiberCoupledNLS(double*)

See also:
void SetPMD_Fiber(double PMD_Fiber2)
void PropagateFiberCoupledNLS(double*)
MeanDGD_SqrtLength,PMD_Fiber, dBeta_2dOmega

References dBeta_2dOmega, DeltaZ, MeanDGD_SqrtLength, PMD_Fiber, RMS2meanDGD(), and sqrt().

Referenced by ReduOptFiber().

void ReduOptFiber::SetPMD_Fiber ( double  PMD_Fiber2  ) 

Sets the PMD strength as given by PMD_Fiber.

Parameters:
PMD_Fiber2,: The PMD of the fiber in ps/sqrt(km)

Initializes PMD_Fiber, MeanDGD_SqrtLength, and dBeta_2dOmega. The parameter dBeta_2dOmega is used in void PropagateFiberCoupledNLS(double*)

See also:
void SetPMD_Fiber(double PMD_Fiber2)
void PropagateFiberCoupledNLS(double*)
MeanDGD_SqrtLength,PMD_Fiber, dBeta_2dOmega

References dBeta_2dOmega, DeltaZ, MeanDGD_SqrtLength, PMD_Fiber, RMS2meanDGD(), and sqrt().

Referenced by ReduOptFiber().

void ReduOptFiber::SetScatteringMatrix ( int  IndexStep,
double  Theta,
double  Phi,
double  Psi 
)

Sets the Euler angles of the fiber section indexed by IndexStep.

Parameters:
IndexStep,: The index of the fiber section. An int in range 0,..,qtSteps-1.
Theta,: The value of the Euler angle theta
Phi,: The value of the Euler angle phi
Psi,: The value of the Euler angle psi.
See also:
void GetFiberSampleCompleteScatt(void)

References GotFiberSample, phi, psi, qtStepsFiber, theta, and TypeScattering.

void ReduOptFiber::SetStepLengths ( double *  usrStepLength  ) 

Allows user to set the lengths of the fiber sections.

Parameters:
double *usrStepLength: A user-defined double array of size qtSteps that contains the values of the lengths of the fiber sections to be set.

The method also computes the sum of the new section lengths and sets FiberLength to this sum.

References LengthFiber, qtStepsFiber, and StepLength.

void ReduOptFiber::SyncNumChannels ( void   ) 

Used to set the class paramter NumChannels to be the same as oReduOptSignal->NumChannels.

References ReduOptSignal::GetNumChannels(), NumChannels, and oReduOptSignal.

void ReduOptFiber::WriteFileSteps ( char *  outFile  ) 

Writes the lengths of the fiber sections to a file.

Parameters:
outFile,: the name of the output file (char*)

The output is given in two columns. Column 1 is the integer index of the fiber section, starting from 0. Column 2 is the length of the section as stored in StepLength[ii].

References qtStepsFiber, and StepLength.


Member Data Documentation

The fiber loss in dB/meter.

See also:
NepperAttenuationFiber

Referenced by CopyReduOptFiber(), GetAttenuationFiber(), ReduOptFiber(), and SetAttenuation().

double ReduOptFiber::CenterFreq [private]

The central frequency of the reduced optical signal.

Referenced by CopyReduOptFiber(), and ReduOptFiber().

double ReduOptFiber::cvtheta [private]
double ReduOptFiber::cvtheta2 [private]
double ReduOptFiber::dBeta_2dOmega [private]

The birefringence strength.

Used in PropagateFiberCoupledNLS(double *).
Obtained from PMD_Fiber by dBeta_2dOmega = PMD_Fiber/sqrt(4*DeltaZ) where DeltaZ is the length of each the fiber section

See also:
PMD_Fiber,
PropagateFiberCoupledNLS(double *)

Referenced by AllocateMemory_SetLengthSteps(), GetInPD_Vector(), GetOutPD_Vector(), PropagateFiberCoupledNLS(), SetMeanDGD_SqrtLength(), and SetPMD_Fiber().

double ReduOptFiber::DeltaZ [private]

The length of the fiber section in meters when the fiber sections all have the same length

Referenced by AllocateMemory_SetLengthSteps(), GetRandomStepLengths(), ReduOptFiber(), SetMeanDGD_SqrtLength(), and SetPMD_Fiber().

double ReduOptFiber::LengthFiber [private]

The mean DGD per square root length of the fiber in s/sqrt(m).

Related to PMD_Fiber by MeanDGD_SqrtLength = sqrt(8./(3.*pi))*PMD_Fiber

See also:
PMD_Fiber, dBeta_2dOmega

Referenced by GetMeanDGD_SqrtLength(), ReduOptFiber(), SetMeanDGD_SqrtLength(), and SetPMD_Fiber().

The fiber loss expressed on a log scale in base e = 2.718...

NepperAttenuationFiber = AttenuationFiber*log(10.0)/10.0;

See also:
AttenuationFiber

Referenced by PropagateFiberCoupledNLS(), ReduOptFiber(), and SetAttenuation().

Number of channels of oReduOptSignal.

Referenced by CopyReduOptFiber(), PropagateFiberCoupledNLS(), ReduOptFiber(), and SyncNumChannels().

Pointer to a random number generator object which is used to obtain the random realization of the fiber in the method GetFiberSampleCompleteScatt()

Referenced by CopyReduOptFiber(), GetFiberSampleCompleteScatt(), GetRandomStepLengths(), and ReduOptFiber().

Pointer to the reduced optical signal object which is to be propagated through the reduced fiber.

Referenced by CopyReduOptFiber(), PropagateFiberCoupledNLS(), ReduOptFiber(), and SyncNumChannels().

double * ReduOptFiber::phi [private]
double ReduOptFiber::PMD_Fiber [private]

The PMD of the fiber in s/sqrt(m).

Related to the mean DGD per squre root length by PMD_Fiber = sqrt(3.*pi/8.)*MeanDGD_SqrtLength

See also:
MeanDGD_SqrtLength, dBeta_2dOmega

Referenced by AllocateMemory_SetLengthSteps(), CopyReduOptFiber(), GetExpectedDGD(), GetFiberSampleCompleteScatt(), ReduOptFiber(), SetMeanDGD_SqrtLength(), and SetPMD_Fiber().

double* ReduOptFiber::psi [private]

The number of sections in the fiber.

In the corase step method for simulating propagation in a fiber with PMD we randomly rotate the axes of birefringence qtStepsFiber times.

Referenced by AllocateMemory_SetLengthSteps(), CopyReduOptFiber(), GetFiberSampleCompleteScatt(), GetInPD_Vector(), GetOutPD_Vector(), GetRandomStepLengths(), PropagateFiberCoupledNLS(), ReduOptFiber(), RMS2meanDGD(), SetLengthFiberAndQtStepsFiber(), SetScatteringMatrix(), SetStepLengths(), and WriteFileSteps().

double* ReduOptFiber::StepLength [private]

Array to store the lengths of the fiber sections.

In this implementation of the reduced model we require that all sections in the same piece of fiber have the same length!!!

Referenced by AllocateMemory_SetLengthSteps(), GetInPD_Vector(), GetOutPD_Vector(), GetRandomStepLengths(), PropagateFiberCoupledNLS(), ReduOptFiber(), ReleaseMemory(), SetStepLengths(), and WriteFileSteps().

double * ReduOptFiber::theta [private]

The type of random rotation.

In the reduced model this is set to 1 to indicate that we use the Euler angles

Referenced by GetFiberSampleCompleteScatt(), GetInPD_Vector(), GetOutPD_Vector(), PropagateFiberCoupledNLS(), ReduOptFiber(), and SetScatteringMatrix().


The documentation for this class was generated from the following files: