4343from cv2 import filter2D
4444from numba import njit , types as numba_types
4545from matplotlib .path import Path as BasePath
46- from pyproj import Proj
4746from pint import UnitRegistry
4847from ..observations .observation import EddiesObservations
4948from ..eddy_feature import Amplitude , Contours
@@ -701,28 +700,29 @@ def eddy_identification(
701700 for contour in contour_paths :
702701 if contour .used :
703702 continue
704- centlon_e , centlat_e , eddy_radius_e , aerr = contour .fit_circle ()
703+ # FIXME : center could be not in contour and fit on raw sampling
704+ _ , _ , _ , aerr = contour .fit_circle ()
705705
706706 # Filter for shape
707707 if aerr < 0 or aerr > shape_error or isnan (aerr ):
708708 continue
709- # Get indices of centroid
710- # Give only 1D array of lon and lat not 2D data
711- i_x , i_y = self .nearest_grd_indice (centlon_e , centlat_e )
712- i_x = self .normalize_x_indice (i_x )
713-
714- # Check if centroid is on define value
715- if data .mask [i_x , i_y ]:
716- continue
717- # Test to know cyclone or anticyclone
718- acyc_not_cyc = data [i_x , i_y ] >= cvalues
719- if anticyclonic_search != acyc_not_cyc :
720- continue
721709
722710 # Find all pixels in the contour
723711 i_x_in , i_y_in = contour .pixels_in (self )
724712
725- # Maybe limit max must be replace with a maximum of surface
713+ # Check if pixels in contour are masked
714+ if data .mask [i_x_in , i_y_in ].any ():
715+ continue
716+
717+ # Test to know cyclone or anticyclone
718+ if anticyclonic_search :
719+ if (data [i_x_in , i_y_in ] < cvalues ).any ():
720+ continue
721+ else :
722+ if (data [i_x_in , i_y_in ] > cvalues ).any ():
723+ continue
724+
725+ # FIXME : Maybe limit max must be replace with a maximum of surface
726726 if (
727727 contour .nb_pixel < pixel_limit [0 ]
728728 or contour .nb_pixel > pixel_limit [1 ]
@@ -775,30 +775,7 @@ def eddy_identification(
775775 pixel_min = pixel_limit [0 ],
776776 )
777777
778- # Use azimuth equal projection for radius
779- proj = Proj (
780- "+proj=aeqd +ellps=WGS84 +lat_0={1} +lon_0={0}" .format (
781- * inner_contour .mean_coordinates
782- )
783- )
784- # First, get position based on innermost
785- # contour
786- centx_i , centy_i , _ , _ = fit_circle (
787- * proj (inner_contour .lon , inner_contour .lat )
788- )
789- centlon_i , centlat_i = proj (centx_i , centy_i , inverse = True )
790- # Second, get speed-based radius based on
791- # contour of max uavg
792- centx_s , centy_s , eddy_radius_s , aerr_s = fit_circle (
793- * proj (speed_contour .lon , speed_contour .lat )
794- )
795- # Computed again to be coherent with speed_radius, we will be compute in same reference
796- _ , _ , eddy_radius_e , aerr_e = fit_circle (
797- * proj (contour .lon , contour .lat )
798- )
799- centlon_s , centlat_s = proj (centx_s , centy_s , inverse = True )
800-
801- # Instantiate new EddyObservation object (high cost need to be review)
778+ # FIXME : Instantiate new EddyObservation object (high cost need to be review)
802779 obs = EddiesObservations (
803780 size = 1 ,
804781 track_extra_variables = track_extra_variables ,
@@ -818,23 +795,38 @@ def eddy_identification(
818795 else :
819796 obs .obs ["uavg_profile" ] = raw_resample (speed_array , sampling )
820797 obs .obs ["amplitude" ] = amp .amplitude
821- obs .obs ["radius_s" ] = eddy_radius_s
822798 obs .obs ["speed_average" ] = max_average_speed
799+ obs .obs ["num_point_e" ] = contour .lon .shape [0 ]
800+ xy_e = uniform_resample (contour .lon , contour .lat , ** out_sampling )
801+ obs .obs ["contour_lon_e" ], obs .obs ["contour_lat_e" ] = xy_e
802+ obs .obs ["num_point_s" ] = speed_contour .lon .shape [0 ]
803+ xy_s = uniform_resample (
804+ speed_contour .lon , speed_contour .lat , ** out_sampling
805+ )
806+ obs .obs ["contour_lon_s" ], obs .obs ["contour_lat_s" ] = xy_s
807+
808+ # FIXME : we use a contour without resampling
809+ # First, get position based on innermost contour
810+ centlon_i , centlat_i , _ , _ = _fit_circle_path (
811+ create_vertice (inner_contour .lon , inner_contour .lat )
812+ )
813+ # Second, get speed-based radius based on contour of max uavg
814+ centlon_s , centlat_s , eddy_radius_s , aerr_s = _fit_circle_path (
815+ create_vertice (* xy_s )
816+ )
817+ # Computed again to use resample contour
818+ _ , _ , eddy_radius_e , aerr_e = _fit_circle_path (
819+ create_vertice (* xy_e )
820+ )
821+
822+ obs .obs ["radius_s" ] = eddy_radius_s
823823 obs .obs ["radius_e" ] = eddy_radius_e
824824 obs .obs ["shape_error_e" ] = aerr_e
825825 obs .obs ["shape_error_s" ] = aerr_s
826826 obs .obs ["lon" ] = centlon_s
827827 obs .obs ["lat" ] = centlat_s
828828 obs .obs ["lon_max" ] = centlon_i
829829 obs .obs ["lat_max" ] = centlat_i
830- obs .obs ["num_point_e" ] = contour .lon .shape [0 ]
831- xy = uniform_resample (contour .lon , contour .lat , ** out_sampling )
832- obs .obs ["contour_lon_e" ], obs .obs ["contour_lat_e" ] = xy
833- obs .obs ["num_point_s" ] = speed_contour .lon .shape [0 ]
834- xy = uniform_resample (
835- speed_contour .lon , speed_contour .lat , ** out_sampling
836- )
837- obs .obs ["contour_lon_s" ], obs .obs ["contour_lat_s" ] = xy
838830 if aerr > 99.9 or aerr_s > 99.9 :
839831 logger .warning (
840832 "Strange shape at this step! shape_error : %f, %f" ,
0 commit comments