Live Code - LIGO

import numpy as np import matplotlib.pyplot as plt from scipy.signal import hilbert from scipy.linalg import hankel, svd from gwpy.timeseries import TimeSeries # ---------------------------------------------------------------------- # 1. Core ESPRIT Function (Your Optimized Implementation) # ---------------------------------------------------------------------- def esprit_optimized(x, fs, order=2, L=None): """ ESPRIT for real damped sinusoids. Optimized for long signals with a limited Hankel matrix size. """ xa = hilbert(x) N = len(xa) if L is None: L = min(2000, N // 2) # <-- optimization: cap L if L < order: return np.nan, np.nan M = N - L + 1 H = np.zeros((L, M), dtype=complex) for i in range(L): H[i, :] = xa[i:i+M] try: U, s, Vh = svd(H, full_matrices=False) except: return np.nan, np.nan if len(s) < order: return np.nan, np.nan U_s = U[:, :order] J1 = np.eye(L-1, L, dtype=int)[:, :L-1] J2 = np.eye(L-1, L, dtype=int)[:, 1:] A = J1 @ U_s B = J2 @ U_s Psi = np.linalg.lstsq(A, B, rcond=None)[0] eigvals = np.linalg.eigvals(Psi) freqs = np.abs(np.angle(eigvals)) * fs / (2 * np.pi) dampings = -np.log(np.abs(eigvals)) * fs pos = freqs > 0 if not np.any(pos): return np.nan, np.nan return freqs[pos][0], dampings[pos][0] # ---------------------------------------------------------------------- # 2. Download and Preprocess LIGO Data (GW150914 as example) # ---------------------------------------------------------------------- # Step 2.1: Download the data from GWOSC print("Downloading LIGO data from GWOSC...") # Use the known GPS times for the GW150914 event (centered at gps 1126259462.4) start_gps = 1126259446 # 16 seconds before the event end_gps = 1126259478 # 16 seconds after the event strain = TimeSeries.fetch_open_data('H1', start_gps, end_gps) # Step 2.2: Preprocess the data fs = strain.sample_rate.value print(f"Sampling frequency: {fs:.0f} Hz") print(f"Data range: {strain.t0.value} to {strain.t0.value + strain.duration.value} seconds GPS") # Apply a bandpass filter to focus on the frequency range of interest (your invariant) low_freq = 35 high_freq = 400 filtered_strain = strain.bandpass(low_freq, high_freq) # Whiten the data to flatten the noise spectrum whitened_strain = filtered_strain.whiten(fftlength=4) # Convert the whitened GWpy TimeSeries object to a NumPy array for ESPRIT signal_array = whitened_strain.value print(f"Signal length: {len(signal_array)} samples") # ---------------------------------------------------------------------- # 3. Run ESPRIT Analysis on the LIGO Data # ---------------------------------------------------------------------- print("\nRunning ESPRIT on LIGO strain data...") f_est, d_est = esprit_optimized(signal_array, fs, order=2) if not np.isnan(f_est): print(f"ESPRIT Results on LIGO Data:") print(f" Frequency: {f_est:.4f} Hz") print(f" Damping: {d_est:.4f} s⁻¹") else: print("ESPRIT could not detect a clear dominant frequency in this segment.") # ---------------------------------------------------------------------- # 4. Visualization: Compare the PSD to the ESPRIT Result # ---------------------------------------------------------------------- fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 8)) # Top plot: Spectrum of the filtered data spec = filtered_strain.spectrum(fftlength=4) ax1.loglog(spec.frequencies.value, spec.value, label='Filtered LIGO Strain (H1)') ax1.set_xlim(low_freq, high_freq) ax1.set_ylabel('Amplitude Spectral Density (1/√Hz)') ax1.set_title('LIGO Strain Data Spectrum (GW150914)') ax1.legend() if not np.isnan(f_est): ax1.axvline(f_est, color='r', linestyle='--', label=f'ESPRIT Frequency = {f_est:.2f} Hz') ax1.legend() ax1.grid(True, alpha=0.3) # Bottom plot: The whitened time series around the event time_axis = strain.times.value - strain.t0.value # seconds relative to segment start event_time_offset = 1126259462.4 - strain.t0.value # Center of GW150914 t_window = 1.0 # seconds t_mask = np.abs(time_axis - event_time_offset) < t_window ax2.plot(time_axis[t_mask], whitened_strain[t_mask], label='Whitened LIGO Strain') ax2.set_xlabel('Time (seconds from segment start)') ax2.set_ylabel('Whitened Amplitude') ax2.set_title('Whitened LIGO Strain around GW150914') ax2.legend() ax2.grid(True, alpha=0.3) plt.tight_layout() plt.show()

