• Publicidad

Interpolación bilineal en Perl

Así que programas sin strict y las expresiones regulares son otro modo de hablar. Aquí encontrarás respuestas de nivel avanzado, no recomendable para los débiles de corazón.

Interpolación bilineal en Perl

Notapor lis » 2008-05-27 21:54 @954

Hola, mi problema es el siguiente: tengo una matriz (con extensión .dat) de 49x50, que en este caso son datos de temperatura, y debo hacer una interpolación bilineal a cualquier punto (x,y), en donde me muestre los 4 puntos que rodean a (x,y) y además me muestre el resultado de la interpolación.

Si alguien me puede ayudar con este script en Perl se lo agradecería mucho porque llevo en esto casi 3 semanas.

Muchas gracias...
lis
Perlero nuevo
Perlero nuevo
 
Mensajes: 106
Registrado: 2008-05-27 21:43 @946

Publicidad

Notapor explorer » 2008-05-28 06:16 @303

¿Quieres trabajar con toda la imagen de golpe o solo un pixel que elijas tu? ¿Cómo quieres elegirlo, con el ratón o con el teclado? ¿Y la salida, cómo debe ser?

Actualización: Ya encontré la expresión matemática de la interpolación bilineal.
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14480
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Notapor explorer » 2008-05-28 08:03 @377

Este es un ejemplo de cálculo de interpolación bilineal:
Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
#!/usr/bin/perl
use warnings;
use strict;

# Toma de datos.
# Suponemos que el origen de coordenadas está en
# la coordenada superior izquierda de la tabla.
my $temperaturas = [
    [ 0, 0, 0, 0, 0, 0, 0, ],
    [ 0, 5, 4, 4, 6, 5, 0, ],
    [ 1, 4, 5, 5, 7, 6, 0, ],
    [ 3, 3, 4, 5, 4, 4, 0, ],
    [ 1, 2, 2, 3, 2, 1, 0, ],
    [ 0, 1, 0, 1, 1, 0, 0, ],
    [ 0, 0, 0, 0, 0, 0, 0, ],
    [ 0, 0, 0, 0, 0, 0, 0, ],
];

# Dimensiones del campo muestreado
my $alto_tabla_temperaturas  = @{$temperaturas     }-1; # Quitamos 1 para evitar el caso
my $ancho_tabla_temperaturas = @{$temperaturas->[0]}-1; # del fallo en los bordes

# Dimensiones físicas del campo
my $factor_ancho = 10;
my $factor_alto  = 15;
my $ancho = $ancho_tabla_temperaturas * $factor_ancho;
my $alto  = $alto_tabla_temperaturas  * $factor_alto;

# Cálculo de la temperatura interpolada bilinealmente
sub temperatura_interpolada {
    my ( $x, $y ) = @_;         # Posición dentro del campo físico

    return 'ERROR: posición fuera de campo'
        if $x < 0 or $y < 0 or $x >= $ancho or $y >= $alto;

    # Cálculo de las coordenadas de los píxeles que rodean
    # al indicado, dentro de la tabla de muestreo
    my ($a , $b ) = ( $x / $factor_ancho, $y / $factor_alto );
    my ($x1, $y1) = ( int $a            , int $b            );
    my ($x2, $y2) = ( $x1 + 1           , $y1 + 1           );

    # Cálculo de la posición de entrada, referida dentro del cuadro
    my ($xp, $yp) = ( $a - $x1, $b - $y1 );

    # Extraemos los valores muestreados del cuadro que lo rodea
    my ( $f00, $f10, $f01, $f11 ) = (
        $temperaturas->[$y1][$x1],
        $temperaturas->[$y1][$x2],
        $temperaturas->[$y2][$x1],
        $temperaturas->[$y2][$x2],
    );

    print "Coordenadas y valores de los puntos que le rodean, en la tabla de muestreo:\n";
    print "\t[ $x1, $y1 ]:$f00 --- [ $x2, $y1 ]:$f10\n";
    print "\t[ $x1, $y2 ]:$f01 --- [ $x2, $y2 ]:$f11\n";
    print "Posición relativa dentro del cuadro: $xp,$yp\n";

    # Cálculo del valor interpolado
    my ( $b1, $b2, $b3, $b4 ) = ( $f00, $f10-$f00, $f01-$f00, $f00-$f10-$f01+$f11 );

    return $b1  +  $b2 * $xp  +  $b3 * $yp  +  $b4 * $xp * $yp;
}

# Bucle principal
print "Cálculo de temperaturas\n";
print "Dimensiones del campo: $ancho x $alto\n";
print "Dimensiones de la muestra: $ancho_tabla_temperaturas x $alto_tabla_temperaturas\n";
print "Factor: $factor_ancho, $factor_alto\n";

