91 lines
2.6 KiB
PHP
91 lines
2.6 KiB
PHP
<?php
|
|
ini_set('display_errors', 1);
|
|
ini_set('display_startup_errors', 1);
|
|
error_reporting(E_ALL);
|
|
|
|
header('Content-Type: application/json');
|
|
|
|
// --- Configuration ---
|
|
$gfs_url = 'https://nomads.ncep.noaa.gov/pub/data/nccf/com/gfs/prod/gfs.20251013/00/atmos/gfs.t00z.pgrb2.0p25.f000';
|
|
$tmp_dir = '/tmp';
|
|
$grib_file = $tmp_dir . '/' . basename($gfs_url);
|
|
$xyz_file = $tmp_dir . '/gfs.xyz';
|
|
|
|
// --- Download the GFS file ---
|
|
if (!file_exists($grib_file)) {
|
|
$fp = fopen($grib_file, 'w');
|
|
$ch = curl_init();
|
|
curl_setopt($ch, CURLOPT_URL, $gfs_url);
|
|
curl_setopt($ch, CURLOPT_FILE, $fp);
|
|
curl_setopt($ch, CURLOPT_FOLLOWLOCATION, 1);
|
|
curl_exec($ch);
|
|
$http_code = curl_getinfo($ch, CURLINFO_HTTP_CODE);
|
|
curl_close($ch);
|
|
fclose($fp);
|
|
|
|
if ($http_code !== 200) {
|
|
unlink($grib_file); // Clean up failed download
|
|
echo json_encode(['error' => 'Failed to download GFS file', 'http_code' => $http_code]);
|
|
exit;
|
|
}
|
|
}
|
|
|
|
// --- Convert GRIB to XYZ using gdal_translate ---
|
|
// We'll use band 1 for this example. You can use gdalinfo to see available bands.
|
|
$gdal_command = 'gdal_translate -b 1 -of XYZ ' . escapeshellarg($grib_file) . ' ' . escapeshellarg($xyz_file);
|
|
$gdal_output = shell_exec($gdal_command);
|
|
|
|
if (!file_exists($xyz_file)) {
|
|
echo json_encode(['error' => 'Failed to convert GRIB to XYZ', 'gdal_output' => $gdal_output]);
|
|
unlink($grib_file);
|
|
exit;
|
|
}
|
|
|
|
// --- Read the XYZ file and create GeoJSON features ---
|
|
$features = [];
|
|
$handle = fopen($xyz_file, 'r');
|
|
$line_count = 0;
|
|
$sample_rate = 100; // Process 1 in every 100 lines
|
|
|
|
if ($handle) {
|
|
while (($line = fgets($handle)) !== false) {
|
|
$line_count++;
|
|
if ($line_count % $sample_rate !== 0) {
|
|
continue;
|
|
}
|
|
|
|
$parts = preg_split('/\s+/', trim($line));
|
|
if (count($parts) === 3) {
|
|
$lon = floatval($parts[0]);
|
|
$lat = floatval($parts[1]);
|
|
$value = floatval($parts[2]);
|
|
|
|
// Skip points that are exactly zero, often represent no data
|
|
if ($value === 0.0) {
|
|
continue;
|
|
}
|
|
|
|
$features[] = [
|
|
'type' => 'Feature',
|
|
'properties' => ['value' => $value],
|
|
'geometry' => [
|
|
'type' => 'Point',
|
|
'coordinates' => [$lon, $lat]
|
|
]
|
|
];
|
|
}
|
|
}
|
|
fclose($handle);
|
|
}
|
|
|
|
// --- Clean up temporary files ---
|
|
unlink($grib_file);
|
|
unlink($xyz_file);
|
|
|
|
// --- Output the GeoJSON ---
|
|
echo json_encode([
|
|
'type' => 'FeatureCollection',
|
|
'features' => $features
|
|
]);
|
|
|
|
?>
|