anchor.diff:          SPLti = 0.
           RETURN
       ENDIF
  -!! Nafnoise source code version see belo 
  -!     Frequency_cutoff = 10*U/PI/Chord
  -!     Ke = 3.0/(4.0*LTurb) 
  -!     Beta2 = 1-Mach*Mach
  -!     ALPSTAR = AlphaNoise*PI/180.
  -!     DO I=1,size(p%FreqList)
  -!         IF (p%FreqList(I) <= Frequency_cutoff) THEN
  -!            Directivity = DBARL
  -!         ELSE
  -!            Directivity = DBARH 
  -!         ENDIF
  -!         WaveNumber = 2.0*PI*p%FreqList(I)/U
  -!         Kbar = WaveNumber*Chord/2.0
  -!         Khat = WaveNumber/Ke
  -!         mu = Mach*WaveNumber*Chord/2.0/Beta2
  -!         SPLhigh = 10.*LOG10(p%AirDens*p%AirDens*p%SpdSound**4*LTurb*(d/2.)/ &
  -!                  (RObs*RObs)*(Mach**5)*tinooisess*tinooisess*(Khat**3)* &
  -!                  (1+Khat**2)**(-7./3.)*Directivity) + 78.4
  -!!   SPLhigh = 10.*LOG10(LTurb*(d/2.)/ &
  -!!                  (RObs*RObs)*(Mach**5)*tinooisess*tinooisess*(WaveNumber**3) &
  -!!                  *(1+WaveNumber**2)**(-7./3.)*Directivity) + 181.3
  -!         SPLhigh = SPLhigh + 10.*LOG10(1+ 9.0*ALPSTAR*ALPSTAR)
  -!         Sears = 1/(2.*PI*Kbar/Beta2+1/(1+2.4*Kbar/Beta2))
  -!!   Sears = 1/(2.*PI*WaveNumber/Beta2+1/(1+2.4*WaveNumber/Beta2))
  -!         LFC = 10*Sears*Mach*Kbar*Kbar/Beta2
  -!!   LFC = 10*Sears*Mach*WaveNumber*WaveNumber/Beta2
  -!!         IF (mu<(PI/4.0)) THEN
  -!!      SPLti(I) = SPLhigh + 10.*ALOG10(LFC)
  -!!   ELSE
  -!!      SPLti(I) = SPLhigh
  -!!   ENDIF
  -!   SPLti(I) = SPLhigh + 10.*LOG10(LFC/(1+LFC))
  -!      ENDDO
  -!!Wei Jun Zhu et al - !Modeling of Aerodynamically Generated Noise From Wind Turbines 2005 paper
  -     Beta2 = 1.d0-Mach**2; !  corresponding line:  Bsq = 1.d0 - Ma**2;
  +	 
  +	! In the following lines, bibliography will be referenced as:  a) Moriarty, Guidati, Migliore, Recent Improvement of a Semi-Empirical Aeroacoustic 
  +	! Prediction Code for Wind Turbines
  +	! ref b) Lowson, Assessment and Prediction of Wind Turbine Noise
  +	
  +   !*********************************************** Model 1:
  +   !!! Nafnoise source code version see below 
  +       Frequency_cutoff = 10*U/PI/Chord
  +       Ke = 3.0/(4.0*LTurb) 
  +       Beta2 = 1-Mach*Mach
  +       ALPSTAR = AlphaNoise*PI/180.
  +   	
          DO I=1,size(p%FreqList)
  -     WaveNumber = PI*p%FreqList(I)*p%SpdSound/U !corresponding line: K = pi*Freq(i)*c/Vrel;
  -     Sears = (2.d0*PI*WaveNumber/Beta2 + (1.d0+2.4d0*WaveNumber/Beta2)**(-1))**(-1);
  -   ! corresponding line: Ssq = (2.d0*pi*K/Bsq + (1.d0+2.4d0*K/Bsq)**(-1))**(-1);
  -      LFC = 10.d0 * Sears*Mach*WaveNumber**2*Beta2**(-1);
  -    ! corresponding line:  LFC = 10.d0 * Ssq*Ma*K**2*Bsq**(-1);
  -   SPLti(I)=(p%AirDens*p%AirDens*p%SpdSound*p%SpdSound*Lturb*d)/(2*RObs*RObs)
  -!   SPLti(I)=SPLti(I)*(Mach**3)*(MeanVnoise**2)*(tinooisess**2)
  -   SPLti(I)=SPLti(I)*(Mach**3)*(tinooisess**2)
  -!   SPLti(I)=SPLti(I)*(Mach**3)*ufluct**2
  -   SPLti(I)=(SPLti(I)*(WaveNumber**3)) / ((1+WaveNumber**2)**(7/3))
  -   SPLti(I)=SPLti(I)*DBARH
  -   SPLti(I)=10*log10(SPLti(I))+58.4
  -    SPLti(I) = SPLti(I) + 10.*LOG10(LFC/(1+LFC))
  -! SPLti(I)=10.d0*log10(DBARH*p%AirDens**2*p%SpdSound**2*Lturb*d/2.0*Mach**3*tinooisess**2* &
  -!WaveNumber**3*(1.d0+WaveNumber**2)**(-7.d0/3.d0)/RObs**2)+58.4d0 + 10.d0*log10(LFC/(1+LFC))
  - ! corresponding line:    SPLti(i)=10.d0*log10(Di_hi_fr*Density**2*co**2*Tbscale*L/2.0*Ma
  -!     & **3*Tbinten**2*K**3*(1.d0+K**2)**(-7.d0/3.d0)/Distance**2)+58.4d0
  -!     &    + 10.d0*log10(LFC/(1+LFC));  
  -             !% ver2.!
  -!    Kh = 8.d0*pi*p%FreqList(i)*Lturb/(3.d0*U);
  -!          SPLti(i) = 10*log10(DBARH*Lturb*0.5*d*Mach**5*tinooisess**2*Kh**3*(1+Kh**2)**(-7/3)/RObs**2) +&
  -!              10*log10(10**18.13) + 10*log10(DBARH*LFC/(1+LFC));   
  -ENDDO
  +           IF (p%FreqList(I) <= Frequency_cutoff) THEN
  +              Directivity = DBARL
  +           ELSE
  +              Directivity = DBARH 
  +           ENDIF
  +   		
  +           WaveNumber = 2.0*PI*p%FreqList(I)/U
  +           Kbar = WaveNumber*Chord/2.0
  +           Khat = WaveNumber/Ke
  +           ! mu = Mach*WaveNumber*Chord/2.0/Beta2
  +   		
  +           SPLhigh = 10.*LOG10(p%AirDens*p%AirDens*p%SpdSound**4*LTurb*(d/2.)/ &
  +                    (RObs*RObs)*(Mach**5)*tinooisess*tinooisess*(Khat**3)* &
  +                    (1+Khat**2)**(-7./3.)*Directivity) + 78.4   ! ref a)
  +   !!!   SPLhigh = 10.*LOG10(LTurb*(d/2.)/ &
  +   !!!                  (RObs*RObs)*(Mach**5)*tinooisess*tinooisess*(WaveNumber**3) &
  +   !!!                  *(1+WaveNumber**2)**(-7./3.)*Directivity) + 181.3  
  +   
  +           SPLhigh = SPLhigh + 10.*LOG10(1+ 9.0*ALPSTAR*ALPSTAR)  ! Component due to angles of attack, ref a)   
  +   		
  +           Sears = 1/(2.*PI*Kbar/Beta2+1/(1+2.4*Kbar/Beta2))      ! ref a)
  +   		
  +   !!!   Sears = 1/(2.*PI*WaveNumber/Beta2+1/(1+2.4*WaveNumber/Beta2))  ! ref b) 
  +   
  +           LFC = 10*Sears*Mach*Kbar*Kbar/Beta2  ! ref a)
  +   !!!   LFC = 10*Sears*Mach*WaveNumber*WaveNumber/Beta2  ! ref b)
  +   
  +   !!!   IF (mu<(PI/4.0)) THEN                     ! ref b)
  +   !!!      SPLti(I) = SPLhigh + 10.*ALOG10(LFC)   ! ref b)
  +   !!!   ELSE                                      ! ref b)
  +   !!!      SPLti(I) = SPLhigh                     ! ref b)
  +   !!!ENDIF
  +     SPLti(I) = SPLhigh + 10.*LOG10(LFC/(1+LFC))
  +   
  +       ENDDO
  +   !!!*********************************************** end of Model 1
  +
  +!  ! ********************************* Model 2:     
  +!  !Wei Jun Zhu et al - !Modeling of Aerodynamically Generated Noise From Wind Turbines 2005 paper
  +!  		Beta2 = 1.d0-Mach**2; !  corresponding line:  Bsq = 1.d0 - Ma**2;
  +!         DO I=1,size(p%FreqList)
  +!       WaveNumber = PI*p%FreqList(I)*p%SpdSound/U !corresponding line: K = pi*Freq(i)*c/Vrel;  ! CarloS: This is a Mistake, c in this case is the Local Chord
  +!       Sears = (2.d0*PI*WaveNumber/Beta2 + (1.d0+2.4d0*WaveNumber/Beta2)**(-1))**(-1);
  +!     ! corresponding line: Ssq = (2.d0*pi*K/Bsq + (1.d0+2.4d0*K/Bsq)**(-1))**(-1);
  +!        LFC = 10.d0 * Sears*Mach*WaveNumber**2*Beta2**(-1);
  +!      ! corresponding line:  LFC = 10.d0 * Ssq*Ma*K**2*Bsq**(-1);
  +!     SPLti(I)=(p%AirDens*p%AirDens*p%SpdSound*p%SpdSound*Lturb*d)/(2*RObs*RObs)
  +!  !   SPLti(I)=SPLti(I)*(Mach**3)*(MeanVnoise**2)*(tinooisess**2)  
  +!     SPLti(I)=SPLti(I)*(Mach**3)*(tinooisess**2)       			 
  +!  !   SPLti(I)=SPLti(I)*(Mach**3)*ufluct**2
  +!     SPLti(I)=(SPLti(I)*(WaveNumber**3)) / ((1+WaveNumber**2)**(7/3))
  +!     SPLti(I)=SPLti(I)*DBARH
  +!     SPLti(I)=10*log10(SPLti(I))+58.4
  +!      SPLti(I) = SPLti(I) + 10.*LOG10(LFC/(1+LFC))
  +!  ! SPLti(I)=10.d0*log10(DBARH*p%AirDens**2*p%SpdSound**2*Lturb*d/2.0*Mach**3*tinooisess**2* &
  +!  !WaveNumber**3*(1.d0+WaveNumber**2)**(-7.d0/3.d0)/RObs**2)+58.4d0 + 10.d0*log10(LFC/(1+LFC))
  +!  ! corresponding line:    SPLti(i)=10.d0*log10(Di_hi_fr*Density**2*co**2*Tbscale*L/2.0*Ma
  +!  !     & **3*Tbinten**2*K**3*(1.d0+K**2)**(-7.d0/3.d0)/Distance**2)+58.4d0
  +!  !     &    + 10.d0*log10(LFC/(1+LFC));  
  +!  !            !% ver2.!
  +!  !    Kh = 8.d0*pi*p%FreqList(i)*Lturb/(3.d0*U);
  +!  !          SPLti(i) = 10*log10(DBARH*Lturb*0.5*d*Mach**5*tinooisess**2*Kh**3*(1+Kh**2)**(-7/3)/RObs**2) +&
  +!  !              10*log10(10**18.13) + 10*log10(DBARH*LFC/(1+LFC));   
  +!  
  +!  ENDDO
  +!  ! ********************************* End of Model 2/ CarloSucameli: I think this model is wrong 
  +
  +
  +
  +!!!!  ! ********************************* Model 3:     
  +!!!!  ! ref b) Lowson, Assessment and Prediction of Wind Turbine Noise
  +!!!!  	Beta2 = 1.d0-Mach**2; !  corresponding line:  Bsq = 1.d0 - Ma**2;
  +!!!!      DO I=1,size(p%FreqList)
  +!!!!       WaveNumber = PI*p%FreqList(I)*Chord/U !corresponding line: K = pi*Freq(i)*c/Vrel;  
  +!!!!       Sears = (2.d0*PI*WaveNumber/Beta2 + (1.d0+2.4d0*WaveNumber/Beta2)**(-1))**(-1);
  +!!!!     ! corresponding line: Ssq = (2.d0*pi*K/Bsq + (1.d0+2.4d0*K/Bsq)**(-1))**(-1);
  +!!!!        LFC = 10.d0 * Sears*Mach*WaveNumber**2*Beta2**(-1);
  +!!!!      ! corresponding line:  LFC = 10.d0 * Ssq*Ma*K**2*Bsq**(-1);
  +!!!!     SPLti(I)=(p%AirDens*p%AirDens*p%SpdSound*p%SpdSound*Lturb*d)/(2*RObs*RObs)
  +!!!!     SPLti(I)=SPLti(I)*(Mach**3)*(MeanVnoise**2)*(tinooisess**2)  
  +!!!!     SPLti(I)=(SPLti(I)*(WaveNumber**3)) / ((1+WaveNumber**2)**(7./3.))
  +!!!!     SPLti(I)=SPLti(I)*DBARH
  +!!!!     SPLti(I)=10*log10(SPLti(I))+58.4
  +!!!!     SPLti(I) = SPLti(I) + 10.*LOG10(LFC/(1+LFC))
  +!!!!  
  +!!!!
  +!!!!  	ENDDO
  +!!!!  ! ********************************* End of Model 3 
  +  
   !!Buck&Oerlamans&Palo - !Experimental validation of a wind turbine turbulent inflow noise prediction code 2016 paper
   !DO I=1,size(p%FreqList)
           ! IF (p%FreqList(I) <= Frequency_cutoff) THEN
  