while ( 1 ) {

    print 'Introduzca las coordenadas a calcular (x,y): ';
    my $entrada = <>;
    chomp $entrada;
    my ($x,$y) = split q{,}, $entrada;
    $x += 0; $y += 0;
    next if !$x or !$y;

    my $t = temperatura_interpolada($x,$y);
    print "Cálculo del valor interpolado: $t\n";
    print "---------------\n";
}

__END__
Coloreado en 0.004 segundos, usando GeSHi 1.0.8.4

Una posible salida sería:
Código: Seleccionar todo
Cálculo de temperaturas
Dimensiones del campo: 60 x 105
Dimensiones de la muestra: 6 x 7
Factor: 10, 15
Introduzca las coordenadas a calcular (x,y): 37,14
Coordenadas y valores de los puntos que le rodean, en la tabla de muestreo:
        [ 3, 0 ]:0 --- [ 4, 0 ]:0
        [ 3, 1 ]:4 --- [ 4, 1 ]:6
Posición relativa dentro del cuadro: 0.7,0.933333333333333
Cálculo del valor interpolado: 5.04
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14480
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Notapor explorer » 2008-05-28 11:29 @520

Obviamente, es mucho mejor usar PDL, para hacer cosas serias...

Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
#!/usr/bin/perl
use warnings;
use strict;

use PDL;
use PDL::Image2D;

my $temperaturas = pdl([
    [ 0, 0, 0, 0, 0, 0, 0, ],
    [ 0, 5, 4, 4, 6, 5, 0, ],
    [ 1, 4, 5, 5, 7, 6, 0, ],
    [ 3, 3, 4, 5, 4, 4, 0, ],
    [ 1, 2, 2, 3, 2, 1, 0, ],
    [ 0, 1, 0, 1, 1, 0, 0, ],
    [ 0, 0, 0, 0, 0, 0, 0, ],
    [ 0, 0, 0, 0, 0, 0, 0, ],
]);
print $temperaturas;

my $superficie_real = zeroes(14,14);

bilin2d($temperaturas, $superficie_real);

print $superficie_real;

__END__
Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4

Código: Seleccionar todo
[
 [0 0 0 0 0 0 0]
 [0 5 4 4 6 5 0]
 [1 4 5 5 7 6 0]
 [3 3 4 5 4 4 0]
 [1 2 2 3 2 1 0]
 [0 1 0 1 1 0 0]
 [0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0]
]

[
 [             0              0              0              0              0              0              0              0              0              0              0              0              0              0]
 [             0      1.2426036      2.4852071      2.4852071      2.2366864      2.1538462      2.1538462      2.4023669      2.8994083       3.147929      2.8994083      2.4852071      1.2426036 -2.3912496e-15]
 [   0.076923077      2.3136095      4.5502959      4.5976331      4.2071006      4.0769231      4.0769231      4.5384615      5.4615385      5.9230769      5.4615385      4.6863905      2.3431953 -4.5092135e-15]
 [    0.61538462      2.3550296      4.0946746      4.4733728      4.5798817      4.6153846      4.6153846      5.0769231              6      6.4615385              6       5.183432       2.591716 -4.9874634e-15]
 [     1.3076923      2.4792899      3.6508876      4.2307692      4.6923077      4.8934911       4.964497      5.3550296      6.0650888       6.408284      6.0177515      5.2544379      2.6272189 -5.0557849e-15]
 [     2.3846154      2.8106509      3.2366864      3.6923077      4.1538462      4.5207101      4.8402367      4.9822485      4.9467456      4.8757396      4.7337278       4.260355      2.1301775  -4.099285e-15]
 [     2.5384615      2.6449704      2.7514793      3.0650888      3.4201183      3.8461538      4.3076923      4.3076923      3.8461538      3.5029586      3.3964497      3.0532544      1.5266272 -2.9378209e-15]
 [     1.4615385       1.816568      2.1715976      2.3195266      2.4260355      2.7692308      3.2307692      3.2307692      2.7692308      2.3431953      1.9881657      1.5621302     0.78106509 -1.5030712e-15]
 [    0.69230769      1.1538462      1.6153846      1.5739645      1.4319527      1.6923077      2.1538462      2.2248521      1.9053254      1.5384615      1.0769231     0.63905325     0.31952663 -6.1489275e-16]
 [    0.15384615     0.61538462      1.0769231     0.82840237     0.43786982     0.61538462      1.0769231      1.2721893      1.2011834              1     0.53846154     0.14201183    0.071005917 -1.3664283e-16]
 [             0     0.28402367     0.56804734     0.37869822    0.094674556     0.18934911     0.47337278     0.61538462     0.61538462     0.52071006     0.23668639              0              0              0]
 [             0    0.035502959    0.071005917    0.047337278     0.01183432    0.023668639    0.059171598    0.076923077    0.076923077    0.065088757    0.029585799              0              0              0]
 [             0              0              0              0              0              0              0              0              0              0              0              0              0              0]
 [             0              0              0              0              0              0              0              0              0              0              0              0              0              0]
]
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14480
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Notapor lis » 2008-05-28 17:25 @767