Popular posts from this blog

THE GOLDEN BALLROOM/BUNKER

Image
Ben Meiselas reports on the shocking admission by Donald Trump’s DOJ in a court case where the DOJ admits to a secret project underneath the ballroom which they claim is needed to protect Donald Trump’s life for “national security purposes.” "You unlock this door with the key of complicity. Beyond it is another dimension — a dimension of betrayal, of indulgence, of fear. You’re moving into a land of both shadow and substance, of politics and paranoia. You’ve just crossed into… the MAGA Zone." "Tonight’s story: A leader sworn to protect his nation makes a bargain with its enemies. The deal? Silence in the face of nuclear annihilation. No retaliation, no defense — only surrender dressed in secrecy. While citizens live unaware, their president builds a palace beneath the earth, a ballroom of gold, of marble and chandeliers, a masquerade hall for billionaires. But behind the gilded doors lies not music and laughter, but a bomb shelter — a sanctuary for the few, pur...
Read more

Conceptual Summary #2: (∂t2​S−c2∇2S+βS3)=σ(x,t)⋅FR​(C[Ψ])

Image
The tests designed to find the Aether—principally the Michelson-Morley experiment (M-M)—would be largely ineffective at detecting the fundamental S field (Substrate) in the Reactive Substrate Theory (RST) framework, though their null results are highly compatible with RST's structure. The reason lies in the distinct conceptual definition and behavior of the two media. The Test: Michelson-Morley and the Aether. The M-M experiment was designed to detect the "aether wind," a predicted change in the speed of light caused by the Earth's motion relative to a stationary, rigid medium—the classical Luminiferous Aether. Classical Aether: A static, absolute reference frame; an elastic, passive, material-like medium that fills space and serves as a carrier for light waves. Predicted Result: Light traveling parallel to Earth's motion through the Aether should be slightly slower than light traveling perpendicular to it. Actual Result (Null Result): No significant differenc...
Read more

ICE PROUDLY ANNOUNCES NEW “ELITE” TASK FORCE COMMANDER JEREMY DEWITTE

Image
For Immediate Release OFFICE OF STRATEGIC SUPPRESSION & DESPERATE RECRUITMENT U.S. IMMIGRATION AND CUSTOMS ENFORCEMENT (ICE) SUBJECT: ICE PROUDLY ANNOUNCES NEW “ELITE” TASK FORCE COMMANDER JEREMY DEWITTE; AGENCY TO ADOPT FULL “MOTOR ESCORT” DOCTRINE BY Q3 WASHINGTON, D.C. — In a historic moment that experts are calling “a cry for help” and insiders are calling “Tuesday,” U.S. Immigration and Customs Enforcement (ICE) has officially appointed Jeremy Dewitte as National Coordinator of the newly created Tactical High-Visibility Fallujah Squadron & Escort Wing (T.H.V.F.S.E.W.) . This appointment follows a PBS NewsHour report highlighting ICE’s bold new hiring strategy: removing every meaningful hiring standard. After eliminating outdated requirements such as "must not be a professional police impersonator," "must possess a valid driver's license," and "must not be a registered sex ...
Read more