This is something of a non-answer, alas; but also the best possible answer?
That's a bit contradictory, so let me explain.
When the resistor is thermally limited, the problem is to avoid exceeding maximum peak rated temperature at any point within the device.
Already we see a problem here, because resistors aren't rated by internal temperatures (how would you figure that out, anyway?). If we had the means to measure every point within, that would be one thing; but since we can't, what then?
The next best thing is to simulate the process -- make a thermal model, following reasonable assumptions about materials and geometry. If we can extract such a model, we can simulate any arbitrary pulse ourselves, and determine if it's within ratings or not.
It's possible that a manufacturer might possess such data. Unfortunately, I haven't come across any myself, even for familiar-packaged parts (like TO-220 thin films, which might be particularly easy and worthwhile to characterize, if anything?). So, lacking data to start with... this isn't a very good answer.
It's also possible we can take the pulse rating (energy or power) curves, and turn them into expected, or limiting, transient thermal impedance curves; but I don't have a strong feel on this and would have to do some work to come up with a method. Anyway, it's likely that such curves are themselves guardbanded (i.e. derated from ultimate device capability) by some margin, by the manufacturer already, given that they probably don't have transient thermal resistance data themselves. (But, they might, and just keep it locked away for internal use, or privileged customers; you'd have to ask.)
Using Transient Thermal Impedance
If you do have a thermal impedance network -- simply set up your power waveform (typically as a current representing power, for voltage as temperature) into the impedance network, and run the simulation.
This is a standard method for semiconductors, which generally do provide transient impedance characteristics (if not necessarily the network parameters, so you may have to do some curve fitting yourself still), and generally do have a single point at which to measure internal temperature (usually by characterizing VGS(th) or body diode VF vs. T).
Here's an exercise where I derived ZθJA for a D2PAK device:
ITH is the thermal source, current corresponding to power at the transistor; CTHn are the thermal capacitances or masses in the system, and RTHn are the resistances connecting them. Voltage corresponds to temperature, and VTA is the voltage corresponding to a temperature of ambient (which I don't care about here, so it's set to zero).
As a plotting aid, the E source buffers the junction temperature TJ, rescales it, and M1 is a log function, so that TJNL matches the scale on the datasheet plot. (I can set the time axis to display logarithmic in the viewer; I could set the same vertical functions in the viewer, actually, but this shows them off explicitly.) I overlay a capture from the datasheet using OnTopReplica, and compare simulated to measured until they match:
Tweaking the RC values is a bit tedious, but as you can see, the 4th order fit is quite good. Well, it's a bit hard to see what all is going on here, there's transparency everywhere, but the on-top window is set to 50% transparency hence the slight outline of the window, and the datasheet plot extends further leftward but I'm only modeling the rightmost curve (which is ZθJA i.e. the device on a real circuit board in ambient conditions; the other is ZθJC i.e. for the case joined to an infinite heatsink). The background matter, "Zth Sim.sdf", the vertical (V) axis, etc., is the simulation output view in Altium.
This was done to simulate a device under pulsed conditions. Here is the representative schematic:
IC and VCE are measured by E and H sources, multiplied by A1 to get voltage representing power, and converted to a current representing power into the simulation network. An integrator also records total energy (here, a bit over 6J after 400ms).
The simulation results:
Note that TJ reads in units of kelvin, temperature rise above ambient. During this pulse (an automotive load dump transient), peak power greatly exceeds the continuous rating of the device (a couple watts for a D2PAK on board with no additional heatsinking), but falls well within its headline rating (over 100W); somewhere inbetween, there is a pulse amplitude and duration that will survive, and the challenge was to determine where. As you can see, the immediate heating is strong, but not destructive, 45K -- and the application was approved for production.
In comparison, a simulated square pulse, of peak amplitude and half width (38W for 205ms), overestimates the temperature rise as 59.5 instead of 44K. A triangular approximation (10ms rise, 5ms peak dwell at 38W, 400ms linear ramp down to zero) does however get pretty close at 45.5K.
This pulse is a bit different from an exponential capacitor discharge (it's a truncated exponential), but a triangle width/duration of about 1.5 time constants should get pretty close. (Calculating the exact ratio would be a good exercise for the proverbial student. :) )
This may be of interest:
Advanced Electro-Thermal Modeling for Surface Mount Resistors, Negrea & Svasta, U.P.B. Sci. Bull., Series C, Vol. 76, Iss. 1, 2014
Appears to be a student paper, short of a thesis but nonetheless exploring the topic in some detail. Also, basically explains the "yes but no" of my answer -- quite some work is needed if the manufacturer hasn't provided it already. The references look handy as well; several Vishay documents are even included (though I don't think any of them contain transient thermal data for comparison; alas).
A relevant Infineon application note:
AN2015-10: Transient thermal measurements and thermal equivalent circuit models | Infineon
Other manufacturers provide information regarding the use or derivation of these data. (Again, few if any resistor manufacturers, sorry.)