Hola, muchas gracias por responderme.

Yo descargo datos grib, los que transformo a .ctl. Estos datos tienen 27 variables (temperatura, viento, presion, etc.), entonces para trabajar con una de ellas uso el wgrib y en este caso ocupo la temperatura como ejemplo. En Perl hago lo siguiente:

Sintáxis: [ Descargar ] [ Ocultar ]
Using perl Syntax Highlighting
#!/usr/bin/perl

open (TAREA, ">tarea6.dat");
$file_save="tarea6.dat";
$filename = "/home/lis/grads/tarea6/cl_wafsgfs_P_t18z_intdsk12";

open(GUARDAR,"$file_save");
system"(wgrib $filename |grep 'kpds5=11' |grep 'kpds6=105' |grep 'kpds7=2' |wgrib -i $filename -text -o $file_save)";
close(GUARDAR);
Coloreado en 0.001 segundos, usando GeSHi 1.0.8.4


Y esto me genera un archivo llamado tarea.dat que contiene una columna de 2450 datos, que en este caso son de temperatura. Lo que tengo que hacer es ingresar una coordenada (latitud, longitud) y me tiene que mostrar los 4 números que rodean a esta coordenada y luego ingresarlos a la ecuación de interpolación lineal y que me muestre resultado de la temperatura interpolada.

Gracias por ayudarme, y estoy muy contenta por la pronta respuesta.
lis
Perlero nuevo
Perlero nuevo
 
Mensajes: 106
Registrado: 2008-05-27 21:43 @946

Notapor explorer » 2008-05-31 19:37 @859

¿2450 valores? Pero si a mi me salen 3447...

¿De qué directorio ftp del NOAA te bajas los ficheros?
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14480
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España

Notapor lis » 2008-06-04 01:16 @095

Hola explorer. Gracias por tu ayuda. Me ha sido de gran utilidad.

Mis datos eran 2450 porque yo acoto los rangos a mi medida, es decir, mi grilla es: -70° latitud sur a -8.75° latitud sur y -120° longitud oeste a -60° longitud oeste.

Gracias a ti ya resolví esa parte, ahora estoy haciendo un script en Perl con la interpolacióon de punto más cercano o vecino más próximo (nearest neighbor), y mi duda es la siguiente:

Con la interpolación bilineal (hecha en el script anterior) yo ingresaba dos datos que por ejemplo son: (x.y) = (longitud, latitud) = (-70.73, 30.25) y me muestra los 4 puntos que rodean al punto ingresado (x,y), con sus respectivas coordenadas y como yo estaba ocupando la temperatura, también me da la temperatura de esos 4 puntos y por último me da la temperatura interpolada bilinealmente.

Ahora para esa misma coordenada ingresada (-70.73, -30.25), en Perl tengo que aplicarle la interpolación del punto más cercano, y como resultado me tiene que mostrar las coordenadas del punto más cercano a (-70.73, -30.25) y la temperatura ahora interpolada por la interpolación del punto más cercano.

Recurro a ti nuevamente porque para esta nueva interpolación no hay una ecuación, y quiero saber si hay algún comando en Perl que haga esto de encontrar el punto más cercano de cualquier variable.

Cualquier ayuda te lo agradecería mucho, Explorer.
lis
Perlero nuevo
Perlero nuevo
 
Mensajes: 106
Registrado: 2008-05-27 21:43 @946

Notapor explorer » 2008-06-05 06:30 @312

Bien, tu rejilla es más pequeña, pero habrás tenido en cuenta que lo que se guarda en un grib no es una rejilla cartesiana, sino que en cada latitud se guarda un número diferente de muestras (menos cuanto más te acercas al polo). Hay que tener cuidado con eso cuando se extraen los datos.

En cuanto a lo del comando Perl, no existe algo parecido. Hay que realizar las cuentas a mano. Ya lo has visto en otro hilo.

Un detalle importante. No es lo mismo hacer esto para uno o dos puntos, que para millones de puntos. Perl no es el lenguaje más adecuado para realizar grandes procesamientos de cálculo, ya que entonces es mejor usar soluciones compiladas. En ese caso se hace imprescindible usar una librería matemática, como la PDL.
JF^D Perl programming & Raku programming. Grupo en Telegram: https://t.me/Perl_ES
Avatar de Usuario
explorer
Administrador
Administrador
 
Mensajes: 14480
Registrado: 2005-07-24 18:12 @800
Ubicación: Valladolid, España


Volver a Avanzado

¿Quién está conectado?

Usuarios navegando por este Foro: No hay usuarios registrados visitando el Foro y 2 invitados