positive.code:     ENDDO
  END SUBROUTINE TipNois
  !==================================================================================================================================!
  SUBROUTINE InflowNoise(AlphaNoise,Chord,U,THETA,PHI,d,RObs,MeanVNoise,TINoise,LE_Location,dissip,ufluct,p,SPLti,errStat,errMsg)
    REAL(ReKi),                                 INTENT(IN   ) :: AlphaNoise        ! AOA
    REAL(ReKi),                                 INTENT(IN   ) :: Chord          ! Chord Length
    REAL(ReKi),                                 INTENT(IN   ) :: U              !
    REAL(ReKi),                                 INTENT(IN   ) :: d              ! element span
    REAL(ReKi),                                 INTENT(IN   ) :: RObs           ! distance to observer
    REAL(ReKi),                                 INTENT(IN   ) :: THETA          !
    REAL(ReKi),                                 INTENT(IN   ) :: PHI            ! Spanwise directivity angle
    REAL(ReKi),                                 INTENT(IN   ) :: MeanVNoise     !
    REAL(ReKi),                                 INTENT(IN   ) :: TINoise        !
    REAL(ReKi),                                 INTENT(IN   ) :: LE_Location    !
    REAL(ReKi),                                 INTENT(IN   ) :: dissip         !
    REAL(ReKi),                                 INTENT(IN   ) :: ufluct         !
    TYPE(AA_ParameterType),                     INTENT(IN   ) :: p              ! Parameters
    REAL(ReKi),DIMENSION(size(p%FreqList)),     INTENT(  OUT) :: SPLti         !
    INTEGER(IntKi),                             INTENT(  OUT) :: errStat        ! Error status of the operation
  
negative.code:             ENDIF
              AlphaNoise= u%AoANoise(J,I) * R2D_D 
  
  
              !--------Read in Boundary Layer Data-------------------------!
              IF (p%X_BLMethod .EQ. 2) THEN
                  call BL_Param_Interp(p,m,Unoise,AlphaNoise,p%BlChord(J,I),p%BlAFID(J,I), errStat2, errMsg2)
  