GNU Radio 3.4.0 C++ API
volk_32f_s32f_calc_spectral_noise_floor_32f_a16.h
Go to the documentation of this file.
00001 #ifndef INCLUDED_volk_32f_s32f_calc_spectral_noise_floor_32f_a16_H
00002 #define INCLUDED_volk_32f_s32f_calc_spectral_noise_floor_32f_a16_H
00003 
00004 #include <inttypes.h>
00005 #include <stdio.h>
00006 
00007 #if LV_HAVE_SSE
00008 #include <xmmintrin.h>
00009 /*!
00010   \brief Calculates the spectral noise floor of an input power spectrum
00011 
00012   Calculates the spectral noise floor of an input power spectrum by determining the mean of the input power spectrum, then recalculating the mean excluding any power spectrum values that exceed the mean by the spectralExclusionValue (in dB).  Provides a rough estimation of the signal noise floor.
00013 
00014   \param realDataPoints The input power spectrum
00015   \param num_points The number of data points in the input power spectrum vector
00016   \param spectralExclusionValue The number of dB above the noise floor that a data point must be to be excluded from the noise floor calculation - default value is 20
00017   \param noiseFloorAmplitude The noise floor of the input spectrum, in dB
00018 */
00019 static inline void volk_32f_s32f_calc_spectral_noise_floor_32f_a16_sse(float* noiseFloorAmplitude, const float* realDataPoints, const float spectralExclusionValue, const unsigned int num_points){
00020   unsigned int number = 0;
00021   const unsigned int quarterPoints = num_points / 4;
00022 
00023   const float* dataPointsPtr = realDataPoints;
00024   float avgPointsVector[4] __attribute__((aligned(128)));
00025     
00026   __m128 dataPointsVal;
00027   __m128 avgPointsVal = _mm_setzero_ps();
00028   // Calculate the sum (for mean) for all points
00029   for(; number < quarterPoints; number++){
00030 
00031     dataPointsVal = _mm_load_ps(dataPointsPtr);
00032 
00033     dataPointsPtr += 4;
00034 
00035     avgPointsVal = _mm_add_ps(avgPointsVal, dataPointsVal);
00036   }
00037 
00038   _mm_store_ps(avgPointsVector, avgPointsVal);
00039 
00040   float sumMean = 0.0;
00041   sumMean += avgPointsVector[0];
00042   sumMean += avgPointsVector[1];
00043   sumMean += avgPointsVector[2];
00044   sumMean += avgPointsVector[3];
00045 
00046   number = quarterPoints * 4;
00047   for(;number < num_points; number++){
00048     sumMean += realDataPoints[number];
00049   }
00050 
00051   // calculate the spectral mean
00052   // +20 because for the comparison below we only want to throw out bins
00053   // that are significantly higher (and would, thus, affect the mean more
00054   const float meanAmplitude = (sumMean / ((float)num_points)) + spectralExclusionValue;
00055 
00056   dataPointsPtr = realDataPoints; // Reset the dataPointsPtr
00057   __m128 vMeanAmplitudeVector = _mm_set_ps1(meanAmplitude);
00058   __m128 vOnesVector = _mm_set_ps1(1.0);
00059   __m128 vValidBinCount = _mm_setzero_ps();
00060   avgPointsVal = _mm_setzero_ps();
00061   __m128 compareMask;
00062   number = 0;
00063   // Calculate the sum (for mean) for any points which do NOT exceed the mean amplitude
00064   for(; number < quarterPoints; number++){
00065 
00066     dataPointsVal = _mm_load_ps(dataPointsPtr);
00067 
00068     dataPointsPtr += 4;
00069 
00070     // Identify which items do not exceed the mean amplitude
00071     compareMask = _mm_cmple_ps(dataPointsVal, vMeanAmplitudeVector);
00072 
00073     // Mask off the items that exceed the mean amplitude and add the avg Points that do not exceed the mean amplitude
00074     avgPointsVal = _mm_add_ps(avgPointsVal, _mm_and_ps(compareMask, dataPointsVal));
00075       
00076     // Count the number of bins which do not exceed the mean amplitude
00077     vValidBinCount = _mm_add_ps(vValidBinCount, _mm_and_ps(compareMask, vOnesVector));
00078   }
00079     
00080   // Calculate the mean from the remaining data points
00081   _mm_store_ps(avgPointsVector, avgPointsVal);
00082 
00083   sumMean = 0.0;
00084   sumMean += avgPointsVector[0];
00085   sumMean += avgPointsVector[1];
00086   sumMean += avgPointsVector[2];
00087   sumMean += avgPointsVector[3];
00088 
00089   // Calculate the number of valid bins from the remaning count
00090   float validBinCountVector[4] __attribute__((aligned(128)));
00091   _mm_store_ps(validBinCountVector, vValidBinCount);
00092 
00093   float validBinCount = 0;
00094   validBinCount += validBinCountVector[0];
00095   validBinCount += validBinCountVector[1];
00096   validBinCount += validBinCountVector[2];
00097   validBinCount += validBinCountVector[3];
00098 
00099   number = quarterPoints * 4;
00100   for(;number < num_points; number++){
00101     if(realDataPoints[number] <= meanAmplitude){
00102       sumMean += realDataPoints[number];
00103       validBinCount += 1.0;
00104     }
00105   }
00106     
00107   float localNoiseFloorAmplitude = 0;
00108   if(validBinCount > 0.0){
00109     localNoiseFloorAmplitude = sumMean / validBinCount;
00110   }
00111   else{
00112     localNoiseFloorAmplitude = meanAmplitude; // For the odd case that all the amplitudes are equal...
00113   }
00114 
00115   *noiseFloorAmplitude = localNoiseFloorAmplitude;
00116 }
00117 #endif /* LV_HAVE_SSE */
00118 
00119 #if LV_HAVE_GENERIC
00120 /*!
00121   \brief Calculates the spectral noise floor of an input power spectrum
00122 
00123   Calculates the spectral noise floor of an input power spectrum by determining the mean of the input power spectrum, then recalculating the mean excluding any power spectrum values that exceed the mean by the spectralExclusionValue (in dB).  Provides a rough estimation of the signal noise floor.
00124 
00125   \param realDataPoints The input power spectrum
00126   \param num_points The number of data points in the input power spectrum vector
00127   \param spectralExclusionValue The number of dB above the noise floor that a data point must be to be excluded from the noise floor calculation - default value is 20
00128   \param noiseFloorAmplitude The noise floor of the input spectrum, in dB
00129 */
00130 static inline void volk_32f_s32f_calc_spectral_noise_floor_32f_a16_generic(float* noiseFloorAmplitude, const float* realDataPoints, const float spectralExclusionValue, const unsigned int num_points){
00131   float sumMean = 0.0;
00132   unsigned int number;
00133   // find the sum (for mean), etc
00134   for(number = 0; number < num_points; number++){
00135     // sum (for mean)
00136     sumMean += realDataPoints[number];
00137   }
00138 
00139   // calculate the spectral mean
00140   // +20 because for the comparison below we only want to throw out bins
00141   // that are significantly higher (and would, thus, affect the mean more)
00142   const float meanAmplitude = (sumMean / num_points) + spectralExclusionValue;
00143 
00144   // now throw out any bins higher than the mean
00145   sumMean = 0.0;
00146   unsigned int newNumDataPoints = num_points;
00147   for(number = 0; number < num_points; number++){
00148     if (realDataPoints[number] <= meanAmplitude)
00149       sumMean += realDataPoints[number];
00150     else
00151       newNumDataPoints--;
00152   }
00153 
00154   float localNoiseFloorAmplitude = 0.0;
00155   if (newNumDataPoints == 0)             // in the odd case that all
00156     localNoiseFloorAmplitude = meanAmplitude; // amplitudes are equal!
00157   else
00158     localNoiseFloorAmplitude = sumMean / ((float)newNumDataPoints);
00159 
00160   *noiseFloorAmplitude = localNoiseFloorAmplitude;
00161 }
00162 #endif /* LV_HAVE_GENERIC */
00163 
00164 
00165 
00166 
00167 #endif /* INCLUDED_volk_32f_s32f_calc_spectral_noise_floor_32f_a16_H */